Metodos Numericos En Polymath 5n3z6y

UNIVERSIDAD NACIONAL DEL CALLAO VICE- RECTORADO DE INVESTIGACION

INSTITUTO DE INVESTIGACION DE LA FACULTAD DE INGENIERIA QUIMICA

INFORME FINAL DEL PROYECTO DE INVESTIGACION TITULADO

TEXTO " METODOS NUMERICOS UTILIZANDO POLYMATH, MATHCAD Y MATLAB APLICADOS A INGENIERIA QUIMICA” PRESENTADO POR

ING. JUAN MEDINA COLLANA

(Del 1 de Marzo del 2010 al 29 de Febrero 2012 Resol. N° 1312-2010 R)

1

ÍNDICE I.

RESUMEN

4

II.

INTRODUCCIÓN

5

III.

MARCO TEÓRICO

8

IV.

MATERIALES Y MÉTODOS

8

V.

RESULTADOS

9

VI.

DISCUSIÓN

9

VII.

REFERENCIALES

10

1. ECUACIONES NO LINEALES 1.1 Método de Newton-Raphson 1.2 Método de la secante 1.3 Método de la bisección 1.4 Método de la regla falsa 1.5 Método de iteración del punto fijo 1.6 Problemas de aplicación a la Ingeniería Química 1.7 Aplicaciones en Ingeniería Química

11 11 14 16 19 22 24 32

2. SISTEMAS DE ECUACIONES NO LINEALES 2.1 Método de newton Raphson 2.1.1Aplicación 2.2 Problemas

36 36 37 40

3. INTERPOLACIÓN 3.1 Polinomios de interpolación de newton 3.1.1 Interpolación lineal 3.1.2 Interpolación cuadrática 3.1.3 Diferencias finitas divididas 3.2 Polinomios de interpolación de Lagrange 3.3 Problemas

44 45 46 47 47 48 50

4. REGRESIÓN 4.1 Regresión lineal 4.2 Regresión polinomial 4.3 Regresión lineal múltiple 4.4 Problemas

53 53 54 55 57

5. DIFERENCIACIÓN NUMÉRICA 5.1 Diferenciación mediante método Newton 5.2 Diferenciación de Lagrange: datos discretos

68 67 67

2

5.3

Problemas

69

6. INTEGRACIÓN NUMÉRICA 6.1 Método del trapezoide 6.2 Regla de Simpson 6.2.1 regla de Simpson 1/3 6.2.2Simphson 3/8 6.3 Problemas

76

7. SISTEMA DE ECUACIONES LINEALES 7.1 Métodos de Jacobi 7.2 Métodode Gauss – Seidel 7.3 Problemas

85 85 87 88

8. ECUACIONES DIFERENCIALES NUMÉRICAS 8.1 Método de Euler 8.2 Método de Euler modificado 8.3 Métodos de Runge-Kutta (rk) 8.4 Método de diferencias finitas 8.5 Sistema de ecuaciones diferenciales de primer orden 8.6 Problemas

91 92 94 96 100 101 103

9. ECUACIONES DIFERENCIALES PARCIALES (EDP)

110

9.1

Problemas

110

10.

APLICACIONES DE INGENIERÍA QUÍMICA EN

76

77 78 80 81

POLYMATH- MATHCAD Y MATLAB 10.1 Problemas con polymath 10.2 Problemas con mathcad 10.3 Problemas con matlab

112 112 121 131

SILABO

147

3

I.

RESUMEN

La presente Investigación tuvo como propósito la elaboración de un texto universitario titulado “TEXTO: " METODOS NUMERICOS UTILIZANDO POLYMATH, MATHCAD Y MATLAB APLICADOS A INGENIERIA QUIMICA” Se trata de un texto básico que se expone brevemente los fundamentos teóricos, ilustraciones

con problemas

resueltos

y propuestos de

ingeniería química como, equilibrio químico, ecuaciones de estado, transferencia de calor, cinética química y al final se plantea problemas resueltos haciendo uso de software de polymath, mathcad y matlab.. Asi mismo en cada capitulo se prantea problemas propuestos aplicados a la ingeniería química .

4

II.

INTRODUCCIÓN.

El tema de investigación, está referida a la elaboración de un texto Universitario, cuyo propósito es apoyar en la labor de formación de los alumnos, en el curso de métodos numéricos aplicada a ingeniería .Quimica en la universidad Nacional del Callao facultad de Ingeniería Química. Durante mi experiencia en la docencia, en el intento de encontrar textos necesarios para el dictado del curso de

métodos numéricos , se ha

podido constatar que existe poca bibliografía en nuestro medio haciendo uso de mathcad,

polymath y matlab en

un solo texto, mediante ste

trabajo se ha desarrrollado problemas aplicados a la ingeniería química haciendo uso de una técnica numérica, al mismo timpo se ha resuelto problemas con el software mathcad, polymath y mathcad, que se vienen usando con mayor intesnsidad en los últimos años a nivel de ingeniería. En este trabajo se ha resaltado el capitulo de ecuaciones no lineales ecuaciones diferenciales puesto

y

que gran parte de los modelos de

ingenieira química se encuentran en esta area.

5

2.1

A.

PLANTEAMIENTO DEL PROBLEMA DE INVESTIGACION

DESCRIPCIÓN Y ANÁLISIS DEL TEMA El presente trabajo de investigación es una propuesta para la

elaboración de un texto Universitario titulado: TEXTO:

“METODOS

NUMERICOS

UTILIZANDO

POLYMATH,

MATHCAD Y MATLAB APLICADOS A LA INGENIERIA QUIMICA.”

Dirigido a estudiantes de pre – grado de Ingeniería Química y otras especialidades afines, que presente de una manera didáctica los principios fundamentales y el uso adecuado de los software de polymath, mathcad y matlab, lo que permitirá cumplir con los propósitos de una adecuada enseñanza y formación profesional.

El texto contendrá una base de teoría apropiada y práctica que van a permitir desarrollar criterios y habilidades, que resultara muy valioso para los propósitos de este texto.

B.

PLANTEAMIENTO DEL PROBLEMA

¿Como

elaborar

un

TEXTO:

“METODOS

NUMERICOS

UTILIZANDO POLYMATH, MATHCAD Y MATLAB APLICADOS A LA INGENIERIA QUIMICA.”, que oriente adecuadamente a los estudiantes de Ingeniería Química?

2.2

2.2.1

OBJETIVO Y ALCANCE DE LA INVESTIGACION

OBJETIVO GENERAL

Desarrollar

un

TEXTO:

“METODOS

NUMERICOS

UTILIZANDO

POLYMATH, MATHCAD Y MATLAB APLICADOS A LA INGENIERIA QUIMICA.”

6

2.2.2.

1.

OBJETIVOS ESPECIFICOS

Recopilar Información básica y actualizada, necesaria para iniciar el desarrollo del texto.

2.

Analizar y procesar la información para iniciar el desarrollo del texto.

3.

Desarrollar los capítulos del texto referido a los fundamentos teóricos.

4.

Desarrollar los capítulos del texto referido a las aplicaciones prácticas aplicada la ingeniería química

5.

Desarrollar los capítulos del texto referido a las aplicaciones prácticas aplicada la ingeniería química haciendo uso del uso polymath , mathcad y matlab

2.3.

ALCANCE DE LA INVESTIGACION

El presente trabajo de investigación de acuerdo a la naturaleza del problema se puede

manifestar que es una investigación básica y

Aplicada, dado los modelos planteados para la solución numérica proviene de fenómenos químicos, cuya solución analítica sería demasiada compleja. El aporte del trabajo de investigación

estará orientado al sector

académico conformado por los profesores, estudiantes y egresados de la Facultad de ingeniería química de la Universidad Nacional del Callao y otras Universidades del país. Por otro lado, este texto también podría ser utilizado por estudiantes de especialidades afines tales como ingeniería de Alimentos, ingeniería ambiental e ingeniería Industrial.

2.4.

IMPORTANCIA Y JUSTIFICACION DE LA INVESTIGACION

2.4.1 IMPORTANCIA Al desarrollar el texto propuesto se facilitará el proceso de enseñanza – aprendizaje en la formación profesional del estudiante universitario

a nivel de pre-grado

y pueda facilitar los cálculos

laboriosos de ingeniería haciendo uso de software. 2.4.2 JUSTIFICACION 7

La contribución del presente trabajo estará orientada a la preparación y el entrenamiento de los alumnos de ingeniería química en el curso de Métodos numéricos, adquiriendo fundamentos teóricos y la parte práctica que consiste en efectuar cálculos de balances de materia, energía, termodinámica, reacciones química,

transferencia de calor

entre otros.

III.

MARCO TEÓRICO En la presente investigación se ha incorporado la teoría resumida y simplificada para nueve capítulos del presente texto. Así por ejemplo, en el capítulo I se hace referencia a los diferentes métodos de

solución de ecuaciones no lineales con

sus respectivos

ejemplos y problemas. En el capítulo

III , se hace referencia de las técnicas de

interpolación , resaltando el método de diferencias por newton . Asimismo en el capítulo VIII de ecuaciones diferenciales , se hace referencia de los diferentes ordenes de solución de Runge Kutta y con sus respectivos ejemplos. En el capítulo X se presenta problemas resueltos con polymath , mathcad y Matlab. IV.

MATERIALES Y MÉTODOS Materiales: 

Materiales De oficina



Material bibliográfico



Software Polymath



Software Mathcad



Software Matlab



Material de cómputo e impresión

Métodos La elaboración del texto, propósito de la investigación le demandó al suscrito, ordenar la información disponible y complicada en función del Syllabus propuesto del curso de métodos numéricos . 8

La estructuración del texto responde a la experiencia docente en la Facultad de Ingeniería Química de la Universidad Nacional del Callao.

Para la elaboración del texto, se tuvo cuidado en recurrir a la síntesis de los aspectos teóricos,

selección de los problemas

de

ingeneiria química, elaboración de los programas en matnab, asi como revisión de tutoriales de mathcad.

En cuando al planteamiento de problemas, se quiere plasmar la experiencia con el dictado del curso, que me ha permitido revisar una extensa bibliografía sobre la materia, así como volcar en el texto los resultados de mi ejercicio profesional en el campo de Ingeniería Química para satisfacer el propósito de la investigación.

V.

RESULTADOS El resultado de la presente investigación es la elaboración de un texto

universitario

titulado:

TEXTO:

“METODOS

NUMERICOS

UTILIZANDO POLYMATH, MATHCAD Y MATLAB APLICADOS A LA INGENIERIA QUIMICA.”, , el cual se adjunta al presente. El texto contiene 9 capítulos. La teoría desarrollada en el texto, responde a los aspectos básicos de métodos numericos. Los problemas resueltos en el texto, tienen el propósito de dar las pautas de la aplicación de la teoría desarrollada. Se ha logrado un texto base para el curso de métodos numéricos necesario

en la formación universitaria del estudiante de Ingeniería

Química. VI.

DISCUSIÓN El texto universitario titulado ““METODOS NUMERICOS UTILIZANDO POLYMATH, MATHCAD Y MATLAB APLICADOS A LA INGENIERIA QUIMICA.”, aplicados a la Ingeniería Química, que es el resultado de la investigación a que se refiere el presente informe, se caracteriza por presentar

la metodología de calculo de los modelos . Los problemas

resueltos y planteados han sido seleccionados con la intención 9

de

brindar

una mayor claridad

a

los alumnos

y puedan entender los

fundamentos teóricos tratados. Los textos de METODOS NUMERICOS

contienen demasiada

información, muchas veces muy detallado, sin conexión directa con aplicación inmediata. Por eso, el presente texto va a tratar de desenvolver el contenido de tal manera que cada capítulo describa en forma precisa y concreta la teoría y los problemas de aplicación. Sin embargo,

y por ende no

sustituye el uso de la bibliografía de la especialidad, a la que deberá referirse necesariamente quienes pretendan profundizar en conocimientos de temas específicos.

VII. REFERENCIALES

Existen textos que se ocupan de los métodos numéricos en general entre los cuales tenemos 1. A. Constantinides and N. Mosotoufi, Numerical Methods for Chemical Engineers with MATLAB Applications Prentice Hall , Upper Saddle River, 1999. 2. Burden, R. Y Faires J. Análisis Numérico. Edit. Iberoamericana, México, 1985. 3. Carnahan, B. Luther, A. Wilkes Cálculo Numérico, Aplicaciones Editorial Rueda, Madrid, 1979 4. Carrasco Venegas Luis . Métodos Numéricos aplicados a la ingeniera 5. Nieves, A., Domínguez, F. Métodos Numéricos Aplicados a la Ingeniería Química Edit. CECSA, México 1985. 6. Nakamura, S. Métodos Numéricos aplicados con

software. Edit.

Prentice – Hall Hispano Americano, S.A. México, 1992. 7. S.C. Chapra y R. P. Canale,. Métodos Numéricos para Ingenieros. McGraw Hill, México,1999

10

1.ECUACIONES NO LINEALES Los métodos numéricos de resolución de ecuaciones no lineales suelen ser métodos iterativos que producen una sucesión de valores aproximados de la solución, que se espera, que converja a la raíz de la ecuación. Estos métodos van calculando las sucesivas aproximaciones en base a los anteriores, a partir de una o varias aproximaciones iníciales.. Para saber que método debemos aplicar, hay que tener en cuenta la capacidad de separar raíces cercanas, confiabilidad en el alcance de soluciones evitando errores numéricos graves y orden de convergencia. Uno de los problemas que con mayor frecuencia aparece en la ciencia y en la ingeniería es hallar las raíces de una ecuación no lineal de la forma f(x) = 0. Estudiaremos Métodos Iterativos para determinar aproximaciones a raíces reales simples de la ecuación no lineal f(x) = 0.

1.1

MÉTODO DE NEWTON-RAPHSON

Este método, es uno de los más usados, a diferencia de los otros métodos, el método de Newton-Raphson no trabaja sobre un intervalo sino que basa su fórmula en un proceso iterativo. Supongamos que tenemos la aproximación xi a la raíz xr de f ( x ) ,

Figura Nº 1: Demostración método de newton

11

Trazamos la recta tangente a la curva en el punto  xi , f ( xi )  ; ésta cruza al eje x en un punto xi 1 que será nuestra siguiente aproximación a la raíz xr . Para calcular el punto xi 1 , calculamos primero la ecuación de la recta tangente. Sabemos que tiene pendiente m  f   xi 

Y por lo tanto la ecuación de la recta tangente es: y  f  xi   f   xi  x  xi 

Hacemos y  0 :  f  xi   f   xi  x  xi 

Y despejamos x :

x  xi 

f  xi   x  xi 

Que es la fórmula iterativa de Newton-Raphson para calcular la siguiente aproximación:

xi 1  xi 

f  xi  , si f   xi   0 f   xi 

Note que el método de Newton-Raphson no trabaja con intervalos donde nos asegure que encontraremos la raíz, y de hecho no tenemos ninguna garantía de que nos aproximaremos a dicha raíz. Desde luego, existen ejemplos donde este método no converge a la raíz, en cuyo caso se dice que el método diverge. Sin embargo, en los casos donde si converge a la raíz lo hace con una rapidez impresionante, por lo cual es uno de los métodos preferidos por excelencia. También observe que en el caso de que f   xi   0 , el método no se puede aplicar. De hecho, vemos geométricamente que esto significa que la recta tangente es horizontal y por lo tanto no intersecta al eje en ningún punto, a menos que coincida con éste, en cuyo caso xi mismo es una raíz de f ( xi )! Ejemplo Usar el método de Newton-Raphson, para aproximar la raíz de

f ( x)  e x  ln x , comenzando con x0  1 y hasta que a  1% .

12

Solución En este caso, tenemos que f ( x )   e  x 

1 x

De aquí tenemos que:

Iniciamos con x0  1 y obtenemos:

En este caso, el error aproximado es,

Continuamos el proceso hasta reducir el error aproximado hasta donde se pidió. Resumimos los resultados en la siguiente tabla: Tabla Nº1 Aproximación de la raíz y porcentaje de error Aprox. a la raíz

Error aprox.

1 1.268941421

21.19%

1.309108403

3.06%

1.309799389

0.052%

De lo cual concluimos que la aproximación obtenida es:

13

1.2

MÉTODO DE LA SECANTE

Uno de los objetivos de este método es eliminar el problema de la derivada de la función, ya que existen funciones que describen fenómenos físicos y químicos, cuya derivada es muy compleja El método de la secante es muy similar al de Newton con la diferencia principal que en este método de la secante no requiere de la derivada. Este método funciona por medio de dos puntos sobre el eje x, es decir un intervalo (xi-1,xi), los cuales se evalúan en la función para sacar los puntos correspondientes en el eje de la y, los puntos a obtener son f(x i-1) y f(xi), por lo que las coordenadas de los puntos que interceptan a la función son (x i1,f(xi-1))

y el (xi ,f(xi)).

Se debe considerar que los puntos xi-1 y xi deben de contener a la raíz, por lo que el punto xi-1 debe estar a la izquierda y el punto xi a la derecha de la raíz. Estos dos puntos que interceptan a la función, se unen por medio de una recta, la cual al cruzar el eje de la x, genera el siguiente punto de acercamiento xi+1 , el cual quedara ubicado entre el intervalo propuesto, como se muestra en la Figura. El método de la secante parte de dos puntos (y no sólo uno como el método de Newton) y estima la tangente (es decir, la pendiente de la recta) por una aproximación de acuerdo con la expresión: Este método se basa en la fórmula de Newton-Raphson, pero evita el cálculo de la derivada usando la siguiente aproximación:

f   xi  

f  xi1   f  xi  xi1  xi

(Recuérdese la solución numérica al problema del cuerpo en caída libre). Sustituyendo en la fórmula de Newton-Raphson, obtenemos:

xi1  xi 

f  xi  f  xi   xi  f  xi1   f  xi  f   xi  xi1  xi

14

Ejemplo Usar el método de la secante para aproximar la raíz de f  x   e x  x , 2

comenzando con x 0  0 , x 0  1 y hasta que a  1% . Solución Tenemos que f  x0   1 y f  x1   0, 632120558 , que sustituimos en la fórmula de la secante para calcular la aproximación x 2 :

Con un error aproximado de:

Como todavía no se logra el objetivo, continuamos con el proceso. Resumimos los resultados en la siguiente tabla: Tabla Nº2 Aproximación de la raíz y porcentaje de error

Aprox. a la raíz

Error aprox.

0 1

100%

0.612699837

63.2%

0.653442133

6.23%

0.652917265

0.08%

De lo cual concluimos que la aproximación a la raíz es:

x4  0,652917265

15

1.3

MÉTODO DE LA BISECCIÓN

Este método es basado en el teorema de Bolzano, que establece que si una función continua cambia de signo en el intervalo (a,b), es decir, f(a)f(b)<0, entonces, existe al menos una raíz α, α(a,b). El método de bisección se basa en el siguiente teorema de Cálculo: El método de bisección sigue los siguientes pasos: Sea f  x  continua, i)

Encontrar valores iniciales x a , x b tales que f  xa 

y f  xb 

tienen

signos opuestos, es decir, f  xa   f  xb   0

ii)

La primera aproximación a la raíz se toma igual al punto medio entre x a y x b : xr 

iii)

Evaluar

xa  xb 2

f  xr  . Forzosamente debemos caer en uno de los

siguientes casos: 

f  xa   f  xr   0

En este caso,

tenemos que f  xa 

y

f  xr 

tienen signos

opuestos, y por lo tanto la raíz se encuentra en el intervalo  x a ,xb  

f  xa   f  xr   0

En este caso, tenemos que

f  xa 

signo, y de aquí que f  xr 

f  xb  tienen signos opuestos. Por

y

y

lo tanto, la raíz se encuentra en el intervalo

f  xr  tienen el mismo

 xr

, xb  .

f  xa   f  xr   0

En este caso se tiene que f  xr   0 y por lo tanto ya localizamos la raíz. El proceso se vuelve a repetir con el nuevo intervalo, hasta que:

16

a s

es decir, xactual  x previa xactual

 100% s

Ejemplo Aproximar la raíz de f  x   e  x  ln x

hasta que a  1% .

Solución Sabemos por lo visto en el ejemplo 1 de la sección anterior, que la única raíz f  x

de

se localiza en el intervalo 1,1, 5  .

Así que este intervalo es

nuestro punto de partida; sin embargo, para poder aplicar el método de bisección debemos checar que f 1

y f 1, 5  tengan signos opuestos.

En efecto, tenemos que f 1  e 1  ln1  e 1

mientras que f 1, 5   e 1.5  ln(1, 5)  0,18233  0

Cabe mencionar que la función f  x  sí es continua en el intervalo 1,1, 5  . Asípues, tenemos todos los requisitos satisfechos para poder aplicar el método de bisección. Comenzamos: i)

Calculamos el punto medio (que es de hecho nuestra primera aproximación a la raíz): xr1 

1  1, 5  1, 25 2

ii)

Evaluamos f 1, 25   e 1,25  ln(1, 25)  0, 0636  0

iii)

Para identificar mejor en que nuevo intervalo se encuentra la raíz, hacemos la siguiente tabla:

Por lo tanto, vemos que la raíz se encuentra en el intervalo 1.25,1.5  .

17

En este punto, vemos que todavía no podemos calcular ningún error aproximado, puesto que solamente tenemos la primera aproximación. Así, repetimos el proceso con el nuevo intervalo 1.25,1.5  . Calculamos el punto medio (que es nuestra segunda aproximación a la raíz): xr2 

1, 25  1,5  1,375 2

Aquí podemos calcular el primer error aproximado, puesto que contamos ya con la aproximación actual y la aproximación previa:

Puesto que no se ha logrado el objetivo, continuamos con el proceso. Evaluamos f 1, 375   e 1,375  ln(1, 375)  0, 06561  0 , y hacemos la tabla:

Así, vemos que la raíz se encuentra en el intervalo 1.25,1.375  . Calculamos el punto medio,

Y calculamos el nuevo error aproximado:

Resumimos los resultados que se obtienen en la siguiente tabla: Tabla Nº3 Aproximación de la raíz y porcentaje de error Aprox. a la raíz

Error aprox.

1.25 1.375

9.09%

1.3125

4.76%

1.28125

2.43%

1.296875

1.20%

1.3046875

0.59% 18

1.4

MÉTODO DE LA REGLA FALSA

Como mencionamos anteriormente, sería bueno considerar si la raíz de una ecuación está localizada más cerca de alguno de los extremos del intervalo. Consideremos nuevamente una gráfica Como la anterior,

Figura Nº2: Demostración de el método regla falsa Donde hemos agregado la línea recta que une los puntos extremos de la gráfica en el intervalo  a , b  . Es claro que si en lugar de considerar el punto medio del intervalo, tomamos el punto donde cruza al ejexesta recta, nos aproximaremos mucho más rápido a la raíz; ésta es en sí, la idea central del método de la regla falsa y ésta es realmente la única diferencia con el método de bisección, puesto que en todo lo demás los dos métodos son prácticamente idénticos. Supongamos que tenemos una función f ( x ) que es continua en el intervalo

 xa , xb 

y además, f ( x a ) y f ( xb ) tienen signos opuestos.

Calculemos

x

a

la

ecuación

de

la

línea

recta

que

une

los

puntos

, f  xa   ,  xb , f  xb   . Sabemos que la pendiente de esta recta está dada

por:

Por lo tanto la ecuación de la recta es:

19

Para obtener el cruce con el ejex, hacemosy = 0:

Multiplicando por xb  xa nos da:

Finalmente, de aquí despejamos x :

Este punto es el que toma el papel de xr en lugar del punto medio del método de bisección. Así pues, el método de la regla falsa sigue los siguientes pasos: Sea f  x  continua, i)

Encontrar valores iníciales xa , xb tales que f  xa  y f  xb 

tienen

signos opuestos, es decir, f  xa   f  xb   0

ii)

La primera aproximación a la raíz se toma igual a:

 f  xa  xb  xa   xr  xa     f  xb   f  xa   iii)

Evaluar

f  xr  . Forzosamente debemos caer en uno de los

siguientes casos: 

f  xa   f  xr   0

En este caso, tenemos que f  xa  y f  xr  tienen signos opuestos, y por lo tanto la raíz se encuentra en el intervalo  xa , xr  . 

f  xa   f  xr   0


f  xa  y f  xr  tienen el mismo signo, y de aquí

que f  xr  y f  xb  tienen signos opuestos. Por lo tanto, la raíz se encuentra en el intervalo  xr , xb  . 

f  xa   f  xr   0

En este caso se tiene que f  xr   0 y por lo tanto ya localizamos la raíz. 20

El proceso se vuelve a repetir con el nuevo intervalo, hasta que: a s

Ejemplo Usar el método de la regla falsa para aproximar la raíz de f  x   e  x  ln x , comenzando en el intervalo 1, 2 y hasta que ea  1% . Solución Este es el mismo ejemplo 1 del método de la bisección. Así pues, ya sabemos que f  x  es continua en el intervalo dado y que toma signos opuestos en los extremos de dicho intervalo. Por lo tanto podemos aplicar el método de la regla falsa. Calculamos la primera aproximación:

Puesto que solamente tenemos una aproximación, debemos seguir con el proceso. Así pues, evaluamos Y hacemos nuestra tabla de signos:

De donde vemos que la raíz se encuentra en el intervalo

1,1.397410482 .

Con este nuevo intervalo, calculamos la nueva aproximación:

En este momento, podemos calcular el primer error aproximado:

Puesto que no se cumple el objetivo seguimos con el proceso. Evaluamos f  xr 2   f 1, 321130513   0, 011654346  0 , y hacemos la tabla de signos: 21

De donde vemos que la raíz se encuentra en el intervalo 1,1.321130513 , con el cual, podemos calcular la nueva aproximación:

Y el error aproximado:

Como se ha cumplido el objetivo, concluimos que la aproximación buscada es: x r 3  1, 311269556

Observe la rapidez con la cual converge el método de la regla falsa a la raíz, a diferencia de la lentitud del método de la bisección.

1.5

MÉTODO DE ITERACIÓN DEL PUNTO FIJO

Este método se aplica para resolver ecuaciones de la forma x  g ( x)

Si la ecuación es

f ( x)  0 , entonces puede despejarse x ó bien sumar x

en ambos lados de la ecuación para ponerla en la forma adecuada.

Ejemplo Usar el método de iteración del punto fijo para aproximar la raíz de , comenzando con

y hasta que

.

Solución Si despejamos la xdel término lineal, vemos que la ecuación equivale a

de donde,

22


.

g ( x )  1 , para x    1,1 lo que es suficiente para deducir que el método sí

converge a la raíz buscada. Aplicando la fórmula iterativa, tenemos: x1  g ( x 0 )   0, 2

Con un error aproximado del 100%. Aplicando nuevamente la fórmula iterativa, tenemos: x 2  g ( x1 )   0,1557461506

Con un error aproximado igual al 28.41%. En este ejemplo, el método solo necesita de 5 iteraciones para reducir el error menor al 1%. Resumimos los resultados en la siguiente tabla:

Tabla Nº4 Aproximación de la raíz y porcentaje de error

Aprox. a la raíz

Error aprox.

0 -0.2

100%

-0.1557461506

28.41%

-0.1663039075

6.34%

-0.163826372

1.51%

-0.164410064

0.35%

De donde vemos que la aproximación buscada es: x5   0,164410064

Use el método de punto fijo para resolver f ( x)  x2  2x  3 en el intervalode : x=-1 y

x=3

f ( x)  x2  2x  3  x  x 2  3  g ( x)  x 2  3 23

x

3 3  g ( x)  x2 x2

x

x2  3 x2  3  g ( x)  2 2

x  x 2  x  3  g ( x)  x2  x  3

1.6 PROBLEMAS DE APLICACIÓN A LA INGENIERÍA QUÍMICA

1. A un reactor ingresa una mezcla de gases de 8% de dióxido de azufre y el 12% de oxigeno y 80%

de nitrógeno, y se desarrolla la siguiente

reacción.

1 SO2( g )  O2( g )  SO3( g ) 2 Calcular la composición en el equilibrio a presión constante de 2 atm y la constante de equilibrio KP es de 160.atm1/2 Solución : Asumiendo 100 moles de mezcla SO2 = O2

8 moles

=

12 moles

SO3 = N2

0

=

80 moles

Moles enl equilibrio SO2 = O2 SO3=

8 -

= 12

x - 0,5x

x

Moles totales = 100

-

0,5x

Luego hallamos el Kp nSO3 x ( ) Pt ( ) PSO 3 nt 100  0.5 x Kp    nSO2 nO PSO2 x ( PO2 )1/ 2 ( ) Pt ( 2 )1/ 2 Pt 1/ 2 ( 8  x ) ( 12  0.5 x ) Pt nt nt 100  0.5 x 100  0.5 x Para Kp=160 y Pt=2atm

La ecuación queda de la siguiente manera:

24

f ( x)  160 

x(100  0.5 x)0.5 (8  x)(12  0.5 x)0.5 2

Utilizando el método de la bisección y tomando como referencia los valores de : X1=7,87 y X2=7,88 Siguiendo el proceso iterativo se tiene los siguientes valores: Tabla Nº5 Aproximación de la raíz X1

F(x1)

X2

F(x2)

x

F x

7.87

12.604

1.88

-0.3001

7.875

6.2314

7.875

6.2314

7.88

--0.3001

7.8775

3.0323

7.8775

3.0323

7.88

-0.3001

7.87875

1.3833

7.87875

1.3833

7.88

-0.3001

7.879375

0.5460

Entonces X=7.878375 moles La composición en el equilibrio seria Moles SO2 Moles O2 Moles SO3

=

0.120625

= =

8.0606125 7.879375

2. La ecuación de estado Redlich- Kwong es :

  a p  (V  b)  RT T . V ( V  b )   Donde a = 17,19344 y b = 0,0221141 para el oxigeno molecular si T = 373 K y P= 30 atm Calcular el volumen molar por el método de la secante Solución: Utilizando como referencia la ecuación del gas ideal para obtener el primer valor 25

V 

RT 0, 082  373   1.0195 P 30

Tomando dos valores Xi-1 = 1

; f(Xi-1 )

Xi

; f(Xi)

=1,5

= -0,3977 = 14,3268

Remplazando en la ecuación

14,3268(1,5  1) 14,3268  (0,3977) X (i  1)  1, 0135 f ( Xi  1)  305520 E  3 X (i  1)  1,5 

Haciendo ahora : X(i-1)=1,5 y

X(i)=1,0135

Remplazando X (i  1)  1, 0135 

3, 5520  10 3 (1, 0135  1, 5)  1, 0136 3, 5520  10 3  14, 3268

X (i  1)  1, 0136

 X (i  1)  V  1.0136 L E  1.0136  1.01351  104 3. Determinar volumen molar del oxigeno mediante la ecuación del VAN DER WAALS a    P  2  V  b   RT V  

P = 100 Atm.,

T = 700 K para un gas que tiene a = 1,36

 PV

2

b = 0,0318

 a  V  b   V 2 RT

PV 3  aV  PV 2b  ab  V 2 RT PV 3  PV 2b  V 2 RT  ab  aV  0 f V   100 V 3  3,18 V 2  57, 4 V 2  1.36 V  0, 043248 f V   100 V 3  60,58V 2  1,36 V  0, 043248 f ' V   300 V 2  121,16 V  1,36

26

Método de Newton Raspón de Primer Orden

v1  v0 

f v f ' v

v1  0,58 

 0,122248  32, 0072

v1  0,5838 Tabla Nº6 Aproximación de la raíz y error Vi

f(V)

f'(V)

Vi+1

E

0.5740

-0.3100

30.6560

0.5841

0.0101

0.5841

0.0115

32.9470

0.5837

3.5 x 10-4

0.5837

1.4 x 10-5

32.8670

0.5837

4.2 x 10-7

0.5837

2.0 x 10-11

32.8670

0.5837

6.0 x 10-13

0.5837

-1.0 x 10-13

32.8670

0.5837

3.0 x 10-15

El volumen seria: V=0.5837mol/L

4. Un gas se encuentra a una presión absoluta de 13.76 bar y una temperatura de 333 K. Encontrar el volumen molar ocupa el gas empleando la ecuación de estado de Redlich-Kwong. P

RT a  1/ 2 V  b T V (V  b )

Para este compuesto las constantes son: P = 13.76 atm T = 333ºk a = 1.5614 x 108 (cm6 bar/(g mol)2 k1/2) b = 44.897 (cm3/ g mol) R = 83.4 (cm3 bar /g mol k) Solución Despejando la ecuación Nº 1 Tenemos: 27

A  ab  PV 3  RTV 2  V  Pb 2  RTb  1/ 2   1/ 2  0 T  T  Reemplazando: F ( v )  13, 76 v 3  27685, 62 v 2  7300696, 52 v  384831290, 3

Aplicando el método de Bisección Tomamos: x1  370

;

f  x1   

v

f (v)

365

v

–

 f  v    11674093.22  

+

f (v)

370  367.5  368.75 2



v 370 367.5

v

–

 f  v    2.313174.844  

+

f (v)



v 370 368.75

v

– +

f (v)

368.75

– +

370  368.75  369.375 2

 f  v    2381900.225  



v 369.375

f  x2   

370  365  367.5 2



v 370

;

x2  365

369.375  368.75  369.06025 2

 f  v    33146.79984  

28

v

f (v)

369.0625  368.75  368.90625 2



v 369.0625

–

368.75

v

 f  v    1140318.152  

+

f (v)



v 369.0625 368.90625

v

–

369.0625  368.90625  368.984375 2

 f  v    553661.689  

+

f (v)



v 369.0625 368.984375

v

– +

f (v)

 f  v    260276.4453  



v 369.0625 368.0234375

v

– +

f (v)

369.0429688



v 369.0625

– +

f (v) –

369.0625  369.0234375  369.0429688 2

 f  v    113569.1971  

v 369.0625

369.0625  368.984375  369.0234375 2

369.0625  369.0429688  369.0527344 2

 f  v    40212.38605  



v

369.0625  369.0527344  369.0576172 2

 f  v    3533.08997   29

369.0527344

v

+

f (v)

369.0625  369.0576172  369.0600586 2



v 369.0625

–

369.0576172

v

 f  v    14806.78073  

+

f (v)

369.0600586  369.0576172  369.0600586 2



v 369.0600586 369.0576172

–

 f  v    5636.82683  

+

v  369, 0588379 cm 3 gmol E  v2  v1  2, 4414  10 3

5. El factor de fricción (f) para el flujo turbulento en una tubería está dado por la correlación de Colebrook:  D 1 2, 54   0,86 ln    f  3, 4 Re f 

Donde Re = es el número de Reynolds (adimensional) , es la aspereza o rugosidad de la tubería (unidad de longitud) D, es el diámetro de la tubería (unidad de longitud) Obtener el factor de fricción para un fluido con un Reynolds de 3E4 que fluye en una tubería con un diámetro de 0,1 m y una rugosidad de 0,0025m. Solución Despejamos la ecuación (1) 

1

Re   0,86 R f   e 2, 51

R f  

f  11952,19124 e

  1    0,86 f   

  f 



Re  1  2, 51  3, 7  D f

 80, 75814889 f  1 30

Utilizando el Método de Bisección Valores iniciales:

f

–

0.10

+

f

F ( f0 )  

f1  0,1



F ( f1 )  

f 

0.05  0.075  0.0625 2 F  f    7.347809254

F (f)

0.05

–

0.075

+

f



0.05  0.1  0.075 2 F  f    23.76508585

F (f)

0.05

f 0  0, 05

f 

0.0625  0.05  0.05625 2 F  f    0.8973609056

F (f)

0.05

–

0.0625

+

f

f 

0.05  0.05625  0.053125 2 F  f    1.866717794

F (f)

0.05

–

0.05625

+

f

f 

0.05625  0.053125  0.0546875 2 F  f    0.5236505211

F (f)

0.053125

–

0.05625

+

f

f 

F (f)

0.0546875

–

0.05625

+

f

F (f)

0.546875

–

0.05546875

+

0.05625  0.0546875  0.05546875 2 F  f    0.1771493539 f 

0.05546875  0.0546875  0.055078125 2 F  f    0.1756818303 f 

31

El factor de fricción es: f  0, 055078125en  e  0, 0546875; 0, 05546875 

Error: Error  E  x2  x1  7,8125  10 4

1.7 APLICACIONES EN INGENIERÍA QUÍMICA 1. En un Proyecto de Ingeniería Química se requiere que se determine exactamente el volumen molar y factor de comprensibilidad del amoniaco a una presión de 120 atm y 500 °K mediante la siguiente ecuación:

a    P  2  (V  b)  RT V   a

27  R 2 T 2    64  PC 

b

RTC 8 PC

Datos: TC = 405,5 °K

PC = 111,3 atm.

R = 0,82

2. La concentración c de una bacteria contaminante en un lago decrece según la expresión: c(t) = 80e-2t + 20e-0.5t Siendo t el tiempo en hs. Determinar qué tiempo se necesita para que el nº de bacterias sea  7. Utilizar un valor inicial t0 =0. (Newton Raphson) 3. Para el diseño hidrodinámico de un proceso aparece la ecuación 0, 242 CF

  log( R.CF )

Donde R es numero de Reynolds evaluar el valor de CF para ( R= 106 y R=103) 32

4. Calculo de la presión de vapor Una de las propiedades de una sustancia pura que más comúnmente se utiliza en cálculos de Termodinámicas es la presión de vapor o presión de saturación. Esta se define como la presión a la cual existen en equilibrio una fase liquida y una fase vapor. Si la presión de vapor iguala a la presión atmosférica, el líquido entrara en ebullición. Solo depende de la temperatura. Existen diversas ecuaciones para calcular. Una de las precisas es la ecuación de Frost-KalkwarfLn Pvap = A – B/T + ClnT + DPvap/T2 Donde: A, B, C, D: Constantes empíricas que dependen de cada sustancia. Calculemos la Pvap del etilbenceno a una temperatura de 347,25 ºK. Los valores de las constante son: A=58,1, B=6792,54, C=-5,802, D=5,75. Las unidades de T son ºK y las de Pvap mmHg.

1 y  y2  y3 5. En Química aparece una ecuación de la forma z  . (1  y) 3 ¿Cuánto vale y si z=0.6789?.

6. Para determinar la constante de nacimientos de una población dada, se necesita calcular el valor de la constante , (0.1,0.9) de la siguiente ecuación: 1 .564  10 6  10 6 e  

0 .435  10 6



( e   1) , para esto utilice el

método de Newton-Raphson hasta que el error sea del orden del 0.01% (0.10249)

7. La constante de disociación del agua pura es función de la temperatura como se muestra a continuación. Log(Kw) =

4470,99 + 6,0875 - 0,01706T T

T en K y Kw en M2

33

Calcular La temperatura en °C cuando el pH del agua condensada es de de 6,2 8. La temperatura media logarítmica

de un intercambiador de calor a

contracorriente esta dado por: LMTD 

T 3  T 2   T 4  T 1 T3 T2  ln    T 4  T1 

Figura Nº 3: Interrcambiadores de calor Se sabe que este sistema la temperatura media logarítmica debe ser de 52ªC( LMTD) . El fluido caliente se alimenta al sistema a 100ªC y sale a 40ªC, mientras que el fluido frió se alimenta a 8ªC ¿A qué temperatura sale del intercambiador el fluido frió? 9. A un reactor ingresa una mezcla de gases de 8% de dióxido de azufre, y 12 % de oxigeno y la diferencia de nitrógeno, y se desarrolla la siguiente reacción. SO2

+

1 O2  2

SO3

Calcular la composición en equilibrio a presión constante de 2 atm. Y la constante de equilibrio (Kp) es de 160. 10. La ecuación de estado R-K

P

RT a  V  b V (V  b ) T

 R 2Tc5 / 2   a  0.42747 P c  

34

 RT  b  0.8664 c   Pc  P = Presión en atm V = volumen molar en L/g-mol T = temperature en K R = gas constant (R = 0.08206 (atm·liter/g-mol·K)) Tc = critical temperature in K Pc = critical pressure in atm The compressibility factor is given by:

z

PV RT

(4)

Calcule el factor de comprensibilidad y el volumen molar a 50 atm t 500 ºC. 11. En un reactor químico ingresaron cierta cantidad de moles con cierta cantidad de impureza como se indica en la tabla de registro a una Temperatura de 400ºC:

Tabla Nº7 Composicion de alimentación al reactor SUSTANCIA Nitrógeno Hidrogeno Amoniaco Dióxido de carbono

MOLES QUE INGRESARON 1663,12 Kmol/h 4990,03 Kmol/h 986,45 Kmol/h 0,076 Kmol/h

Metano

334,16 Kmol/h

a partir de la ecuación de Larson que proporciona el valor de la constante de equilibrio de la reacción.

Log ( Kp )  2,113 

2074,8  2, 49 Log (T )  1, 256  10 4 T  1,85  10 7 T 2 T

Donde T en ªC y Kp en bares Calcular la composición a 300 ºC y 160 bares

35

2. SISTEMAS DE ECUACIONES NO LINEALES Un sistema de ecuaciones no lineales que se expresa como un sistema no lineal igualando a cero de la forma F(x) =0Este sistema puede presentar múltiples soluciones matemáticamente posibles y su solución numérica debe proporcionar la solución físicamente correcta, uno de los métodos es el de newton Raphson estándar 2.1 MÉTODO DE NEWTON RAPHSON 

Este procedimiento corresponde, analíticamente, a extender el uso de la derivada, ahora para calcular la intersección entre dos funciones no lineales.



Al igual que para una sola ecuación, el cálculo se basa en la expansión de la serie de Taylor de primer orden, ahora de múltiples variables, para considerar la contribución de más de una variable independiente en la determinación de la raíz.



Para dos variables, la serie de Taylor de primer orden se escribe, para cada ecuación no lineal:



ui 

ui u u u xi1  i xi  i yi1  i yi  0 x x y y

vi 

vi v v v xi1  i xi  i yi1  i yi  0 x x y y

Pero ui+1 = vi+1 = 0 : u u u u ui  i xi1  i xi  i yi1  i yi  0 x x y y

vi  

vi v v v xi1  i xi  i yi1  i yi  0 x x y y

Que reescribiendo en el orden conveniente:

Y cuya solución es:

ui u u u xi1  i yi1  ui  i xi  i yi x y x y v i v v v xi1  i yi1  v i  i xi  i yi x y x y

36

ui xi1  xi 

vi u u v  vi i v i i  ui i y y x yi1  yi  x J J

Donde J es el determinante jacobiano del sistema

 ui  x J  ui  y 

vi  x  vi  y 

2.1.1 APLICACIÓN 1. Se tiene un sistema de 3 reactores continuos tipo tanque en donde se lleva a cabo la reacción A  P operando isotérmicamente manteniendo los volúmenes constantes. Calcular la concentración de A en régimen permanente en cada reactor si la reacción es de segundo orden

FA0

R1

R2

R2

Figura Nº4: Sistema de reactores en serie

FA0 = 10L/min., FR = 5L/min. CA0 = 1mol/l. V1 = 100L, V2 = 50L, V3 = 50L, K = 0,1 RESOLUCIÓN: A partir de la siguiente ecuación tenemos: ACUMULACIÓN = ENTRADA-SALIDA-REACCIONA 0= ENTRADA-SALIDA-REACCIONA Las concentraciones en cada reactor son las siguientes: Primer reactor: CA1 Segundo reactor: CA2 Tercer reactor: CA3 37

FR

ECUACIÓN del primer reactor: F1 CA1 , CA2 , CA3   10  5CA3  15CA1  10C 2 A1 df 1  15  20CA1 dx df 1 0 dy df 1 5 dz

ECUACIÓN del segundo reactor: F2  CA1 , CA2 , CA3   15CA1  10C 2 A1

df 2  15 dx df 2  15  10CA2 dy df 2 0 dz ECUACIÓN del tercer reactor: F3 CA1 , CA2 CA3   15CA2  15CA3  5C 2 A3 df 3 0 dx df 3  15 dy df 3  15  10CA3 dz

APLICANDO EL METODO DE NEWTON RAPHSON ESTANDAR TENDREMOS:

df1 , dCA1

df1 , dCA2

df 2 , dCA1

df 2 , dCA2

df3 , dCA1

df3 , dCA2

df1 dCA3

h f1 ( CA1 , CA2 , CA3 ) df 2  i   f 2 ( CA1 , CA2 , CA3 ) dCA3 j f3 ( CA1 , CA2 , CA3 ) df3 dCA3 38

Reemplazaremos las derivadas y las funciones en la ecuación matricial anterior:  15  20CA1

0

h

0

i

  15 CA1  15 CA2  5 C 2 A2

 15  10CA3

j

15 CA2  15 CA3  5 C 2 A3

 15  10CA2

15 0

15

10  5CA3  15 CA1  10 CA 2 1

5

A partir de estas ecuaciones tenemos: CA1K 1  CA1K  h CA2K 1  CA2K  i CA3K 1  CA3K  j

Para k=0 CA11  CA10  h CA21  CA20  i CA31  CA31  j

PRIMERA ITERACIÓN

CA

0 1



, CA 20 , CA 30   0 . 7 , 0 . 3 , 0 . 1 

Reemplazamos en la matriz: 29

0

5

h

15

18

0

i

0

15

16

j

4,9  

5,55 2,95

h  0.1 i  0.2242 j  0.3945 Reemplazamos estos valores en las ecuaciones anteriores

39

CA11  CA10  h  0.7  h  0.7  0.1  0.6 CA21  CA20  i  0.3  i  0.3  0.2242  0.5242 CA31  CA31  j  0.1  0.3945  0.4945

Luego de cinco iteraciones tendremos Tabla Nº8 Aproximación de la raíz en cada reactor CA1

CA2

CA3

0.7

0.3

0.1

0.599

0.5242

0.4945

0.5848

0.5012

0.4382

0.5848

0.5009

0.4372

0.5845

0.5009

0.4372

Finalmente las concentraciones en cada reactor serian: Primer reactor: CA1=0.5845mol/L Segundo reactor: CA2=0.5009mol/L Tercer reactor: CA3=0.4372mol/L 2.2 PROBLEMAS 1. La reacción irreversible en fase liquida se lleva a cabo en tres reactores en serie como se muestra en la figura. El volumen de cada reactor es de 200L y la constante de la velocidad de reacción es de 0,5 l/mol-g min.

A B  R Calcular la concentración a l salida de cada reactor asumiendo un comportamiento estable.

40

Q = 6L/min. CA= 2 mol-g/l R2

R1

R3

Q = 6L/min. CA= 2 mol-g/l

2. A alta temperatura y baja presión

H2S y SO2 experimentan las

siguientes reacciones H2S( g)



H2(g) +

H2S( g) + SO2 (g)



1 S2 (g) 2

Kp1 =0,45atm1/2

2H2O(v) +

3 S2(g) Kp2 =28,5 atm1/2 2

La mezcla inicial contiene 45% de H2S( g) 25% de SO2 (g)

y gas

inerte nitrógeno con una presión total de 1,2 atm. Calcular las fracciones molares en equilibrio de todos los componentes 3. Una reacción

tiene lugar en una serie de cuatro reactores

continuamente agitados Como se muestra en la Figura. La reacción es irreversible y de primer orden

A B Las constantes de velocidad ki y los volúmenes Vi tienen los siguientes valores en cada reactor Tabla Nº9 Volumen de cada reactor y su constante de velocidad Reactor

Vi(L)

ki

1

1000

0.1

2

1500

0.2

3

100

0.4

4

500

0.3

41

Se supone: c0= 2mol/l El sistema está en estado estacionario.

Figura Nº5: Sistema de reactores en serie

4. Las siguientes reacciones se lleva a cabo en un reactor a 1000°C y presión de 1 atm.

C3 H8  3H 2O  3CO  7 H 2 C3 H8  6H 2O  3CO2  10H 2 A 1000 °C KP1 = 0,13x1012 y KP2 = 0,33X1012 Calcular los moles en equilibrio de cada sustancia para una alimentación de 1 mol de gas propano y 10 moles de vapor de agua.

5. En una reacción química A  B desde un balance de masa y energía se obtuvieron en las siguientes funciones 120X -75K(1-X) =0 -X(873-T) + 11(T-300) =0 K = 0,12 exp( 12581 ( T-298/298T) Donde x es la conversión, y T temperatura del sistema, y K constante de velocidad en base a esta información determine la conversión, temperatura y la constate de velocidad. 6. Para un sistema binario Liquido – Liquido el modelo de Wilson permite evaluar el coeficiente de actividad mediante la siguiente ecuación. A12 A21 Ln(  1 )  Ln ( x1  A12 x2 )  x2 (  ) x1  A12 x2 x2  A21 x1

42

Ln(  2 )  Ln ( x2  A21 x1 )  x1 ( Donde:

A12 A21  ) x1  A12 x2 x2  A21 x1

x 1 , x 2 fracciones molares

1 ,  2 coeficientes de actividad A12 , A21 Parámetros para cada sistema Evalué los parámetros para el sistema metanol agua cuando x 1 = 0,2 x 2 = 0,8 1 = 1,2 ;  2 = 1,6 7. La siguiente ecuación muestra el efecto de las variables x e y sobre el costo total para una operación particular: 11900 C = 2,33x + xy + 1,86y + 10 T

Determine los valores de x e y que den el costo total mínimo

8. En el sistema mostrado se lleva a cabo una reacción irreversible isotérmica de orden 1,8 respecto al reactante con los datos que a continuación se indica, calcular la concentración del reactante en los reactores 1 y 2. AB FA0 F

R1

R2

FR

Figura Nº6: Reactores en serie.

F=25 L/min, CA0 = 1mol/L FR =100 L/min K = 0,2 V1 =80L V2 = 20L.

43

3. INTERPOLACIÓN Introducción

Se denomina interpolación a la construcción de nuevos puntos partiendo del conocimiento de un conjunto discreto de puntos. En ingeniería y otras ciencias es frecuente disponer de un cierto número de puntos obtenidos por muestreo o experimento y pretender construir una función que los ajuste. Otro problema estrechamente ligado con el de la interpolación es la aproximación de una función complicada por una más simple. Si tenemos una función cuyo cálculo resulta costoso, podemos partir de un cierto número de sus valores e interpolar dichos datos construyendo una función más simple. En general, por supuesto, no obtendremos los mismos valores evaluando la función obtenida que si evaluásemos la función original, si bien dependiendo de las características del problema y del método de interpolación usado la ganancia en eficiencia puede compensar el error cometido. En la grafica se observa diferentes puntos ,

Figura Nº7: Puntos de variables x,y.

44

3.1. POLINOMIOS DE INTERPOLACIÓN DE NEWTON Uno de estas formas de interpolación se denomina Polinomios de Interpolación de Newton, que trabaja directamente en la tabla obtenida mediante el proceso de Diferencias Divididas; 3.1.1 Interpolación Lineal La forma más simple de interpolar es la de conectar dos puntos con una línea recta. Este método, llamado interpolación lineal, se muestra en la figura

Figura Nº8: Interpolación lineal

Usando triángulos semejantes, se tiene:

f 1( x)  f ( x0 ) f ( x1 )  f ( x0 )  x  x0 x1  x0 que se puede reordenar como :

f 1( x)  f ( x0 ) 

f ( x1 )  f ( x0 ) ( x  x0 ) x  x0

La cual es una fórmula de interpolación lineal. La notación f 1(X) indica que se trata de un polinomio de interpolación de primer orden. Nótese que además de representar la pendiente de la linera que conecta los dos puntos, el termino 45

f ( x1 )  f ( x0 ) x1  x0 Es una aproximación de diferencias divididas finitas a la primera derivada. En general, entre más pequeño sea el intervalo entre dos puntos, más exacta será la aproximación. 3.1.2 Interpolación Cuadrática Una estrategia que mejora la aproximación es la introducir cierta curvatura en la línea que conecta a los puntos. Si se dispone de tres datos, lo anterior se puede llevar a cabo con un polinomio de segundo orden (llamado también polinomio cuadrático o parábola). Una manera conveniente para este caso es: f 2  x   b0  b1 ( x  x0 )  b1 ( x  x0 )  b2 ( x  x0 )( x  x1 )

Nótese que aunque la ecuación parezca diferente de la ecuación general de un polinomio: f  x   a0  a1 x  a2 x 2  ...  an x n

Las dos ecuaciones son equivalentes. Se puede usar un procedimiento simple para determinar los valores de los coeficientes. Para bo , se usa la ecuación con X=X0 y se obtiene.

b0  f ( x0 ) Sustituyendo la ecuación y evaluando en X=X1 se obtiene

b1 

f ( x1 )  f ( x0 ) x1  x0

Y por ultimo las ecuaciones se sustituyen en la ecuación y se evalúa está en X=X2 y se obtiene

f ( x2 )  f ( x1 ) f ( x1 )  f ( x0 )  x2  x1 x1  x0 b2  x2  x0 Nótese que, al igual que en el caso de interpolación lineal, b1 aun representa la pendiente de la línea que une los puntos X0 y X1. Por lo tanto, los primeros dos términos de la ecuación son equivalentes a la 46

interpolación de X0 a X1. El ultimo termino b2(X-X0) (X-X1), introduce la curvatura de segundo orden de la formula. Tabla de diferencias divididas ---------------------------------------------------------------------------------------------

x

f  x

diferencia dividida

diferencia dividida

primer orden

segundo orden

---------------------------------------------------------------------------------------------x0

f  x0  f [ x0 , x1 ] 

x1

x1  x0

f  x1 

f [ x0 , x1 , x2 ]  f [ x1 , x2 ] 

x2

f  x1   f  x0  f  x1 , x2   f  x0 , x1 

f  x2   f  x1 

f  x2 

x2  x0

x2  x1

------------------------------------------------------------------------------------------------

3.1.3 Diferencias Finitas Divididas Sea F una función de valor real definida sobre xk , xk 1 , ..... xk  n no necesariamente equidistante. Se define: F  xk , xk 1  

F ( xk )  F ( xk 1 ) xk  xk 1

F  xk , xk 1 , xk  2  

 Primer grado

F ( xk , xk 1 )  F ( xk  2 ) xk  xk  2

F  xk , xk 1 ,......., xk  n  

 2 grado

F ( xk , xk 1 ,.......xk  n 1 )  F ( xk  n )  ngrado xn  xk  n

Ejemplo

f  x 0 ,x1 , x2  

f  x1 , x 2   f  x 0 , x1  x2  x0

47

En la práctica, los cálculos se disponen en una tabla de diferencias divididas, colocando en la primera columna los valores de la función o diferencias divididas de orden 0, en la segunda columna las diferencias divididas de primer orden, en la tercera columna las de orden 2, y así sucesivamente. Usando cuatro puntos ( x0 , y0 ), ( x1 , y1 ), ( x2 , y2 ), y ( x3 , y3 ),

f3 ( x)  f [ x0 ]  f [ x1 , x0 ]( x  x0 )  f [ x2 , x1 , x0 ]( x  x0 )( x  x1 )  f [ x3 , x2 , x1 , x0 ]( x  x0 )( x  x1 )( x  x2 ) b0

3.2.

x0

f x0 

x1

f x1 

x2

f x2 

x3

f x3 

b1 b2

f x1,x0 

b3

f x2 , x1, x0  f x3, x2 , x1, x0 

f x2 , x1 

f x3 ,x2 

f x3 , x2 , x1 

POLINOMIOS DE INTERPOLACIÓN DE LAGRANGE

Presentamos ahora una forma alternativa del polinomio de interpolación P(x) asociado con una tabla de datos (xi , yi) con 0≤ i ≤ n. Es importante entender que existe uno y solo un polinomio de interpolación de grado ≤ n asociado con los datos (suponiendo, claro está, que las n+1 abscisas xi son distintas). Sin embargo, existe ciertamente la posibilidad de expresar este polinomio de maneras distintas y de llegar a él a través de distintos algoritmos. El problema al utilizar Polinomio de Newton para aproximar es que se debe tener la derivada y muchas veces este lado no se tiene; una forma de evitar esto es trabajar con una interpolación de Lagrange, que es una reformulación del Polinomio de Newton que evita las diferencias divididas y se representa como: 48

x  xj

n

Li ( x )  

xi  x j

j0 ii

; Fn ( x ) 

n

 L ( x)F ( x )  R i0

i

i

n

Para obtener el polinomio de grado uno (lineal) reemplazamos n=1 1

F1 ( x)   Li ( x) F ( xi ) i 0

 L0 ( x) F ( x0 )  L1 ( x) F ( x1 )

j 0 i j

x  xj xi  x j



x  x1 x0  x1

x  x0 x1  x 0

L1 ( x ) 

F1 ( x) 

n 1

L0 ( x)  

x  x0 x  x1 f ( x0 )  f ( x1 ) x0  x1 x1  x0

Ahora calcularemos el polinomio de interpolación de grado dos (Cuadrático), haciendo n=2: 2

F2 ( x)   Li ( x) F ( xi ) i 0

 L0 ( x) F ( x0 )  L1 ( x) F ( x1 )  L2 ( x) F ( x2 ) F2 ( x ) 

x  x0 x  x 2 x  x0 x  x1 x  x1 x  x2 f ( x0 )  f ( x1 )  f ( x2 ) x0  x1 x0  x2 x1  x0 x1  x2 x2  x0 x2  x1

La aproximación del polinomio cúbico es: 3

F3 ( x)   Li ( x) F ( x i ) i 0

 L0 ( x) F ( x 0 )  L1 ( x) F ( x1 )  L2 ( x) F ( x 2 )  L3 ( x) F ( x 3 ) F3 ( x) 

x  x0 x  x2 x  x3 x  x1 x  x2 x  x3 f ( x0 )  f ( x1 ) x0  x1 x0  x2 x0  x3 x1  x0 x1  x2 x1  x3 

x  x0 x  x1 x  x3 x  x0 x  x1 x  x2 f ( x2 )  f ( x3 ) x2  x0 x2  x1 x2  x3 x3  x0 x3  x1 x3  x2

49

Ejemplo Con un polinomio de interpolación de Lagrange de primero, segundo y tercer grado evalué 1,5 ; basándose en los datos dados a continuación: x

1

2

3

4

F(x)

1

1,4142

1,732

2

Solución: Primero hallamos el polinomio lineal:

F1 ( x) 

x2 x 1 (1)  (1,414213)  0.414213x  0.585787  F1 (1.5)  1.2071065 1 2 2 1

Ahora hallamos el polinomio cuadrático:

x 2 x 3 x 1 x  3 x 1 x  2 (1)  (1.414213)  (1.732050)  F2 (1.5)  1.219153825 1 2 1 3 2 1 2  3 3 1 3  2 Finalmente el polinomio cúbico es: F2 ( x) 

x2 x3 x4 x 1 x  3 x  4 x 1 x  2 x  4 (1)  (1.414213)  (1.732050) 1 2 1 3 1 4 2 1 2  3 2  4 3 1 3  2 3  4 x 1 x  2 x  3  (2) 4 1 4  2 4  3 F3 ( x ) 

F3 ( x )  1,22205934

3.3 PROBLEMAS 1.

A 400ºC y con concentraciones:

CO  

NO 2   0 ,10

mol/L

se

obtuvieron los siguientes datos para la reacción dada: CO g 





NO 2  g 



CO 2  g 



NO  g 

Tiempo (s)

0

10

20

30

CO  mol L

0,1

0,067

0,05

0,04

Cuál es la concentración del monóxido de carbono para el tiempo de 25 s

50



Usando un polinomio de segundo grado mediante el método de LaGrange ¿cuáles son las constantes ( a0 , a1 ; a2 ; )

f(x) = a 0

+ a1x + a2x2 . 2.

En una prueba experimental de una sedimentación discontinua se a determinado los siguientes datos: t ( s)

Zi ( m )

0

0,35

300

0,21

600

0,129

900

0,094

1200

0,083

Determine el polinomio interpolante de grado 3 3. Con una función f las diferencias divididas progresivas están dadas por

x0=0,

f[x0]

0 f[x0,x1] x1=0, 04

f[x0,x1,x2]=

f[x1]

50 7

f[x1,x2]=10 x2=0,

f[x2]=

7

6

Determine los datos que faltan en la tabla. 4. Dados los nodos (0, -5), (1, -3), (2, 1), (3, 13), averiguar el polinomio interpolador mediante: La resolución de un sistema de ecuaciones. El método de Lagrange. El método de las diferencias divididas de Newton. 3 2 (Solución: P( x)  x  2 x  3 x  5 )

51

5. En un proceso de desalinización usando un equipo de osmosis inversa se ha registrado los siguientes datos p(bar)

5

10

15

20

25

Flujo( mL/h

103

202

287

386

501

¿Cuál es el polinomio usando todos los puntos mediante el método de Newton? 6. Cuál es el flujo para una caída de presión de 23 bar usando un polinomio de grado 3 r (cm)

0

5

10

15

20

25

30

35

40

45 47 50

V(m/s) 50 49,5 49

48

46,5

45

43

40,5

37,5 34

25

0

Determine el polinomio de grado 3 para r= 30 y r= 45cm 7. La tabla siguiente enumera la población de los Estados Unidos a partir de 1940 a 1990. años

1940

1950

Poblacion

132,165

151,326

1960

1970

1980

1990

179,323

203,302

226,542

249,633

(ien millares)

Encuentre el polinomio de Lagrange del grado 4 , y utilice este polinomio para estimar la población en los años 1965 8. Aproximar f (0.05) , mediante el polinomio interpolación de Newton usando la tabla siguiente x

0

0.2

0.4

0.6

0.8

f (x )

2.00001

2.13456

2.56789

2.67845

2.78954

52

4. REGRESIÓN Las técnicas de regresión permiten hacer predicciones sobre los valores de cierta variable Y (dependiente), a partir de los de otra X (independiente), entre las que intuimos que existe una relación. 4.1. REGRESIÓN LINEAL Consiste en ajustar a una línea recta un conjunto de parejas de datos observados: (X1, Y1), (X2, Y2),………, (Xn, Yn) y para ello utilizaremos el método de mínimos cuadrados. La expresión matemática de la función lineal es g(x)=a 0 + a1x donde a0 y a1 son constantes a determinar. La desviación de la recta con respecto a cada dato se define como: di = yi – g(xi)

y como g(xi) = a0 + a1 xi

Entonces di = yi – a0 – a1xi

donde i = 1,2………,n siendo n número de

datos. El cuadrado total de las desviaciones para obtener un buen ajuste está dado por: n

n

D    d i    yi  a0  a1 xi  i 1

2

2

i 1

cuando este valor es mínimo. Para determinar los valores de las constantes a0 y a1 y hacer mínimo a D, se deriva la ecuación anterior con relación a cada uno de los coeficientes (derivados parciales) y se igualan a cero n

Considerando que

a i 1

0

na0 , entonces las anteriores ecuaciones se pueden

expresar Como un conjunto de 2 ecuaciones lineales con 2 incógnitas:

na0  a1xi yi

a0 xi a1xi2 xi yi Al expresar estas 2 ecuaciones en forma matricial, tenemos:

53

 n   xi

 x  a     y   x   a   x y  i 2

i

i

0

1

i

i

Que se resuelven por los métodos conocidos ( gauss-seidel)

4.2. REGRESIÓN POLINOMIAL El método de mínimos cuadrados puede extenderse para ajustar un polinomio de cualquier valor a los datos de una medición. g(x) = a0 + a1x + a2x2 + amxm

di=yi – g(xi)

Desviación de la curva con respecto a cada Polinomio. La suma de los cuadrados de la desviación así: D =  (yi – a0 – a1 xi – a2xi2 + ……… amxim)2

Siguiendo el mismo procedimiento hecho para la regresión lineal, tenemos: derivados parciales con respecto a los coeficientes al polinomio.

D  2  (y i  a 0  a 1 x i  a 2 x 2i ..........a m x mi )  a0 D  2 x i ( y i  a 0  a 1 x i  a 2 x 2 .......a m x m2 )  a1 D  2 x 2i ( y i  a 0  a 1 x i  a 2 x 2i .......a m x mi )  a2

D  2 x mi ( y i  a 0  a 1 x i  a 2 x 2i .......a m x mi )  am Al igualar estas ecuaciones a cero y reordenando, se tiene el siguiente conjunto de ecuaciones. a0n + a1xi + a2xi2 + ……………+amxim = yi a0xi + aixi2 + a2x23 + …………+ amxim+1 = xiyi a0xi2 + a1xi3 + a2xi4 +…………+ amxim+2 = xi2 yi

54

a0xim + a1xim+1 + a2xim+2 +……+ amxi2m = ximyi (Todos los  van desde i =1 hasta n) el anterior sistema lo podemos expresar así:

x x x

 n   xi   xi2     m m 1  x   i  xi

i 2 i 3 i

  

  m2  xi 

x x x

m i m 1 i m2 i

  xi2 m

  a0    yi        a1    xi y    a2     x12 y             a   xm y    m   1 i 

El sistema de ecuaciones se resuelve por los métodos vistos para resolver sistema de ecuación. Ojo: si el polinomio es de grado m => necesitamos m + 1 ecuaciones lineales con m+1 incógnitas.

4.3. REGRESIÓN LINEAL MÚLTIPLE Una extensión útil en la regresión lineal es el caso en que y es una función lineal de dos o más variables. Por ejemplo, y pudiera ser una funciónlineal de x1 y x2, de la forma: y = a0 + a1x1 + a2x2 Tal ecuación es útil particularmente cuando se ajustan datos experimentales en donde la variable que está analizando, a menudo es función de otras dos variables. En este caso bidimensional, la “línea” de regresión viene a ser un “plano”. Como con los casos anteriores, los “mejores” valores de los coeficientes se determinan agrupando la suma de los cuadrados de los residuos: n

D    yi  a0  a1 x1  a2 x2 

2

i j

y derivando con respecto a cada uno de los coeficientes: D1  2  yi  a0  a1 x1  a2 x2  a0 D1  2x1  yi  a0  a1x1  a2 x2  (1) a0

55

D1  2 x2  yi  a0  a1 x1  a2 x2  a 2

Los coeficientes que generan la suma mínima de los cuadrados de los residuos se obtienen igualando cada una de las derivadas parciales a cero y expresando la ecuación (1) como un conjunto de ecuaciones lineales simultáneas, de la forma: nao + ∑ x1 a1 + ∑ x2 a2 = ∑ y1 ∑ x1ao + ∑ x12 a1 + ∑ x1 x2 a2 = ∑ x1 y1 ∑ x2ao + ∑ x1 x2 a1 + ∑ x22 a2 = ∑ x2 y1 o como una matriz:  n    x1   x2 

x x x x x x x x 1 2 1

1 2

  a0    y1      1 2   a1     x1 y1  2      x2 y1  2   a3    2

56

4.4 PROBLEMAS

1. Una ecuación para la variación de viscosidad del líquido con la temperatura

es:

b     a  e T 

Donde a, b son constantes empíricas y T, temperatura en K. T (ºC)

20

30

40

μ (N.s/m2)

1x10-3

8,08x10-3

50

6,53x10-3

5,33x10-3

Estime la viscosidad a 25 y 35ºC utilizando la ecuación empírica Solución: Sacando “ln” a ambos lados de la ecuación:  1  1  ln    ln a       b  T 

T (K)

μ (N.s/m2)

1/T

ln μ

293,15

1x10-3

3,411x10-3 -6,907

303,15

8,08x10-3

3,298x10-3 -7,12

313,15

6,53x10-3

3,193x10-3 -7,333

323,15

5,33x10-3

3,094x10-3 -7,536

Obtenemos una recta, de la forma: Cálculo de m: m 

y  c  mx

1 b

1 1 4    ln        ln   4  -0,093797  0, 0129978  -28,89962 T  T  m  2 2 4  4, 422915  105  1, 68943  104  1    1  4          T    T  m  1989,369



m

1 b

Cálculo de c:



b  5,0267 10-4

c  ln a 57

c  y  a1 x  -7,22491  1989, 369  3, 24945  10 3

c  -13,689

 c  ln a  -13,689



a  1,134106

1 ln    -13,689  1989 ,369   T

Por lo tanto el modelo es:  5,0267104      T    1,134  106 e

De la ecuación, calculamos μ para 25 y 35 ºC: T = 25ºC = 298,15 K

 298,5K

 5,026710  4      298 , 15   1,134  10  6 e

 8,966  10  4 N.s/m 2

T = 35ºC = 308,15 K

308,15K

 5,026710  4      308 , 15   1,134  10  6 e

 7,22  10  4 N.s/m 2

2. Se desea determinar la energía de activación del proceso de degradación de un contaminante a un producto inocuo en fase líquida. En el estudio de la 58

cinética de la reacción a diferentes temperaturas se obtuvieron los siguientes datos: 4,04x10-5

k (L/mol.s) T (ºC)

7,72x10-5

0

1,29x10-5

2,50x10-5

15

25

7

Determine la energía de activación Solución: De la ecuación de Arrhenius, determinamos el modelo:  Ea    RT  k  A  e

Sacando “ln” a ambos lados de la ecuación:

 E  1  ln k  ln A   a    RT  T  T (K)

k

1/T

ln k

273,15

4,04x10-5

3,66x10-3

-10,116

280,15

7,72x10-5

3,569x10-3

-9,469

288,15

1,29x10-5

3,47x10-3

-8,955

298,15

2,50x10-5

3,354x10-3

-8,294

Obtenemos una recta, de la forma: y  a0  a1 x Cálculo de a1: a1 

Ea R

1 1 4    ln k       ln k  4  -0,1297  0, 01405  -36,8355 T  T  a1   2 2 4  4,9437 105  1,9754  104  1    1  4          T    T  a 1  - 5850 ,1838



Ea  a1  R  - 5850,1838 K  8,314

Ea  48,638  103

J mol.K

J mol

Cálculo de a0: a0  ln A a0  y  a1 x  -9,208875   -5850,1838   3, 513734  10 3 59

a 0  11, 3471 

a 0  ln A  11, 3471



A  84720, 26

1 ln k   11,3471  5850 ,1838    T

Por lo tanto el modelo es:

k

 - 48,638103     RT     84720,26  e

3. Se cuenta con los siguientes datos de temperatura de ebullición de la acetona a diferentes presiones: P (atm)

1

T (ºC)

56,5

2 78,6

5 113,0

10

20

30

40

144,5

181,0

205,0

214,5

Si se grafican estos valores se observa que la temperatura no tiene un comportamiento lineal con respecto a la presión, por lo que se propone el modelo: T  AP B Determine las constantes A y B Solución: Sacamos el “ln” a ambos lados de la ecuación y obtenemos: 60

ln T   ln A   B  ln P 

P (atm)

T (ºC)

ln P

ln T

1

56,5

0

4,0342

2 5 10 20

78,6 113 144,5 181

0,6931 1,6094 2,3025 2,9957

4,3643 4,7273 4,9732 5,1984

30 205 3,4011 40 214,5 3,6888 Obtenemos una recta, de la forma: y  a0  a1 x

5,323 5,3683

Cálculo de a1: a1  B a1 

7   ln T  ln P     ln T   ln P  7   ln P      ln P   2

2



7  75, 565  14, 6909  33, 989 7  41,8263  215,8248

a1  0, 362

 B  0,362

Cálculo de a0: a0  ln A

a0  y  a1 x  4,8555  0,362  2,0987 a 0  - 4,095 

a 0  ln A  - 4, 095



A  60,039

ln T   4 , 095  0 , 362  ln  P 

61


T  60,039 P0,362 4. Experimentos de la cinética enzimática fueron llevados a cabo en el laboratorio y se registraron los siguientes datos:

  reactantes v rx

(s)(g/L)

(r)(g/L-h)

1

2

3

4

5

1

4

5

6

6

Determine las constantes a y b, si la expresión se ajusta al siguiente modelo: r

as bs

Solución: Invertimos la ecuación y obtenemos: 1 1 b 1      r a a s

  (s)

v rx

1

1/r

1/s

1

1

1

2

4

1/4

1/2

3

5

1/5

1/3

4

6

1/6

1/4

5

6

1/6

1/5

(r)

Obtenemos una recta, de la forma: Cálculo de m: m 

y  k  mx

b a

 1  1  1 1  19   137   107  5           5      r  s  r s 15   60   60     m  2 2 2  137   1    1  5  1, 463611111  5            60   s    s  m  1,074

Cálculo de k: k 

1 a 62

 107   137  k  y  mx     1,074    300   300  k  - 0,134  k

1  - 0,134 a



a  -7,459

Reemplazamos en “m”: m

b  1,074 a



b  1,074  -7,459  -8,015

1 1    - 0,134  1,074    r s


r

- 7,459 s - 8,015  s

PROBLEMAS 1.

La densidad del aire es función de la altitud como se muestra en los siguientes datos. Altitud( Km)

0,32

0,64

1,28

1,6

Densidad( Km/m3)

1,15

1,10

1,05

0,95

63

Se propone un modelo de regresión :   k1e k2 h Determine la densidad a 1,5 Km 2.

En una celda electrolítica se ha determinado la densidad de corriente enfunción de la caída de la tensión como se muestra a continuación.

Voltaje(V)

0,03

0,07

0,11

0,21

0,51

0,85

I(mA/cm2)

4

8

12

16

20

24

¿Cuál es la densidad de corriente a 0,25 V usando un polinomio de tercer grado? ¿Cuáles son las constantes ( a0 , a1 ; a2 ; a3 ) f(x) = a0 + a1x + a2x2 + a3x3 3.

Los siguientes conjuntos de mediciones de temperatura tomados de la cabeza de cilindros de un motor nuevo que está probando con miras a usarse en un auto de carrera. TABLA 10 Tiempo vs temperatura Tiempo (s)

Temperatura( ºF)

0

72,5

0,5

78,1

1

86,4

1,5

92,3

2

110,6

2,5

111,5

3

109,3

3,5

110,2

5

110,5

4,5

109,9

5

110,2

Calcule valores de temperatura en los siguientes intervalos usando un polinomio de grado 2 para 4,48

64

4.

Se desea determinar la energía de activación del proceso de degradación de un contaminante a un producto inocuo en fase líquida. En el estudio de la cinética de la reacción a diferentes temperaturas se obtuvieron los resultados que se muestran a continuación: k (L/mol.s) T (ºC)

4,04x10-5

7,72x10-5

1,29x10-4

2,5x10-4

0

7

15

25

Determinar: La energía de activación.

5. La ecuación de Arrenihus : K  Ae

 EA R .T

En la siguiente tabla se muestra la constante de velocidad en función de la temperatura. Temperatura ( K)

K(1/s)

273

7,78x10-7

298

3,46x10-5

318

4,98x10-4

328

1,5x10-3

338

4,87x10-3

Donde K : constante de velocidad de reacción EA: energía de activación R: constante universal de los gases R= 8,314J/ mol.K A: factor de frecuencia Determinar la constante de velocidad a 300 K E a =1.03021 105 joule mole

1

A= 4.00467 1013 s 1

65

6.

Un estudio en mecánica de fluidos indica que el flujo de un fluido a través de una tubería está relacionada con el diámetro de la tubería y la pendiente como se muestra en la siguiente tabla.

Experimento

Diámetro en pies

Pendiente,

Flujo( pie3/s)

1

1

0,001

1,4

2

2

0,001

8,3

3

3

0,001

24,2

4

1

0,01

4,7

5

2

0,01

28,9

6

3

0,01

84

7

1

0,05

11.1

8

2

0,05

69

9

3

0,05

200

Use una regresión lineal múltiple para analizar estos datos, la ecuación de potencia a evaluar es Q = a0Da1Sa2 Donde Q= Flujo volumétrico S: Pendiente D: diámetro de la tubería Después que se a encontrado el modelo, determine el flujo en una tubería de 2,5 pies y una pendiente de 0,025

7.

La densidad del aire es función de la altitud como se muestra en los siguientes datos. Altitud( Km)

0,32 0,64 1,28

1,6

Densidad( Km/m3) 1,15 1,10 1,05 0,95 Se propone un modelo de regresión :   k1e k2 h Determine la densidad a 1,5 Km 66

8.

Experimentos de la cinética enzimática fueron llevados a cabo en el laboratorio y se registraron los siguientes datos concentración del reactante (s) (g/L)

Velocidad de reacción (r) g/L-h

1

1

2

4

3

5

4

6

5

6

Determine la constante a y b Si la expresión de la velocidad se ajusta al siguiente modelo r

9.

a*s bs

Se cuenta con los siguientes datos de temperaturas de ebullición de la acetona a Diferentes presiones: P(atm) T(ºC) 1

56,5

2

78,6

5

113,0

10

144,5

20

181,0

30

205,0

40

214,5

Si se grafican estos valores se observa que la temperatura no tiene un comportamiento Lineal con respecto a la presión, por lo que se propone el modelo: T  AP B …..(1) Determinar los valores de las constantes A y B 67

5. DIFERENCIACIÓN NUMÉRICA Diferenciación numérica es una técnica de análisis numérico para producir una estimación del derivado de la función matemática o función usando valores de la función. Dos métodos son los más usuales a la hora de resolver tal problema: 

Derivar un polinomio de interpolación construido mediante alguno de los métodos estudiados



Calcular directamente la derivada utilizando para ello aproximaciones de la función mediante los polinomios de Taylor. Las formulas obtenidas de esta manera reciben el nombre de formulas.

La diferenciación numérica puede calcularse usando la definición de derivada

f '  x0   lim

f  x0  h   f  x0 

h 0

h Tomando una h pequeña. Si h > 0 se llama fórmula de diferencia progresiva, Fórmulas de diferencias divididas hacia adelante Sea los puntos x

Xi

Xi+1

Xi+2

XI+3

…….

f

f (Xi)

f (Xi+1)

f (Xi+2 )

f (XI+3 )

…….

5.1.DIFERENCIACIÓN MEDIANTE MÉTODO NEWTON Primera derivada

f '  xi  

f  xi 1   f  xi  h

, f '  xi  

 f  xi  2   4 f  xi 1   3 f  xi  2h

Segunda derivada

f ''  xi   f ''  xi  

f  xi 2   2 f  xi 1   f  xi  h2  f  xi 3   4 f  xi  2   5 f  xi 1   2 f  xi  h2

5.2. DIFERENCIACIÓN DE LAGRANGE: datos discretos Construimos el polinomio de Lagrange 68

L  x   L1  x  y1  L2  x  y2  L3  x  y3 

 x  x2  x  x3   x  x1  x  x3   x  x1  x  x2  y1  y2  y  x1  x2  x1  x3   x2  x1  x2  x3   x3  x2  x3  x1  3

Diferenciando la función f  x  L  x  

2 x  x2  x3 y  x1  x2  x1  x3  1

2 x  x1  x3 2 x  x1  x2 y2  y  x2  x1  x2  x3   x3  x2  x3  x1  3

Asumiendo espacio constante f  x 

2 x  x2  x3 2 x  x1  x3 2 x  x1  x2 y1  y2  y3 2 2 2 x x 2 x 2

f  x 


f   x1  

2 x1  x2  x3 2x  x  x 3 y1  4 y2  y3 2x  x  x y1  1 1 2 3 y2  1 1 2 2 y3  2 2 x 2 x x 2 x

f   x2  

2 x2  x2  x3 2 x  x1  x3 y  y1 2 x  x1  x2 y1  2 y2  2 y3  3 2 2 2 2 x 2 x x 2 x

f   x3  

2 x3  x2  x3 2x  x  x 2x  x  x y  4 y2  3 y3 y1  3 1 2 3 y2  3 1 2 2 y3  1 2 2 x 2 x x 2 x

f  x 


f   x  

y  2 y2  y3 1 2 1 y  y2  2 y3  1 2 1 2 x x x x 2

segunda derivada

5.3 PROBLEMAS 1. El flujo de calor en la interfaz suelo-aire puede calcularse con ley dT q  z  0   k  C dz z 0

69

Z (m) T (º C)

0 13,5

1,25 12

1,75 10

Donde q = flujo de calor, k = coeficiente de difusividad térmica (3.5x10 -7),  = la densidad del suelo (1800Kg/m3), C = calor específico del suelo (840).

f '  0   13.5

2  0   1.25  3.75

 0  1.25 0  3.75

 12

2  0   0  3.75

1.25  01.25  3.75

 10

2  0   0  1.25

 3.75  0 3.75  1.25

 1,333

q = -3.5x10-7 . 1800 (-1,333). =70.56 2. A 400ºC y con concentraciones: CO   NO2   0,10 mol/L se obtuvieron los siguientes datos para la reacción dada:

COg   NO2

g  

CO2

g   NOg 

Tiempo (s)

0

10

20

30

CO  mol L

0,1

0,067

0,05

0,04

Determinar el orden de la reacción dada Determinar la constante de velocidad Calcular la velocidad instantánea de la reacción en el instante t = 10 s Solución: De la ecuación: vA  -

dCA  k  CA  dt

n

 dC  Tenemos: ln  - A   ln x  n ln  CA  ecuación lineal  dt  Diferenciación numérica con 4 puntos: x0

x1

x2

x3

Tiempo (s)

0

10

20

30

CO  mol L

0,1

0,067

0,05

0,04

y0

y1

y2

y3

h  10 70

dCA 1   11y0  18 y1  9 y2  2 y3  dt 6h

Para x0:



dCA dt



1 11 0,1  18  0,067  9  0,05  2  0,04  4,4 103 6  10

dCA 1   2 y0  3y1  6 y2  y3  dt 6h

Para x1:



dCA dt



1 2  0,1  3  0,067  6  0,05  0,04  2,35 103 6 10

dCA 1   - y0  6 y1  3y2  2 y3  dt 6h

Para x2:



dCA dt



1 - 0,1  6  0,067  3  0,05  2  0,04  1,2 103 6  10

dCA 1   2y0  9 y1  18 y2  11y3  dt 6h

Para x3:



dCA dt



1 2  0,1  9  0,067  18  0,05  11 0,04  9,5 104 6  10

Llenamos la tabla y procedemos a graficar: TABLA 11 CALCULO DIFERENCIACIÓN

CO  mol L

- dCA

0

0,1

10

 

 dC A  ln dt  

ln CA

4,4 x 10-3

-5,426

-2,302

0,067

2,35 x 10-3

-6,053

-2,703

20

0,05

1,2 x 10-3

-6,725

-2,995

30

0,04

9,5 x 10-4

-6,959

-3,218

Tiempo (s)

dt

71

Del ajuste lineal tenemos que la ecuación es: y  1, 7 x  1, 41 Donde: Orden de reacción es: n  pendiente  1, 7  2 Constante de velocidad es: ln K  -1,41 Para t=10 s: v  k  C At 10 s 



K  0,244

n

v  0,244   0,067 2  1,095  103

PROBLEMAS 1. La primera ley de difusión de Fick establece que el flujo másico Flujo másico = F  D

dc dx

Donde flujo másico es cantidad de masa que pasa a través de una unidad de area por unidad de tiempo ( g/ cm 2/s. D= coeficiente de difusión ( cm2/s). C: concentración X: distancia en cm. Un ing. químico mide la concentración de un contaminante en los sedimentos depositados en un lago obteniéndose los siguientes datos. X (cm)

0

1

2

c. 10-6 g/cm3

0,1

0,4

0,9

72

Use la mejor técnica de derivación numérica para estimar en x =1,5 y calcule el flujo másico del contaminante cuando D= 2x10 -6cm2/s, para un lago de 3x106m2 de sedimentos ¿cuánto contaminante podría ser transportada hacia el lago en un año? 2. Para la descomposición del pentaóxido de dinitrógeno a 45ºC se obtuvieron los siguientes datos:

2N 2 O5 ( g )  O2 ( g )  4NO2 ( g ) Tiempo(s) 0

 N 2 O5 

200

400

600

800

mol 2,5 2,22 1,96 1,73 1,53 L

determina de qué orden será la reacción 3. A 440ºC con y con concentraciones |CO| = |NO2| = 0,10 mol/L se obtuvieron los siguientes datos para la reacción dada:

CO( g )  NO2 ( g )  CO2 ( g )  NO( g ) Tiempo(s)

CO 

mol L

0

10

20

30

40

100

1000

0,1 0,067 0,05 0,04 0,033 0,017 0,002

a) Determina el orden de la reacción dada b) Determine la constante de velocidad c) Calcula la velocidad instantánea de la reacción en el instante t= 10 s.

4. Se ha estudiado la descomposición térmica de arsenamina sobre vidrio, 2 AsH3(g) 2 As(s) + 3H2(g) La presión total del sistema varia con el tiempo, a 350C, según se indica a continuación: t/h

0

4.33

16.00 25.50 37.66 44.75

Pt/cmHg 39.2 40.30 43.65 45.35 48.05 48.80

a) Determine el orden de reacción con respecto a la arsenamina. b) Calcule la constante de velocidad. 73

5. Para poder estudiar el comportamiento físico químico de una solución de bromo en presencia de luz solar, se coloca una pequeña cantidad de solución de bromo en un matraz y se expone al sol, obteniéndose los siguientes resultados: Tiempo( min.)

10

20

30

40

50

60

Ppm Br2 a 25ºC 2,45 1,74 1,23 0,88 0,62 0,44

En base a esta información calcular la constante de velocidad en las unidades adecuadas y orden de la reacción.

6. En el estudio experimental de la velocidad de reacción para la reacción 2 NO (g) + 2 H2(g) → N2(g) + 2 H2O (g) se han obtenido los siguientes datos a 298 K: [NO] mol/L 0,02 0,02 0,03 0,06 Calcular:

[H2] mol/L 0,15 0,45 0,50 0,50

Velocidad mol/l s 1,2 x 10-8 3,6 x 10-8 9,0 x 10-8 36,0 x 10-8

a) La ecuación de velocidad para esta reacción. b) El orden de reacción c) La constante de velocidad

7. Determinar el orden de reacción para la descomposición en fase gaseosa del peróxido de di-tert-butilo (CH3)3COOC(CH3)  C2H6 + 2CH3COCH3 La reacción se lleva a cabo en el laboratorio en un reactor isotérmico batch La presión total de reacción a lo largo de reacción evoluciona según la siguiente tabla

74

tiempo(min) Presión total (MmHg) 0.0 7.5 2.5 10.5 5.0 12.5 10.0 15.8 15.0 17.9 20.0 19.4

8. Determinar el orden de reacción: CH3-Cl(g) + H2O(g)  CH3-OH(g) + HCl(g) usando los datos de la tabla. v  k  [CH3 -Cl]n  [H2 O]m

Experiencia [CH3-Cl] (mol/l) [H2O] (mol/l) v (mol·l–1·s–1) 1

0,25

0,25

2,83

2

0,50

0,25

5,67

3

0,25

0,5

11,35

9. La transferencia de calor está dada por la siguiente ecuación dT q   kA dy donde

J   k = conductividad térmica   smK 

 

A = área m2 T = temperatura K  y = distancia m  k  0.025

J smK

A  3 m2

La temperatura en función de distancia

T  1493y3  2200y 2  1076y  500 La transferencia de calor para y= 1,2m

75

6. INTEGRACIÓN NUMÉRICA En esta lección comenzamos el estudio de métodos numéricos para el cálculo numérico de integrales de la forma.  f  

b

 f  x  dx a

Un método común para aproximar I(f) es reemplazando f(x) con un polinomio de interpolación. Este procedimiento se conoce como las reglas de Cuadratura de Newton. Examinamos los primeros dos casos de este método donde se usan polinomios de interpolación lineales y cuadráticos.

6.1

Método del trapezoide: Sea p1(x) el polinomio lineal que interpola a f(x) en x=a y x=b

f(x) L(x)

x0

x1

x

Figura Nº9: Integración método de trapecio usando dos puntos

Usando la fórmula para el área de un trapezoide o integrando p1(x) directamente se obtiene que



b



b

a

a

f ( x)dx 

h  f ( x0 )  f ( x1 ) 2

f ( x ) dx 

h  f (a )  f (b) 2

Generalizando Generalizando el método de trapecio

76

f(x )

h 

b a n x0

h

x1

h

x2

h

x3

h

x4

x

Figura Nº10: Integración método de trapecio usando n puntos



b

a

x1

x2

xn

x0

x1

x n 1

f(x)dx   f(x)dx   f(x)dx    

f(x)dx

h h h  f(x 0 )  f(x1 )   f(x1 )  f(x 2 )     f(x n 1 )  f(x n ) 2 2 2 h   f(x 0 )  2f(x1 )    2f(x i )    2 f ( xn 1 )  f ( xn )  2 b h a f(x)dx  2 f(x 0 )  2f(x1 )    2f(x i )    2 f ( xn1 )  f ( xn ) 

4

Evaluar la integral :  xe 2 x dx método de trapecio 0

4

I   xe 2 x dx  0

40  f (0)  f (4)  2(0  4e8 )  23847.66 2

6.2 REGLA DE SIMPSON Además de aplicar la regla trapezoidal con segmentos cada vez más finos, otra manera de obtener una estimación más exacta de una integral, es la de usar polinomios de orden superior para conectar los puntos. Por ejemplo, si hay un punto medio extra entre f(a) y f(b), entonces los tres puntos se pueden conectar con un polinomio de tercer orden. A las fórmulas resultantes de calcular la integral bajo estos polinomios se les llaman Reglas de Simpson.

77

6.2.1 REGLA DE SIMPSON 1/3 La Regla de Simpson de 1/3 proporciona una aproximación más precisa, ya que consiste en conectar grupos sucesivos de tres puntos sobre la curva mediante parábolas de segundo grado, y sumar las áreas bajo las parábolas para obtener el área aproximada bajo la curva. Por ejemplo, el área contenida en dos fajas, bajo la curva f(X) en la Figura. 2, se aproxima mediante el área sombreada bajo una parábola que pasa por los tres puntos:

f(x)

x0

h

L(x)

x1

h

x2

x

Figura Nº11: Integración método de Simpson usando tres puntos



b a

f ( x)dx 

2

c i0



i

f ( x i )  c 0 f ( x 0 )  c1 f ( x1 )  c 2 f ( x 2 )

h  f ( x 0 )  4 f ( x1 )  f ( x 2 )  3

Generalizando

78

h

ba n f(x)

…... x0 h x1 h x2 h x3 h

x4

xn-2 xn-1

xn

x

Figura Nº12: Integración método de Simpson usando n puntos



b

a

x2

x4

xn

x0

x2

xn 2

f(x)dx   f(x)dx   f(x)dx    

f(x)dx

h h f(x 0 )  4f(x1 )  f(x 2 )  f(x 2 )  4f(x 3 )  f(x 4 ) 3 3 h     f(x n  2 )  4f(x n 1 )  f(x n )  3 h   f(x 0 )  4f(x1 )  2f(x 2 )  4f(x 3 )  2f(x 4 )   3 4f(x 2i-1 )  2 f ( x2i )  4f(x 2i 1 )   

2 f ( xn  2 )  4 f ( xn 1 )  f ( xn )  Evaluar la integral



4

0

xe 2 x dx 4

I   xe 2 x dx  0



h  f (0)  4 f (2)  f (4) 3

2  0  4(2e 4 )  4e8   8240.411 3

Simpson’s 3/8 Approximate by a cubic polynomial La derivación de la Regla de los Tres Octavos de Simpson es similar a la regla de un tercio, excepto que se determina el área bajo una parábola de tercer grado que conecta 4 puntos sobre una curva dada. La forma general de la parábola de tercer grado es: 79

L(x)

x0

h

f(x)

x1

h

x2

h

x3

x

Figura Nº13: Integración método de Simpson usando cuatro puntos



b

a

3

f ( x) dx   ci f ( xi )  c 0 f(x 0 )  c1f(x1 )  c 2 f(x 2 )  c3 f(x 3 ) i 0



h

3h  f ( x0 )  3 f ( x1 )  3 f ( x2 )  f ( x3 ) 8

b-a 3

6.2.2 SIMPHSON 3/8 4

I   xe 2 x dx  0

3h  4 8  f (0)  3 f ( )  3 f ( )  f (4)  8  3 3

3(4/3) 0  3(19.18922)  3(552.33933)  11923.832  6819.209 8 5216.926  6819.209   30.71% 5216.926 

80

6.3 PROBLEMAS 1.

Los datos que aparecen en la tabla siguiente corresponden a una velocidad transversal de salida en una tubería ¿cuál es el flujo volumétrico en m3/s?

r (cm)

0

5

10

15

V(m/s) 50 49,5 49 48

20

25

30

35

40

45 47 50

46,5

45

43 40,5 37,5 34 25 0

Q  2  V * rdr

2.

La cantidad de calor necesaria para elevar la temperatura de 1000gr de H2O desde –1000C hasta 2000C se puede evaluar de la siguiente forma:  H



 m

T 2 T 1

 T  dT

 2.64 107 T 2  1.56 104 T  0.132 Determinar  empleando la técnica del trapecio y la regla de Simpson

3.

La masa que entra

o sale de un reactor perfectamente agitado en

un periodo especifico se puede determinar mediante: t2

M   Q.cdt ti

Donde t1 y t2 son el tiempo inicial y final respectivamente Q es el flujo volumétrico y c es la concentración.

81

4.

T(,min )

C,mg/m3

0

10

5

22

10

35

15

47

20

55

25

58

30

52

35

40

40

37

45

32

50

34

Los datos que se muestra a continuación corresponde a la reacción de tipo. A  B en un reactor tubular X

(-1/rA)

0

0.0053

0.1

0.0052

0.2

0.005

0.3

0.0045

0.4

0.0033

0.5

0.0025

0.7

0.0018

0.8

0.00125

0.85

0.001

La ecuación de diseño para este tipo de reactor es de: V  FA0  

x

0

1 dx  rA

V:volumen(m3)FA0: flujo de alimentación A ,X: conversión el flujo molar de alimentación es de 2 mol/ min. el 80 % es de A y la diferencia de inertes Determine el volumen del reactor para una conversión de 90% de A 82

5.

Se tiene los siguientes puntos con h constante x

X0

X1

X2

X3

X4

y

Y0

Y1

Y2

Y3

Y4

Demostrar una fórmula de integración

6.

La reacción de saponificación de acetato de etilo en medio alcalino se conduce en fase homogénea de acuerdo con la siguiente ecuación. CH3COOCH2CH3

+ NaOH 

CH3COONa + CH3CH2OH

Donde la velocidad de reacción es de

ra 

5817 dCA  kCACB , k  2, 028 *109 * e T dt

de reacción

Constante de velocidad

el tiempo de reacción para un reactor discontinuo es

calculado mediante la siguiente reacción. CA

t

1 dCA ( ra ) CA0



Donde : CA = CA0(1-XA) ; CB = CB0 – CA0XA Si la concentración inicial de A y B es de 0,8 mol-g/ L ¿ Cual es el tiempo para alcanzar una conversión de 95%? A 500K

7.

La cantidad de calor necesaria para elevar la temperatura de 1 000 g de H2O

desde 100C hasta 900C se puede evaluar de la siguiente

forma:

H  m  T dT T2

T1

Donde:

 2.64 107 T 2  1.56 104 T  0.132 Determinar  empleando la técnica del trapecio y la regla de Simpson 8.

Resolver la siguiente integración.

83

9.

En un reactor tubular continuo se realiza la eliminación de un contaminante en fase gas mediante la reacción: A + B  R + S, a 200ºC y 7 atm. Para ello, se alimenta una corriente con A y B en relación equimolar, en las mismas condiciones de presión y temperatura, con un caudal de 500 L/min. Si se quiere alcanzar una conversión del 85%, ¿qué volumen de reactor es necesario? Datos: (rA) = 10.CA.CB [mol/L.min]

10.

Una reacción A  R, de ecuación cinética r (mol/L.h) = 1,5 CA (a 50ºC), se lleva a cabo en un reactor discontinuo, de forma isoterma, siendo la temperatura de trabajo constante e igual a 50ºC. La mezcla reaccionante tiene una concentración inicial de 10 mol/L en A ¿Cuál es el tiempo de reacción necesario para alcanzar una XA de 0,90?

84

7. SISTEMA DE ECUACIONES LINEALES

INTRODUCCIÓN

La solución de sistemas de ecuaciones lineales es un tema de gran utilidad en diversas ramas del conocimiento como la economía, la biología, física, ingenierías, etc. La resolución de sistemas de cualquier número de ecuaciones (10, 50, 100, 500, etc.) es una realidad hoy en día, gracias a los computadores, lo cual proporciona solución directa. Un gran número de problemas prácticos de ingeniería se reduce a resolver un sistema de ecuaciones.

Un sistema de

ecuaciones lineales con

incógnitas, tiene la forma:

a11x1 + a12x2 + ….. + a1nxn = b1 a21x1 + a22x2 + ….. + a2nxn = b2  am1x1 + am2x2+ …. + amnxn = bm Con notación matricial se escribe así: a11 a  21   a m1

a12 ..... a1n  a 22 ......a 2 n    a m 2 .....a mn 

 x1  b1  x     2   b2           xm  bm 

 AX=B Donde

A es la matriz de coeficientes del sistema X es el vector incógnita B es el vector de términos independientes.

7.1 METODOS DE JACOBI (Método de desplazamiento simultáneos)

a11 x1 + a12 x2 + a13 x3 = b1 a21 x1 + a22 x2 + a23 x3 = b2 a31 x1 + a32 x2 + a33 x3 = b3 85

cona11, a22 y a33 distintos de cero. Se despeja x1de la primera ecuación,

de la segunda y

de la tercera con lo

que se obtiene. =−

=− que en notación matricial queda

=

+

+

−

+

−

=−

⎡ 0 ⎢ = ⎢− ⎢ ⎢− ⎣

−

− −

+

⎡ ⎢ ⎢ +⎢ ⎢ ⎢ ⎣

−

⎤ ⎥ ⎥ − ⎥ 0 ⎥ ⎦

0

⎡ 0 ⎢ = ⎢− ⎢ ⎢− ⎣

− −

0

−

⎤ ⎥ ⎥ − ⎥ 0 ⎥ ⎦

⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎡ ⎢ ⎢ =⎢ ⎢ ⎢ ⎣

⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦

Una vez que se tiene la forma en notación matricial, se propone un vector inicial

( )

que puede ser

solución .

( )

= 0, o algún otro que sea aproximado al vector

Para iterar existe dos variantes. Si

( )

=

es el vector aproximación a la solución

después de

tiene para la siguiente aproximación.

86

iteraciones, entonces se

(

)

1 ( ⎡ ⎢ ⎢ 1 ( =⎢ ⎢ ⎢ 1 ( ⎣

=

−

−

−

−

−

)⎤ ⎥ ⎥ )⎥ ⎥ )⎥ ⎦

−

7.2 METODO DE GAUSS - SEIDEL (método de desplazamientos sucesivos) En este método los valores que se van calculando en la (

+ 1) −ésima

iteración se emplean para calcular los valores faltantes de esa misma iteración; es decir, con

( )

se calcula

(

)

de acuerdo con:

( ⎡ ⎢ 1 ( =⎢ ⎢ 1 ⎢ ( ⎣

−

−

−

−

−

87

−

)

⎤ )⎥ ⎥ ⎥ )⎥ ⎦

7.3 PROBLEMAS 1) Determine los valores de x1 , x2 y x3 usando el método de gaus seidel

2 x1  7 x2  11x3  6 x1  2 x2  x3   5 7 x1  5x2  2 x3  17 2) Se tiene la siguiente torre de separación realizar los balances y calcular los flujos masicos de cada corriente. H2O

0,2% de HCl 99.8% aire

Masa = 2640Kg/h

43,1% de HCl

63,8% H2O

56,9% aire

36,2% de HCl

Figura Nº14:Diagrama de una columna de absorción 3) Un proceso de 5 etapas en equilibrio para una extracción liquido

o

absorción gaseosa puede ser modelado mediante un sistema

de

ecuaciones lineales de acuerdo como se muestra a continuación.

(1  rr ) X 1  rrX 2   F0

X i1  (1  rr) X i  rrXi1  0 para i = 2,3, ………..(n-1) X n 1  (1  rr) X n  F Para : rr = 0,9 ; F0 = 0,05 y F = 0,5 Resolver el sistema de ecuaciones lineales

88

4) Para el sgt.sistema de separación, se conoce el flujo masico(Kg/hr) y las fracciones de masa de cada sustancia de entrada de la columna de separacion1 (F0) y la salida de la columna 2(F3 ,F1Y F4)¿ Cual es el valor de F2 

F1

W 1 = 0.2 W 2 = 0.6 W 3 =0.2

W 1 = 0.2 W 2 = 0.6 W 3 =0.2

F0

W 1 = 0.2 W 2 = 0.6 W 3 =0.2

F3

F2

F4

W 1 = 0.2 W 2 = 0.6 W 3 =0.2

5) Resolver la siguiente sistema de ecuación

4 x1  x 2   1  8 x1  x 2  x 3  1 3 3 x2  2 x3  4 x4   3 x 3  x 4  x 5  2 ,1 2 x 4  6 x 5  3, 4 6) La resistencia de un termistor varia con la temperatura 1 2 3  a 0  a1 ln( R )  a 2 ln R   a 3 ln R  T

T es la temperatura en Kelvin, R es resistencia en ohms, y a0 , a1 , a2 , a3 son constantes R

T

ohm

C

1101.0

25.113

911.3

30.131

636.0

40.120

451.1

50.128

Determine las constantes 89

7) En la figura se muestra la disposición de 5 tanques donde se lleva a cabo una mezcla de corrientes de entrada y salida. En base a la información determine concentración en cada en estado estacionario

2m3/min. C5 3

3m /min. 2m3/min. 5m3/min. 10mg/m3

C1 C1

C2

C4

3

3m /min.

3

1m /min.

1m3/min.

C3 1m3/min.

8m3/min.

8m3/min.20mg/m3 Figura Nº15: Sistema de reactores

90

11m3/min.

8. ECUACIONES DIFERENCIALES NUMÉRICAS

INTRODUCCIÓN Con mucha frecuencia aparecen problemas en ingeniería, física, química, ecología meteorología, sociología, etc., que exigen el manejo de ecuaciones diferenciales, muchas de las cuales no se pueden resolver pos los métodos convencionales, teniéndose que recurrir a métodos numéricos aproximados. Se llama ecuación diferencial aquella ecuación que contiene una variable dependiente y sus derivadas con respecto a una o más variables independientes. Se dividen en dos grandes grupos: Ordinarias, si contienen una sola variable independiente y Parciales, cuando contienen varias variables independientes. Estudiaremos únicamente las ecuaciones diferenciales ordinarias (EDO). Estas se clasifican y estudian según el orden de la mayor derivada que aparece en la respectiva ecuación diferencial. Ejemplo:

dy  ky dt

m

d2 y  ky dt 2

es de 1er orden es de 2º orden

d3 y dy  5  6y  0 es de 3er orden 3 dx dx Los problemas que encierran el uso de ecuaciones diferenciales ordinarias, constarán de: Una ecuación diferencial ordinaria: dy  y1  f(x, y) dx

Las condiciones iniciales. El intervalo para la variable independiente x, para el cual debe determinarse la función y, si existe.

91

8.1.

MÉTODO DE EULER

Consiste en dividir el intervalo de X0 a Xn en n subintervalos de ancho h.

h

xn - xo obteniéndose un conjunto discreto de (n + 1) punto: x0, x1, x2 …… n

xn en el intervalo [x0 , xn]

generándose la sucesión de aproximaciones

siguientes: y1 = y0 + hf (x0 , y0) y2 = y1 + hf (x1 , y1) y3 = y2 + hf (x2 , y2) . yn+1 = yn + hf(xn , yn) con i = 0 hasta n siendo f(xi , yi) la ecuación diferencial evaluada en xi y yi. yi+1 = yi + hf(xi, yi)

Ejemplo 1: Utilice el método de Euler para integrar numéricamente la ecuación: y’=f(x,y) = -2x3 + 12x2 – 20x + 8.5 de x=1.0 hasta x=3.0 considerando 4 intervalos y el valor inicial en x=1 es y=3 x0 = 1

h

f(x0 , y0) = f(1 , 3) 1.0

1.5

x0

x1

2.0 x2

2.5 x3

3.0 x4x4

3 1  0.5 4

yi+1 = yi + hf(xi , yi)

i

xi

yi

0 1

3

F(xi , yi) yi+1=yi+hf(xi ,yi) -1.5

2.25

1 1.5 2.25

-1.25

1.625

2 2.0 1.625

0.50

1.875

3 2.5 1.875

2.25

3.000

4 3.0 3.000

2.5

4.25

92

Ejemplo 2. Dada la ecuación diferencial

dy  f(x, y) = x - y , resuelva por el dx

método de Euler para el intervalo [0, 1] considerando 5 intervalos y siendo y(0)=2 1 0 = 0.2 5

h

Y(1)=?

yi+1 = yi + hf(xi , yi)

i 0

xi 0

x0 0

x1

x2

x3

0.2

0.4

0.6

yi

x4

x5

0.8 00.800.80

F(xi , yi) yi+1=yi+hf(xi ,yi)

2

-2

1.6

1 0.2

1.6

-1.4

1.32

2 0.4

1.32

-0.92

1.136

3 0.6

1.136

-0.536

1.0288

4 0.8

1.0288

-0.2288

0.98304

5 1.0

0.9834

-0.01696

0.986432

El valor exacto es 1.10364 E R 

1.10364  0.986432  100%  10.93% 1.10364

APLICACIÓN Un balón

de acero

de 1200 K es enfriado con corriente de aire

a

temperatura ambiente de 300K. Asumiendo la transferencia de calor en forma de radiación es mediante la ecuación diferencial . dT  2,2067  10 12 T 4  81 108 , T 0   1200 K dt Determine la temperatura para t= 480s usando el método de euler dT  2, 2067  1012 T 4  81  108  dt f  t , T    2, 2067  10 12  T 4  81  108 

T1  T0  f  t0 , T0  h Ti 1  Ti  f  ti , Ti  h  1200  f  0,1200 240





 1200   2.2067  10 12 1200 4  81  10 8  240 93

 1200   4.5579 240

T1  106.09K  0  240

t  t1  t0  h  240 T1  T  240  106,09K

2  1  f t1 ,1  h  106.09  f  240,106.09 240





 106.09  2.2067  10 12 106.09 4  81  108  240

 106.09   0.017595 240  110.32K

t  t2  t1  h  240  240  480

2    480  110.32K

8.2.

MÉTODO DE EULER MODIFICADO

Una fuente fundamental de error en el método de Euler es que la derivada al principio del intervalo se supone que se aplica a través del intervalo entero y para obtener una exactitud razonable se utiliza un intervalo muy pequeño a cambio de un error de redondeo mayor. El método de Euler modificado trata de evitar este problema utilizando un valor promedio de la derivada tomada en los dos extremos del intervalo en lugar de la derivada tomada en un solo extremo. El método de Euler modificado consta de dos pasos básicos:

1.

Se parte de (x0,y0) y se utiliza el método de Euler a fin de calcular el valor de y correspondiente a x1. Este valor de y se denotará como y1 ya que solo es un valor transitorio para y1. Esta parte del proceso se conoce como paso predictor.

94

2.

El segundo paso se llama corrector, pues trata de corregir la predicción. En el nuevo punto obtenido  x1 , y1  se evalúa la derivada de f  x1 , y1  usando la ecuación ordinaria que se está resolviendo. Se

obtiene la media aritmética de esta derivada y la derivada en el punto inicial (x0,y0). 1  f(x 0 , y 0 )  f(x1 , y1 )   derivada promedio 2

se usa la derivada promedio para calcular un nuevo valor de yi con la ecuación de Euler

yn+1 = yi + hf(xi, yi) que deberá ser más exacto que

y1.  y1  y 0 

h  f(x 0 , y 0 )  (f(x1 , y1 ))   2

que se tomará como valor definitivo de y1.

xn  x0   h    n  Este procedimiento se repite hasta llegar a yn. El esquema iterativo para este método quedaría en general así: Usando el paso de predicción resulta:

yi 1  yi  hf(xi ,yi ) Se calcula la derivada f ( x i  1 , y i  1 ) Se establece la derivada promedio (llamémosla B) B

1 f(x i , y i )  f(x i  1 , y i  1 ) 2





Se sustituye f(xi,yi) con este valor promedio en la ecuación de iteración de Euler y se obtiene y i 1  y i 

h  f(x i , y i )  f(x i 1 , y i 1 )   2

o simplificado: yi+1= yi+ hB Ejemplo 1: Resuelva f(x,y) = x - y en el intervalo [0,1] Con y(0)=2

x0=0

y0=2

xn=1

yn=?

y(1)=?

h  n=5

y Con 5 intervalos.

x

n

 x n

0

95



1 0 5

h = 0 .2

0

0.2

0.4

0.6

0.8

x0

x1

x2

x3

x4

1 x5

En forma iterativa, hagamos la siguiente tabla: de i=0 hasta n=5 Recuerde que:

xi+1=xi+h

f(xi,yi) = xi - yi

yi1  yi  hf(xi ,yi ) i

xi

0 0

Yi

2

f(xi,yi)

y i1

B



1 f(xi , y i )  f(xi 1 , y i 1 ) 2

Xi+1

f(x i 1 , y i 1 ) B

 Yi+1=yi+h b

-2

1.6

0.2

-1.4

-1.7

1.66

1 0.2 1.66

-1.46

1.368

0.4

-0.968

-1.214

1.4172

2 0.4 1.4172

-1.0172

1.21376

0.6

-0.61378

-0.81548

1.254104

3 0.6 1.254104

-0.654104

1.12328

0.8

-0.32328

-0.488694

1.156365

4 0.8 1.156365

-0.356365

1.085092 1.0

-0.08509

-0.2207286

1.112220

5 1.0 1.112222

Para i=5

y5=y(1)=1.112222 El valor exacto es 1.10364

1.10364  1.112222  100% 1.10364 E R  0.78%

 ER 

8.3.

MÉTODOS DE RUNGE-KUTTA (RK)

Los métodos de Runge-Kutta tienen la exactitud del esquema de la serie de Taylor sin necesidad de cálculo de derivadas superiores y consisten en obtener un resultado que se obtendría al utilizar un número finito de términos en una serie de Taylor de la forma: yi 1  yi  f ( xi , yi )h 

f ' ( xi , yi )h 2 f " ( xi , yi )h 3 f ' ' ' ( xi , yi )h 4   2! 3! 4!

+.......... con una aproximación en la cual se calcula yi+1 de la fórmula:

Yi 1  yi  h( 0 f ( xi , yi )  1 f ( xi  1h, yi  b1h)  ... ....   p f ( xi   p h, yi  bp h)

96

(2)

donde los

,  ,

b

se determinan de modo que si expandiera

f ( xi   j h, yi  bj h) con i  j  p en series de Taylor alrededor de (xi, yi), se observaría

que los coeficientes de h, h2 , h3 , etc. Coincidirían con los

coeficientes correspondientes de la ecuación (1). Después de efectuar los cálculos convenientes a la ecuación (2) para establecer los valores de  ,  , b , tenemos las siguientes expresiones de Runge-Kutta (RK):

1. Método RK de segundo orden: yi+1 = yi +

h [k1 +k2] 2

donde: k1=f(xi, yi) k2=f(xi+h, yi+hk1) 2. Método RK de tercer orden: yi+1 = yi + donde:

h [k1 +4k2 +k3] 6

k1=f(xi, yi) k2=f(xi+

h h , yi+ k1) 2 2

k3= f(xi+h, yi -hk1+2hk2) 3. MétodoRK de cuarto orden: yi+1 = yi +

h [k1 +2k2 +2k3 +k4] 6

donde: k1=f(xi, yi) k2=f(xi+

h h , yi+ k1) 2 2

k3=f(xi+

h h , yi+ k2) 2 2

k4=f(xi+h, yi+hk3) otra forma de expresar

97

k1  hf ( x n , y n ) k h k 2  hf ( x n  , y n  1 ) 2 2 k h k 3  hf ( x n  , y n  2 ) 2 2 k 4  hf ( x n  h, y n  k 3 ) 1 y n 1  y n  (k1  2k 2  2k 3  k 4 ) 6



El método más utilizado es el de cuarto orden ya que coincide con los primeros cinco términos de la serie de Taylor lo cual significa gran exactitud sin cálculo de derivadas aunque haya que evaluar la función f(x,y) cuatro veces en cada subintervalo.

Ejemplo: Aplíquese el método de Runge-Kutta de cuarto orden para f(x,y)=x-y donde x=0 hasta x=1 con 5 intervalos y siendo para x=0, y=2.

Ejemplo: 1.

Aplíquese el método de Runge-Kutta de cuarto orden para f(x,y)= x – y desde x=0 hasta x=1 Con 5 intervalos siendo para x=0 y=2. Solución: h

(x0, y0)= (0, 2)

1 0  0, 2 5

para x=1 y=? x0=0

x1=0.2

x2=0.4

Aplicamos: yi+1 = yi +

x3=0.6

x4=0.8

x5=1.0

h [k1 +2k2 +2k3 +k4] 6

Se calculan los valores de k1, k2, k3, k4 en cada iteración y se halla yi+1

i xi

Yi

K1

K2

K3

K4

Yi+1

0 0

2

-2

-1.7

-1.73

-1.454

1.6562

1 0.2 1.6562

-1.4562

-121058 -1.235142 -1.009172 1.410973

2 0.4 1.410973 -1.010973 -0.809876 -0.829985 -0.644976 1.246451 3 0.6 1.246451 -0.746451 -0.471806 -0.49927 -0.346797 1.148004 4 0.8 1.148004 -0.348004 -0.213204 -0.226684 -0.10268 1.103656 5 1.0 1.103656 98

para x=1.0, y=1.103656 Que al comparar con el valor exacto que es 1.10364, tenemos: ER 

2.

1.10364  1.103656  100%  0.001% 1.10364

Resolver F ( X i , Yi ) 

1 1  X i  Yi 2 aplicando el método de Runge-Kutta. 2

Solución De la condición inicial del problema se tiene que X = 0, y Y = 1; además, h = 0.1. Sustituyendo estos valores en se obtiene:

Llevando estos valores a (16) y el resultante a (12) se obtiene que para X = 0.1 la solución del problema es

Los valores de las kipara este punto obtenido de la solución, son:

luego

Continuando de la misma forma se obtiene la solución que se muestra en la siguiente tabla: 99

8.4.

X

Y

k1

k2

k3

k4

0.0

1.0000

0.5000

0.5516

0.5544

0.6127

0.1

1.0554

0.6126

0.6782

0.6823

0.7575

0.2

1.1236

0.7575

0.8431

0.8494

0.9494

0.3

1.2085

0.9492

1.0647

1.0745

1.2121

0.4

1.3158

1.2119

1.3735

1.3896

1.5872

0.5

1.4545

1.5868

1.8234

1.8517

2.1509

MÉTODO DE DIFERENCIAS FINITAS

Las diferencias

divididas finitas se sustituyen por las derivadas en la

ecuación original. Así una ecuación diferencial se transforma en un conjunto de ecuaciones algebraicas lineales simultaneas. Primera derivada:

dy yi1  yi  dx x

Segunda derivada:

d 2 y yi 1  2 yi  yi 1  dx 2 x 2

Resolver la ecuación mediante diferencias finitas

d2 y  2  106 y  7.5  107  (75  x) 2 dx i 1

i 1

i

d 2 y yi 1  2 yi  yi 1  dx 2 (x) 2

yi 1  2 yi  yi 1  2  106 yi  7.5  107 xi (75  xi ) 2 (x) i 1 x0

i2

i 3

i4

x  25

x  50

x  75

x  25 , x  0 a x  75 con x  25 .

100

x0  0 x1  x0  x  0  25  25 x 2  x1  x  25  25  50

x3  x2  x  50  25  75 x  0 y1  0

y3  2 y 2  y1  2  10 6 y 2  7.5  10 7 x2 (75  x2 ) 2 (25) 0.0016 y1  0.003202 y 2  0.0016 y 3  7.5  10 7 ( 25)(75  25) 0.0016 y1  0.003202 y 2  0.0016 y 3  9.375  10 4

y 4  2 y3  y 2  2  10 6 y3  7.5  10 7 x3 (75  x3 ) 2 (25) 0.0016 y 2  0.003202 y 3  0.0016 y 3  7.5  10 7 (50 )(75  50 ) 0.0016 y 2  0.003202 y 3  0.0016 y 3  9.375  10 4

x  75

y4  0  0 0 0   y1  0  1  y   0.0016  0.003202  0.0016 0   2  9.375  10 4     0 0.0016  0.003202 0.0016  y3  9.375  10 4       0 0 1   y 4  0  0   y1  0  y     2    0.5852  y 3   0.5852       y 4  0

8.5.

SISTEMA DE ECUACIONES DIFERENCIALES DE PRIMER ORDEN

Ya sea que estemos afrontando un problema de ingeniería cuya solución implique la resolución de una ecuación diferencial de orden n, o uno que implique la resolución de un sistema de ecuaciones diferenciales ordinarias de primer orden, nos enfrentaremos a la necesidad de resolver un sistema que podremos expresar como sigue:

101

dy1  f1  x , y1 , y 2 ,............, y n  dx dy 2  f 2  x , y1 , y 2 ,............, y n  dx . . dy n  f n  x , y1 , y 2 ,............, y n  dx

Por supuesto, la solución de tal sistema requiere de que se conozcan las n condiciones iniciales en el valor inicial del intervalo correspondiente a x. Todos los métodos vistos anteriormente para simples ecuaciones pueden extenderse a la resolución de sistemas como el anterior. El procedimiento para resolver un sistema de ecuaciones simplemente involucra aplicar las técnicas conocidas para cada ecuación en cada paso, antes de proceder con el siguiente. Esto quedará claramente ilustrado con el siguiente ejemplo donde hemos aplicado el método de Euler.

Ejemplo 1:Resolución de un sistema de EDOs mediante el método de Euler. Resolver el siguiente conjunto de ecuaciones diferenciales ordinarias: dy1  - 0, 5  y1 dx dy 2  4 - 0, 3  y 2 - 0,1  y1 dx

Resolveremos el sistema en el intervalo entre x = 0 y x = 2, con condiciones iniciales en x = 0, y1 = 4 , y2 = 6. Utilizaremos un paso h = 0.5. Se implementa el método de Euler para cada variable mediante la ya conocida expresión:

yi 1  yi  f ( xi , yi )  h Primero calculamos las pendientes: dy1  f1 ( 0 , 4, 6)  - 0.5  4  - 2 dx dy 2  f 2 ( 0 , 4, 6)  4 - 0.3  6 - 0.1  4  1.8 dx

Y luego los valores de la función para el primer paso:

102

y1  0.5  y1  0  f1 ( 0 , 4 , 6)  h  4   -2 0.5  3 y2  0.5  y2  0  f 2 ( 0 , 4 , 6)  h  6  1.8 0.5  6,9 Para un segundo paso volvemos a calcular las pendientes: dy1  f1 ( 0.5 , 3, 6.9)  - 0.5  3  - 1.5 dx dy 2  f 2 ( 0.5 , 3, 6.9)  4 - 0.3  6.9 - 0.1  3  1.63 dx

Y luego los valores de la función para el segundo paso:

y1 1.0  y1  0.5  f1 ( 0.5 , 3, 6.9)  h  3   -1.5  0.5  2.25 y2 1.0  y2  0.5  f2 ( 0.5 , 3, 6.9)  h  6.9  1.63  0.5  7.715 Y así continúa el cálculo hasta el final. Los resultados se resumen en la siguiente tabla: x

y1

y2

0.0

4.000000

6.000000

0.5

3.000000

6.900000

1.0

2.250000

7.715000

1.5

1.687500

8.445250

2.0

1.265625

9.094870

8.6 PROBLEMAS

1. Dada la ecuación diferencial f(x,y)=yx2-y con x=0 a x=2

siendo y(0)=1 y

n=4 intervalos resuelva por el método sencillo de Euler 2. Dada y’=f(x,y) = 2x3 –3x2 entre x=0 y x=1 siendo f(0)=1 y con n=2, evalúe utilizando a)el método sencillo de Euler. 3. Dada y’=f(x,y)= -2x3 + 12x2 –20x + 8.5 resuelva utilizando el método RK de cuarto orden desde x=1 hasta x=3 considerando 4 intervalos y siendo y(1)=3. 4. Utilice el método de Euler Modificado para resolver: a) dy/dx=2x3 –2xy con y(0)=0

y(2.5)=? 103

h=0.5

b) y’=2x2 –3y2 con y(1)=0.5

y(2)=?

h=0.2

5. Resolver mediante el método euler dy  100000( y  e  x )  e  x dx

y (0)  0 , use the implicit Euler method to obtain the solution from x  0 to 0,3

using the step-size of 0.1 6. Un tanque de 600 galones de capacidad contiene inicialmente 200 galones de salmuera con 25 lb. de sal. Salmuera con 2 lb. por galón entra al tanque a un flujo de 13 galones/ s. La salmuera mezclada en el tanque fluye hacia fuera a un flujo volumétrico de 8 galones/ s. ¿Cuál es la cantidad de sal y la concentración cuando el tanque se Encuentra lleno? Utilizar el método de runge kutta de segundo orden. 7. Las plantas de desalinización se usan para purificar el agua de mar , para que pueda beberse . El agua de mar contiene disuelto 8 g de sal /kg. Y es bombeada hacia un tanque de mezcla a razón de 0,5 kg./ min. Debido a una falla en el diseño del equipo, el agua se evapora del tanque a razón de 0,5kg/min. .La solución salina sale del tanque a razón de 10 kg/min. Supóngase que el balance de la disolución es agua pura. Si el tanque se llena inicialmente con 1 000 kg. de solución a la entrada. ¿ determine la concentración a la salida del tanque en función del tiempo hasta que el tanque quede vació.? 8. Las siguientes ecuaciones definen las concentraciones de tres reactores

dCA  20CACC  2CB dt dCB  20CACC  2CB dt dCC  20CACC  2CB  0,2CC dt Si las condiciones iniciales son CA = 500, CB = 0 y Cc = 500 halle las Concentraciones para los tiempos que van desde 0 a 30s usando un h = 3

9. La descarga de un fluido por gravedad por un orificio por el fon deo del tanque puede ser modelado mediante el siguiente sistema d e ecuaciones diferenciales. 104

dv  0,0010 y  0,00205 v 2 dt dy 1  0 , 0624 (1  ) dt 100  v

donde v es la velocidad del fluido, y nivel del liquido las condiciones iniciales son de v0 = 2,5 ; y0 = 20 cual son los valores de v, y para t = 5 con h igual 0,5 10.

Par un circuito eléctrico los resistores no obedecen la ley de Ohm, y el

circuito dinámico se describe por la relación siguiente l

di i 3 i  R ( )   0 dt I  I

i = corriente eléctrica ; I = corriente de referencia igual 1; R resistencia eléctrica 2 ; Resolver la ecuación diferencial si i(0) = 0,6 A para

t=

30s con h = 3s 11.

Determine usado técnicas numéricas

la ecuación diferencial ordinaria

de segundo orden

d2x dx  3  10 x  0 2 dt dt

La condición inicial para la ecuación para t = 0

 x =1

Determine los valores de x entre 1,2 usando h = 0,1 12.

Un tanque

de 2 000 L contiene, en el inicio 400L de agua pura .

Comenzando en t = 0 una solución acusa que contiene 1 g/ L de cloruro de potasio fluye hacia el tanque a razón de 8L/s y al mismo tiempo, comienza a fluir una corriente de salida a razón de 4L/s .El contenido del tanque esta mezclado perfectamente y la densidad de la corriente de alimentación y salida puede considerarse constante.

Calcule la concentración de cloruro

de potasio en el tanque en el momento en que este rebosa.

13.

Un Reactor CSTR como se muestra en la figura, la concentración de

alimentación (CA0) of 0,5 mole/m3.

Determine la concentración al cabo de

0,5h

105

rA  V

 k1C A (1  k2 C A )

dC A  F (C A0  C A )  (V )(rA ) dt

V  2, 0m3

F  1, 0m3 / h

k1  1, 0h 1

k2  1, 0m3 / mole

14. Se tiene un intercambiador de calor de tubos concéntricos en contracorriente y sin cambio de fase .Las ecuaciones que describen el intercambiador de calor en ciertas condiciones de operación son: dTB  0.03 x TS  TB  dx dTS  0.04 x TS  TB  dx

CONDICIONES

TB (0)  20 0 C TS (0)  1000 C

TB (3)  ¿? TS (3)  ¿? Evalué la temperatura para x= 3 m

15.

Resolver la siguiente ecuación mediante diferencias finitas

d 2u 1 d u u   2  0, u(5)  0,08731, u(8)  0.0030769 dr 2 r dr r

106

16. El siguiente diagrama ( Figura.18) muestra tres tanques continuos con agitación conectados en serie

Se ha disuelto1500g de Na 2SO4 en el

primer tanque, llenando los

otros tanques con solvente puro, e

iniciando después un flujo de 40L/ a través del sistema. Calcular el tiempo para alcanzar una concentración de 0,01g/l en el tanque tres.

Fig Nº16: Sistema de tanques de mezclado.

17.

Determine el valor de y(1) resolviendo y t   0,05 y t   0,15 y t   0, y 0   0, y0  1

utilizando el método de Runge-Kutta de segundo orden, con h = 0,5

18. El circuito que se muestra en la figura E9.5 tiene una autoinductancia de L = 50H, una resistencia de R = 20 ohms y una fuente de voltaje de V = 10 volts. Si el interruptor se cierra en el instante t = 0, la corriente I(t) satisface la ecuación: L

d I t   RI t   E , dt

I 0   0

Determine el valor de la corriente para 0 < t  10 segundos, mediante el método de Runge-Kutta de segundo orden, con h = 0.1. L SW

i

E

R

Figura 19: Circuito electrico

107

19. El agua que sale del tanque entra a otro de 20 galones, en cual también se vierte agua pura a razón de 3 galones/minuto y se mezcla bien. La concentración de sal en el segundo tanque satisface y 2 t   

3 2 y 2 t   y1 , 20 20

y 2 0   10

donde y1(t) es la concentración de sal del tanque de 50 galones del problema anterior. Utilice el método de Runge–Kutta de segundo orden para determinar cuando alcanza su máximo la concentración de sal en el tanque de 20 galones. Suponga que el segundo tanque tiene agua pura en el instante t = 0.

20. Un tanque de 50 galones de agua contiene sal con una concentración de 10 onzas/galón. Con el fin de diluir el contenido de sal, se suministra agua pura a razón de 2 galones/minuto. Si el depósito tiene una mezcla uniforme y la misma cantidad de agua que entra sale del depósito cada minuto, la concentración de sal satisface y1 t   

2 y1 , 50

y1 0   10

donde y1(t) es la concentración de sal en onzas/galón y t es el tiempo en minutos. Utilice el método de Runge–Kutta de segundo orden con h = 1 minuto para determinar cuánto tiempo debe transcurrir para que la concentración de la sal sea 1/10 de su valor inicial. b) El agua que sale del tanque entra a otro de 20 galones, en cual tambien se vierte agua pura a razon de 3 galones/minuto y se mezcla bien. La concentracion de sal en el segundo tanque satisface y 2 t   

3 2 y 2 t   y1 , 20 20

y 2 0   10

donde y1(t) es la concentración de sal del tanque de 50 galones del problema anterior. Utilice el método de Runge–Kutta de segundo orden para determinar cuando alcanza su máximo la concentración de sal en el tanque de 20 galones. Suponga que el segundo tanque tiene agua pura en el instante t = 0.

108

21. Resolver el sistema de EDOs mediante el método de Runge – Kutta

dy1  - 0.5  y1 dx clásico de cuarto orden.

dy 2  4 - 0.3  y 2 - 0.1  y1 dx Resolveremos el sistema en el intervalo entre x = 0 y x = 2, con condiciones iniciales en x = 0, y1 = 4 , y2 = 6. Utilizaremos un paso h = 0.5.

109

9. ECUACIONES DIFERENCIALES PARCIALES (EDP) INTRODUCCIÓN Una EDP es una ecuación que tiene como incógnita a una función de dos o más variables y que involucra a una o más de sus derivadas parciales. El orden de una EDP es el de la derivada con mayor orden en la ecuación. La linealidad de las ecuaciones se establece como sigue: Si los coeficientes dependen sólo de las variables independientes entonces a la ecuación se le denomina lineal. Si además dependen de la propia función o de alguna de sus derivadas parciales entonces la ecuación es no lineal. EJEMPLO 

Ecuación diferencial parcial de transferencia de calor unidimensional en estado no estacionario.

T T 2  2 t x T: Temperatura en K , t es el tiempo en s,  es la difusivilidad térmica en m2/s



Ecuación de transferencia de calor en estado estacionario bidimensional T 2 T 2  2 x 2 y T es la temperatura en K , x,y ejes coordenados

9.1 PROBLEMAS 1. resuelva la ecuación diferencial parcial de transferencia de calor en dos dimensiones.

 2T  2T  0 X 2 Y 2

0  X 2 ,

0  Y  2 , X  Y  0,5

T( 0, Y) = 100 , T( X ,0 ) = 60 , T( 2,Y) = 100 , T(X, 2) = 60

2. Resuelva la siguiente ecuación diferencial de transferencia de calor dada por 110

 2T  2T  0 X 2 Y 2 T = 40 .

0  y 1 ,

x=0

T = 120

0 y  1

x =4

T = 40

0 X  4

y=0

T= 50

0 X  4

y=1

Utilizando y  0,25 y x  1

111

10. APLICACIONES DE INGENIERÍA QUÍMICA EN POLYMATHMATHCAD Y MATLAB 10.1 PROBLEMAS CON POLYMATH

1. Determine el volumen Molar del gas dióxido de carbono mediante la ecuación de VAN DER WAALS . a   P  V  b   RT V *V  

a= 3,592

b= 0,04267

R = 0,082

PARA T= 300K , P =100 ATM f(v) = (P + a/(V2))*(V-b) –R*T = 0 POLYMATH Results NLE Solution Variable V

Value 0.0793969

a

3.592

b

0.04267

T

300

R

0.082

P

100

f(x) -1.776E-13

Ini Guess 0.2505

2. En la batería de reactores de mezcla completa que se muestra en el diagrama, se lleva a cabo una reacción isotermica de Segundo orden a volumen constante si el flujo volumétrico a cada reactor es constante de 25 l/min. calcular la concentración de salida de cada reactor

112

Figura Nº17: Bateria reactores

Balance de materia en cada reactor: Reactor N°1: q*CA0 –q*CA1 –K*CA12*V1 =0 Reactor N°2 q*CA1-q*CA2 –K*CA22*V2 =0 Reactor N°3 q*CA2-q*CA3 –K*CA32*V3 =0 F(CA1) = q*(CA0 –CA1) –K*CA12*V1 =0 F(CA2) = q*(CA1-q*CA2) –K*CA22*V2 = F(CA3) = q*(CA2-q*CA3) –K*CA32*V3 =0 Variable

Value

f(x)

CA1

0.6031317

-7.105E-15

2

CA2

0.2869484

-1.332E-14

1

CA3

0.1725794

-2.944E-13

0.8

V

1200

K

0.08

CA0

2

Q

3. En

s

25

un reactor químico se obtiene

glucosa a partir de la hidrólisis de

almidón. La ecuación diferencial que gobierna el proceso para la concentración del almidón es .

V*Dca/dt = q*CA0 - q*CA + ra ra = -k*CA2V CA0: Concentración de la alimentación a la entrada = 6,3 mol/m3 113

CA: Concentración del almidon a la salida V : volumen del reactor = 500m3 Q : caudal volumetrico 100 m3/h ra : velocidad de transformación de almidon en

glucosa

3

K= 0.02 m /h*mol . Calcular la composición de salida para t = 20h ODE Report (RKF45) Differential equations as entered by the [1] d(CA)/d(t) = q*CA0 /V - q *CA/V -K*CA2 Explicit equations as entered by the [1] q = 100 [2] V = 500 [3] CA0 = 6.3 [4] K = 0.02 Independent variable variable name : t initial value : 0 final value : 20 Variableinitial valueminimal valuemaximal valuefinal value t

0

CA

6.3

q

20

20

4.38179

6.3

4.38179

100

100

100

100

V

500

500

500

500

CA0

6.3

6.3

6.3

6.3

0.02

0.02

0.02

0.02

K

0

4. La concentración saturada de oxigeno disuelto en agua como función de temperatura y de la concentración de cloro se lista en la siguiente tabla. *Determine la ecuación que correlacione C=f(T) de grado 2 y 3 

Evaluar la concentración para T = 20°C

114

Temperatura

Cloro( 10 000mg/l)

5

11,6

10

10.3

15

9,1

20

8,2

25

7,4

30

6,8

REGRESION LINEAL Model: C = a0 + a1*T + a2*T2 Variable Value 95% confidence a0

13.11

0.1595539

a1

-0.3195

0.0208769

a2

0.0036429

5.839E-04

REGRESION POLINOMICAt Model: C = a0 + a1*T + a2*T2 + a3*T3 Variable Value 95% confidence a0

13.133333

0.5263102

a1

-0.3253704

0.1199362

a2

0.0040317

0.0076758

a3

-7.407E-06

1.451E-04

5. La siguiente secuencia de reacciones se llevan a cabo en estado no estacionario según la reacción:

Determinar la ecuación de cada componente de la reacción así como los productos en función al tiempo. Construya la tabla tiempo vs concentración en un paso de 1 a 1 hasta 10 minutos, asuma que las concentraciones iniciales son A (0)=1 mol/L, B (0)=1 mol/L, C (0)=1 mol/L. Datos: 115

K1=0.5 min-1 K2=1.0 min-1 Solución: formando las ecuaciones diferenciales en función de la reacción. C1   k1C1 t C2  k1C1  k2C2 t C2  k 2 C2 t

FIGURA 18: Datos de ingreso en Polymath

FIGURA 19: Gráfica concentración vs tiempo 116

6. Resolver el siguiente sistema: y   k1 x  k2 y t x  k1 x t z  k2 y t

Con valores iníciales x (0)=1, y (0)=0, z(0)=0 con t=0 a t=3 y k 1=1, k2=2 Solución:

FIGURA 20: Ingreso de Datos en Polymath

117

FIGURA 21: Resultados con el Software Polymath

FIGURA 22: Grafica de resultados de la Ecuación diferencial. 118

7. Resolver las siguientes ecuaciones:

dC1  0,12C1  0, 02C3  1 dt dC2  0,15C1  0,15C2 dt dC3  0, 025C2  0, 225C3  4 dt dC4  0,1C3  0,1375C4  0, 025C5 dt dC5  0, 01C1  0, 01C2  0, 04C5 dt

Con C1 (0)=1, C2 (0)=1, C3 (0)=1, C4 (0)=1 y C5 (0)=1 de t=0 a t=1 Solución:

FIGURA 23: Ingreso de Datos en Polymath

119

FIGURA 24: Resultados Polymath

FIGURA 25: Grafico Polymath

120

10.2 PROBLEMAS CON MATHCAD

1.

Para flujo turbulento para todas las tuberías, el Instituto de Hidráulica de Estados Unidos y la mayoría de ingenieros consideran la ecuación de Colebrook como la más aceptable para calcular f. La ecuación es:

 e 0, 251  0,86 ln    0,36d Re f f 

1

  

Donde: Re

: l número de Reynolds (adimensional)

e

: aspereza o rugosidad de la tubería (unidad de longitud)

d

:, diámetro de la tubería (unidad de longitud)

Obtener el factor de fricción para un fluido con un Reynolds de 3E4 que fluye en una tubería con un diámetro de 0,1 m y una rugosidad de 2,x10-3 m. Solución con mathcad Ingresando datos: E = 0,002 D = 0,1 Re = 3104 f = 0,001 Given

 E 0, 251  0,86 ln    f  0,36 D Re f

1

  

Calculo de la fricción F = Find (f) f = 0,055 de la figura efectuar un balance de masa y determine los flujos

121

Figura Nº26: Sistema de columna de separación

0.07x1 +0.18x2 + 0.15x3 +0.24x4

= 10.5

0.04x1 + 0.24x2 + 0.1x3 +0.65x4 =17.5 0.54x1 + 0.42x2

+ 0.54x3 + 0.10x4 = 28

0.35x1 + 0.16x2 + 0.21x3 + 0.01x4 = 14  0.07 0.04 A    0.54  0.35 

0.18 0.15 0.24  0.24



0.1 0.65 

  0.160 0.21 0.01  0.42 0.54 0.1

122

 10.5     26.25  17.5     B   17.5  solucion   28   8.75   14   17.5      2.

La ecuación de Antonie para el cálculo de la presión de vapor (P*) de un compuesto esta dado por: Donde: P*: presión de vapor (mmHg) T, temperatura (ºC)

Las constantes para el benceno y el hexano son:

Compuesto

A

B

C

benceno

6.89745

1260.350

220.237

Hexano

6.87773

1171.530

224.366

Determinar la temperatura de ebullición de la mezcla a 1 y a 5 atm. Solución

log ( P*)

A

B C T

sean las fracciones molares xa = 0,5 xb = 0,5 constantes para el benceno A1  6.89745

B1  1206.350 C1  220.237

constantes para el hexano A2  6.87773

B2  1171.530 C2  224.366

para una presión de 1atm en mmHg 123

Pt  760

asumiendo una temperatura inicial de 30ºC T  30 Given

 A1 B1   A2 B2      C1  T  C2  T    Pt Xa 10  Xb 10 T  Find(T) T  73.975

3.

Resolver el siguiente sistema de ecuaciones no lineales con mathcad

2

k1

k2

x y

2

( 1  x)  ( 1  x  y ) y ( x  y) ( x  y)  ( 1  x  y)

Solución

Especifique las constantes de equilibrio Tome dos valores iniciales Given 2

k1

k2

x  y

2

( 1  x)  ( 1  x  y ) y (x  y ) (x  y )(1  x  y )

vec  Find ( x  y ) vec 

 0.652     0.198  124

k1  0.6

k2  0.2

x  0.8

y  0.5

4..

Mediante la ecuacion de Vand Der Waals calcular el volumen molar a P = 60 atm y T = 500K.

 P  a   ( V  b ) R T  2 V  

ecuacion de Van der Waals

Definiendo las constantes R  0.082 Para el amoniaco

Tc  405

Pc  111.3

Constante a y b 27  R  Tc 2

a   64 

2

 

Pc

b 

R Tc 8 Pc

Especificando las constantes R  0.082

Para el amoniaco

Constante a y b a 

Tc  405

 R 2  Tc  64  Pc 27

Evaluando las constantes a y b

2

a 

Especifique la presion y temperatura

b 

R  Tc 8  Pc b 

T  500

Asumiendo un valor inicial( comportamiento ideal) R T V  V  P Dado

 P  a  (V  b )  2  V   calculo del volumen Molar V  Find

  

Pc  111.3

(V)

V 

125

R T

P  60

5.

Encontrar el volumen molar que ocupa el gas empleando la ecuación de Redlinh Kwong RT a P  V  b V (V  b) T  R 2TC5 / 2   RTc  a  0, 42747   b  0, 08664    Pc   PC 

 Atm  L    mol  gK 

P en Atm V en L/mplg T en K R  0, 082 

Pc: presión critica de amoniaco 111,3 atm. TC: temperatura critica de amoniaco 405K ECUACION DE ESTADO DE R-K Ingresando las contantes: Ingresando los datos:

R  0.082 P  60

Pc  111.3

Tc  405

T  500

Calculo de a yb: 5   R 2  Tc 2 a  0.42747   Pc

  

b  0.08664 



a 

b 

Asumiendo comportamiento ideal

V  V

Dado P

R T (V  b)



V  Find ( V)

a V ( V  b )  T V

126

R T P

R  Tc Pc

  

6.

A un sistema de separación FLASH ( Fig. 21) ingresa 100 moles/h de una mezcla

VAPOR(V)

2

yi

1 Alimentacion F zi

3

xi

Liquido(L) Figura Nº27:Sistema de separación de componentes Alimentación al Flash( Mol/h) Zi: Fraccion molar de los componentes de la mezcla i = 1,2,3, …NC V: Caudal de salida del vapor ( mol/h) yi : Fraccion molar de los componentes en la fase vapor L: Caudal De salida de la fase liquida

( mol/h)

xi: Fraccion molar de los componentes en la fase liquida Balance de materia global: F = V + L Balance de componentes : Zi*F = xi*L + yi*V Equilibrio Liquido vapor : yi = Ki *xi V F Remplazando en las anteriores ecuaciones expresiones

Fraccion vaporizada :

x



zi  ( Ki  1)  1 127

Se obtiene las siguientes

yi 

Ki * zi  ( Ki  1)  1

Ni

Zi( Ki  1)

  (Ki  1)  1  0 i 1

Teniendo en cuenta los datos de la siguiente tabla Calcular la fracción vaporizada así como los flujos en cada corriente con sus respectivas composiciones COMPONENTES

FRACCION MOLAR

C1

0.05

16.25

C2

0.15

5.25

C3

0.25

1.99

C4

0.2

0.75

C5

0.35

0.29

Ingresando los datos:

F  100 k3  1.99

k1  16.25 k4  0.75

z1  0.05 z3  0.25

K

k2  5.25 k5  0.29

z2  0.15 z4  0.2

Asumiendo un valor:

z5  0.35

  0.5

Dado z1 ( k1  1)  ( k1  1)  1



  Find( )

z2 ( k2  1)  ( k2  1)  1



z3 ( k3  1)  ( k3  1)  1



128



z4 ( k4  1)  ( k4  1)  1



z5 ( k5  1)  ( k5  1)  1

0

7.

Encontrar el volumen molar que ocupa el gas empleando la ecuación de estado de Redlinh-Kwong.

A una presión de 14 bar y una temperatura de 333K

RT a  V  b V (V  b) T

P

Ingresando los datos: T  333 a b

f ( V)

P T

b  44.897

8

a  1.5641410 

R  83.14

P  14

  b  2

b  R T



P



  V  R T  V2  V3  P P T  a

Se observa un plinomio de grado 3 Ingrese los coeficientes en una matriz de 4*1

     2  b  V        

  P T  b  R T a    P P T   R T    P  1  a b

S  polyroots ( V)

8.

 71.299    229.996 S   3  1.676  10 

La temperatura media logarítmica de un intercambiador de calor a contracorriente esta dado por:

LMDT 

T1  t2   T2  t1  T t  ln  1 2   T2  t1 

129

Para este sistema la temperatura media logarítmica debe ser de 50°C y el fluido caliente se alimenta a 100°C y sale a 60 °C, mientras que el fluido frió se alimenta a 15°C ¿A qué temperatura sale del intercambiador el fluido frió?

Ingrese los datos : T1  15

T3  100

T4  60

Asumiendo un valor inicial:

LMTD  50

T2  50

Dado LMTD

( T3  T2)  ( T4  T1) ln

T3  T2 

  T4  T1 

T2  Find( T2) T2 

9.

El factor de fricción para diferentes números de Reynolds empleando la ecuación de Colebrook

  1 2,51   0,86 ln    f  3, 7d Re f  Obtener el factor de friccion para Re de 10 000 , e = 0,002cm y diámetro de 2cm CALCULO DEL FACTOR DE FRICCION

Ingresando datos:

Re  10000

Asumiendo para f : f  0.001 Dado 0.86 ln

1

e

 3.7 d

f



  Re f  2.51

f  Find( f ) f 130

e  0.002

d  2

10.3. PROBLEMAS CON MATLAB 1. A un reactor ingresa una mezcla de gases de 8% de dióxido de azufre y el 12% de oxigeno y la diferencia de nitrógeno, y se desarrolla la siguiente reacción. 1 SO2  O2  SO3 2

Calcular la composición en el equilibrio a presión constante de 2 atm y la constante de equilibrio KP es de 160. Algorritmo de Solucion en Matlab: clear all, close all, clc disp('H2SO4'); SO2=input('Ingrese moles de SO2: '); O2=input('Ingrese moles de O2: '); SO3=input('Ingrese moles de SO3: '); N2=input('Ingrese moles de N2: '); P=input('Ingrese Presion Total (atm): '); KP=input('Ingrese Constante de Equilibrio: '); syms x; disp('Moles en Equilibrio: '); ESO2=SO2-x; EO2=O2-0.5*x; ESO3=SO3+x; disp(ESO2); disp(EO2); disp(ESO3); ET=ESO2+EO2+ESO3+N2; f=KP-((ESO3/ET)*P)/((ESO2/ET)*((EO2/ET)^0.5)*P^1.5); 131

xai=input('Ingrese Limite Inferior: '); xbi=input('Ingrese Limite Superior: '); tol=input('Ingrese Tolerancia: '); f=inline(f); i=1; ea(1)=100; if f(xai)*f(xbi) < 0 xa(1)=xai; xb(1)=xbi; f1=f(xai); xr(1)=(xa(1)+xb(1))/2; fprintf('Iter.

x1

f1

x2

Error aprox \n');

fprintf('%2d \t %11.7f \t %11.7f \t %11.7f \n',i,xa(i),f1(i),xr(i)); while abs(ea(i)) >= tol, if f(xa(i))*f(xr(i))< 0 xa(i+1)=xa(i); xb(i+1)=xr(i); end if f(xa(i))*f(xr(i))> 0 xa(i+1)=xr(i); xb(i+1)=xb(i); end xr(i+1)=(xa(i+1)+xb(i+1))/2; ea(i+1)=abs((xr(i+1)-xr(i))/(xr(i+1))*100); fprintf('%2d \t %11.7f \t %11.7f \t %11.7f \t %7.3f \n',... i+1,xa(i+1),xr(i+1),xb(i+1),ea(i+1)); i=i+1; 132

end else fprintf('No existe una raíz en este intervalo'); end Resultados: TABLA 12 Resultados en Matlab

2. Se sabe que la temperatura media logarítmica de un intercambiador de calor contracorriente está dada por:

 ΔT  AQ

T = LMTD =

ent cal

sald - T frio  - Tcalsald -T frioentr 

sald  Tcalent - T frio ln  sald entr  T -T frio  cal

133

  

Se tiene el siguiente intercambiador: entr = 15ºC T frio

INTERCAMBIADOR DE CALOR

sald = T frio

Tcalsald= 60ºC

t

Tcalent = 100ºC

Se sabe que en este sistema la temperatura media logarítmica debe ser de 50ºC. El fluido caliente se alimenta al sistema a 100ºC, mientras que el fluido frío se alimenta a 15ºC,. ¿a que temperatura sale del intercambiador de fluido frío? Algoritmo de Solución en Matlab: clear all, close all, clc disp ('Temperatura en un Intercambiador de Calor') xo=input('Temperatura Inicial ='); n=input ('Numero de Iteraciones='); salida=ones(n,3); % matiz de salida de datos for i=1:n x1=xo-[(45*exp((55-xo)/50)+xo-100)]/[1+(-45/50)*(exp((55-xo)/50))]; vsal=[xo;x1]; ea=[[abs((x1-xo)/x1)]]; % error xo=x1; salida(i,1)=i; salida(i,2)=x1; salida(i,3)=ea; end disp('Iter.

T(i)

Error');

disp(num2str(salida)); Resultados:

134

TABLA 13 Resultados Temperatura en Matlab

3. Un gas se encuentra a una presión absoluta de 13.76 bar y una temperatura de 333 K. Encontrar el volumen molar ocupa el gas empleando la ecuación de estado de Redlich-Kwong.

P=

RT a - 1/2 V - b T V V + b 

Para este compuesto las constantes son: P = 13.76 atm T = 333ºk a = 1.5614 x 108 (cm6 bar/(g mol)2 k1/2) b = 44.897 (cm3/ g mol) R = 83.4 (cm3 bar /g mol k) Algoritmo de Solución en Matlab: clear all, close all, clc disp('Ecuacion de Redlich-kwong 2'); a=input('Ingrese a: '); b=input('Ingrese b: '); T=input('Ingrese T(K): '); P=input('Ingrese P(bar): '); R=input('Ingrese R(cm^3.bar/K.mol): '); syms x 135

f=R*T/(x-b)-a/((sqrt(T))*x*(x+b))-P; x0=(R*T)/(P); disp('Primer punto inicial: '); disp(x0); xai=x0; xbi=input('Ingrese Limite Superior: '); tol=input('Ingrese Tolerancia: '); f=inline(f); i=1; ea(1)=100; if f(xai)*f(xbi) < 0 xa(1)=xai; xb(1)=xbi; f1=f(xai); xr(1)=(xa(1)+xb(1))/2; fprintf('Iter.

V1

F1

V2

Error aprox \n');

fprintf('%2d \t %11.7f \t %11.7f \t %11.7f \n',i,xa(i),f1(i),xr(i)); while abs(ea(i)) >= tol, if f(xa(i))*f(xr(i))< 0 xa(i+1)=xa(i); xb(i+1)=xr(i); end if f(xa(i))*f(xr(i))> 0 xa(i+1)=xr(i); xb(i+1)=xb(i); end xr(i+1)=(xa(i+1)+xb(i+1))/2; 136

ea(i+1)=abs((xr(i+1)-xr(i))/(xr(i+1))*100); fprintf('%2d \t %11.7f \t %11.7f \t %11.7f \t %7.3f \n',... i+1,xa(i+1),xr(i+1),xb(i+1),ea(i+1)); i=i+1; end else fprintf('No existe una raíz en este intervalo'); end Resultados: TABLA 14 Resultados volumen en Matlab

137

4. El factor de fricción (f) para el flujo turbulento en una tubería está dado por la correlación de Colebrook:

ε/D

1 = - 0.86 ln f

2.54

 +   3.4 Re

  f 

Donde Re = es el número de Reynolds (adimensional) , es la aspereza o rugosidad de la tubería (unidad de longitud) D, es el diámetro de la tubería (unidad de longitud) Obtener el factor de fricción para un fluido con un Reynolds de 3E4 que fluye en una tubería con un diámetro de 0.1 m y una rugosidad de 0.0025m. Algoritmo de Solución en Matlab: clear all, close all, clc disp('Friccion en un Flujo Turbulento'); Re=input('Ingrese el Numero de Reynolds: '); e=input('Ingrese la Rugosidad: '); D=input('Ingrese el Diametro de Tuberia: '); syms x; f=-sqrt(1/x)-0.86*log((e/D)/3.4+2.54/(Re*sqrt(x))); xai=input('Ingrese Limite Inferior del Intervalo: '); xbi=input('Ingrese Limite Superior del Intervalo: '); tol=input('Ingrese Tolerancia: '); f=inline(f); i=1; ea(1)=100; if f(xai)*f(xbi) < 0 xa(1)=xai; xb(1)=xbi; 138

f1=f(xai); xr(1)=(xa(1)+xb(1))/2; fprintf('Iter.

f1

F1

f2

Error aprox \n');

fprintf('%2d \t %11.7f \t %11.7f \t %11.7f \n',i,xa(i),f1(i),xr(i)); while abs(ea(i)) >= tol, if f(xa(i))*f(xr(i))< 0 xa(i+1)=xa(i); xb(i+1)=xr(i); end if f(xa(i))*f(xr(i))> 0 xa(i+1)=xr(i); xb(i+1)=xb(i); end xr(i+1)=(xa(i+1)+xb(i+1))/2; ea(i+1)=abs((xr(i+1)-xr(i))/(xr(i+1))*100); fprintf('%2d \t %11.7f \t %11.7f \t %11.7f \t %7.3f \n',... i+1,xa(i+1),xr(i+1),xb(i+1),ea(i+1)); i=i+1; end else fprintf('No existe una raíz en este intervalo'); end Resultados: TABLA 15 Resultados volumen en Matlab

139

5. Determinar volumen molar del oxigeno mediante la ecuación del VAN DER WAALS

a    P  2 v  b   RT v   p = 100 Atm., T = 700 K para un gas que tiene a = 1,36 b = 0,0318 Algoritmo de Solución en Matlab: clear all, close all, clc disp ('Ecuacion de Van der Walls') T=input('Ingrese T(K): '); P=input('Ingrese P(atm): '); R=input('Ingrese R(L.atm/K.mol): '); 140

disp('Volumen Inicial: '); xo=R*T/P; disp(xo); n=input ('Numero de Iteraciones='); salida=ones(n,3); % matiz de salida de datos for i=1:n x1=xo-(100*xo^3-60.58*xo^2+1.36*xo-0.043248)/(300*xo^2121.16*xo+1.36); vsal=[xo;x1]; ea=[[abs((x1-xo)/x1)]]; % error xo=x1; salida(i,1)=i; salida(i,2)=x1; salida(i,3)=ea; end disp('Iter. V(i)

Error');

disp(num2str(salida));

Resultados:

TABLA 16 Resultados volumen en Matlab

141

6. La siguiente secuencia de reacciones se llevan a cabo en estado no estacionario según la reacción: k1 k2 A   B  C Determinar la ecuación de cada componente de la reacción así como los productos en función al tiempo. Construya la tabla tiempo vs concentración en un paso de 1 a 1 hasta 10 minutos, asuma que las concentraciones iníciales son A (0)=1 mol/L, B (0)=1 mol/L, C (0)=1 mol/L. Datos: K1=0,5 min-1 K2=1,0 min-1 Solución: formando las ecuaciones diferenciales en función de la reacción. C1   k1C1 t C2  k1C1  k2C2 t C2  k 2 C2 t

Metafile: function[f]=rxnserie(t,C); 142

%constantes K(1)=0.5; K(2)=1; f(1)=-K(1)*C(1); f(2)=K(1)*C(1)- K(2)*C(2); f(3)=K(2)*C(2); f=f'; file: tr=[0:1:10];C0=[1,0,1]; [t,C]=ode45('rxnserie',tr,C0); [t,C] TABLA 17 Resustado de concentracion en Matlab

143

plot(t,C); xlabel('tiempo (min)'); ylabel('concentración (mol/L)');

FIGURA 28 Resultados volumen en Matlab

144

UNIVERSIDAD NACIONAL DEL CALLAO FACUTLAD DE INGENIERIA QUIMICA ESCUELA PROFESIONAL DE ING. QUIMICA

SILABO

I.

II.

DATOS GENERALES 1.1 Asignatura

:

Métodos Numéricos

1.2 Código de la Asignatura

:

IG-301

1.3 Semestre Académico

:

2010-A

1.4 Ciclo

:

V

1.5 Créditos

:

04

1.6 Horas de teoría

:

06

1.7 Horas de práctica

:

06

1.8 Duración

:

16 semanas

1.9 Pre-requisito

:

IG 102, FM 202

1.10 Profesor de Teoría

:

Ing. Juan Medina Collana

SUMILLA 2.1 Naturaleza de la Asignatura El curso de Métodos Numéricos es un curso teórico, práctico de carácter obligatorio que sirve como base para las asignaturas de ingeniería, debido a que permite realizar cálculos complejos como resultado del modelamiento de fenómenos físicos y químicos.

2.2 Síntesis del Contenido Ecuaciones algebraicas no lineales. Sistema de Ecuaciones algebraicas no lineales. Interpolación Polinómica. Diferenciación. Integración. Ecuaciones Diferenciales Ordinarias. Sistema de ecuaciones Diferenciales Ordinarias. Sistema de Ecuaciones algebraicas lineales. Ecuaciones Diferenciales parciales. Análisis de regresión lineal, no lineal y multivariable.

145

III.

OBJETIVOS 3.1 Objetivo General Al término del curso, el alumno será capaz de resolver modelos matemáticos complejos que resulten del modelamiento de los fenómenos físicos y químicos, aplicando técnicas numéricas.

3.2

Objetivos Específicos  Realizar el estudio comparativo entre los métodos analíticos y métodos numéricos.  Adquirir destreza en el dominio de las técnicas de solución numérica.  Aplicar las técnicas numéricas para resolver problemas diversos de Ingeniería Química que resultan del modelameinto de fenómenos físicos y químicos.

IV.

PROGRAMA DE CONTENIDOS

PRIMERA SEMANA Ecuaciones algebraicas no lineales. Método de la bisección. Método de Newton Raspón de primer y segundo orden. Aplicaciones a la Ing. Química. Método de la secante. Método del punto fijo. Método de falsa posición. Métodos cuasi Newton. Métodos de la convergencia. Aplicaciones a la Ingeniería Química.

SEGUNDA SEMANA Sistema de ecuaciones no lineales. Método del punto fijo. Método de Newton Raspón estándar. Método de Newton Raspón modificado. Métodos cuasi Newton. Aplicaciones a la Ingeniería Química.

TERCERA SEMANA Interpolación polinómica. Ecuaciones en diferencias progresivas, regresiva y central de Newton. Aplicaciones a la Ingeniería Química.

146

Generación de polinomios mediante las diferencias progresiva, regresiva y central de Gauss. Interpolación con estaciones no equidistantes.

CUARTA SEMANA Análisis de Regresión Lineal y no lineal. Mínimos cuadrados. Regresión multivariable. Aplicaciones a la Ingeniería Química

QUINTA SEMANA Diferenciación numérica. Aproximación por diferencias. Diferenciación de los polinomios de Newton. Aplicaciones a la Ingeniería Química.

SEXTA SEMANA Integración numérica. Regla del Trapecio. Reglas de Simpson 1/3 y 3/8 y

SEPTIMA SEMANA Integración por cuadratura Integración múltiple. Aplicaciones a la Ingeniería Química.

OCTAVA SEMANA EXAMEN PARCIAL

NOVENA SEMANA Ecuaciones Diferenciales Ordinarias. Método de Euler. Método de EulerGauss. Método de Taylor.

DECIMA SEMANA Método de Runge-Kutta. Aplicaciones a la Ingeniería Química. Métodos Predictor-Corrector. Métodos de paso variable. Aplicaciones a la Ingeniería Química.

147

DECIMA PRIMERA SEMANA Sistema de Ecuaciones Diferenciales Ordinarias. Ecuaciones Diferenciales de Orden Superior: métodos de Ruge Kutta y matricial, explícita en la variable dependiente. Aplicaciones a la Ingeniería Química. DECIMA SEGUNDA SEMANA Sistema de ecuaciones algebraicas lineales. Métodos de Jacobi y GaussSeidal. Matrices y vectores. Sistemas especiales. Aplicaciones a la Ingeniería Química. DECIMA TERCERA SEMANA Ecuaciones Diferenciales Parciales. Ecuaciones en Diferenciales. Métodos explícito e implícito. Aplicaciones a la Ingeniería Química. Aplicaciones a la Ingeniería Química.

DECIMA CUARTA SEMANA Practica calificada

EXAMEN FINAL

EXAMEN SUSTITUTORIO

V.

PROCEDIMIENTO DIDÁCTICO Las clases teóricas se desarrollan en forma expositiva, y al mismo tiempo resolviendo ejercicios base. Con estos conocimientos se resolverán ejercicios de mayor complejidad con participación de los estudiantes. El profesor de práctica se encargará de la resolución de problemas aplicados a la Ingeniería Química con la finalidad de que los estudiantes se familiaricen con los tópicos que se estudian en Ingeniería Química.

VI.

EQUIPOS Y MATERIALES Para el desarrollo adecuado del curso, se requiere lo siguiente:  Calculadora científica.  Computadora con un software de programación de alto nivel.  Tiza blanca y de color. 148

VII. SISTEMA DE EVALUACIÓN El promedio final, se obtendrá del siguiente modo:

Peso  Examen Parcial

1.0

 Examen Final

1.0

 Promedio de Prácticas Calificadas

1.0

PF 

EP  EF  PP  11 3

Al final se tomará un Examen de Recuperación que sustituirá a la nota mas baja del Examen Parcial o Final. El examen de recuperación se tomará en base al total del avance del curso.

149

VIII. BIBLIOGRAFÍA

1. Carnahan, B. Luther, A. Wilkes

Cálculo

Numérico,

Aplicaciones

Editorial Rueda, Madrid, 1979.

2. Burden, R. Y Faires J.

Análisis Numérico. Edit. Iberoamericana, México, 1985.

3. Nieves, A., Domínguez, F.

Métodos Numéricos Aplicados a la Ingeniería Química. Edit. CECSA, México 1985.

4. Valderrama, R.

Métodos Numéricos. Edit. Trillas, UNAM, México, 1985.

5. Nakamura, S.

Métodos Numéricos aplicados con software. Edit. Prentice – Hall Hispano Americano, S.A. México, 1992.

6. Gerald, C.

Análisis Numérico Edit. Alfa – Omega, México, 1991.

7. Maron, M. y López R.

Análisis

Numérico.

Un

práctico. Edit. CECSA, México, 1998.

150

enfoque

Metodos Numericos En Polymath 5n3z6y

Overview 4q3b3c

More details 26j3b

Related Documents 171j1w

Metodos Numericos En Polymath 5n3z6y

Metodos Numericos En Matlab 5e2749

Metodos Numericos En Javascript 4uy68

Metodos Numericos n5r11

Metodos Numericos n5r11

Metodos Numericos n5r11

More Documents from "Edgar De La Cruz" 4i3z4m

2c82t

Cableado Vertical 1qj49

2da Ley De Newton 39422f

Ejercicios De Muestreo De Aceptacion 2pa59

Prueba De Tukey k3k4u

02 Modelos Arma 1sw4