METODOS NUMERICOS

PRESENTADO POR:
Angulo Heyner
Lambis Tiffany
Luna Maria

PRESENTADO A:

Dina Mendoza

Facultad de Ingenieria

Ignenieria Quimica

UNIVERSIDAD DE CARTAGENA

CARTAGENA DE INDIAS D.T.Y.C

=2018=
 TALLER DE METODOS NUMERICOS

1. Realizar aproximaciones por el método de mínimos cuadrados, a

partir de los datos experimentales de P, X1, Y1, datos de equilibrio
1-propanol (1) Agua (2) a 45°C, utilizando los modelos de Margules y
el de Van Laar. Escoger el que se ajusta con mejor exactitud a los
datos experimentales.

2. Del ejemplo 10.1 que se encuentran en la parte final de las

diapositivas de métodos numérico; deben primeramente entender
el ejemplo, para que puedan llenar la matriz paso a paso y obtener
la matriz tridiagonal resultante mostrada en la diapositiva, y
generen un informe detallado en Word de cómo se llega a cada
elemento de la matriz, partiendo de los datos de especificaciones
de la torre y propiedades de los componentes procesados.
Finalmente resuelva la matriz tridiagonal usando el código de
Matlab que se encuentra también en las diapositivas enviadas;
modificándolo para que imprima los cálculos intermedios que
ustedes consideren necesarios, además de los nombres de los
integrantes del grupo.
DESARROLLO:

1.

Para el desarrollo de nuestro ejercicio de una mezcla binaria

inicialmente se ha calculado los valores de X 2 y Y2, para ello fue
necesario el uso de las siguientes formulas:

X 2 =1−X 1
Y 2=1−Y 1

Ahora se procese a calcular los coeficientes de actividad, Los cuales

fueron obtenidos a partir de la ecuación de Raoult modificada, como se
muestra a continuación.

Y i P=X i γ i Pi sat
YiP
=γ i
X i P isat

Para hallar las presiones de saturación, se es necesario usar la ecuación

de Antoine, está dada en (mmHg):

B
ln Psat =A−
T ( k )+ c

Sustancia\constantes A B C
1-Propanol 17,54 3166,38 -80,15
Agua 18,30 3816,44 -46,13

Remplazando valores nos queda que las presiones de saturación son:

Psat (1−propanol)= 69, 39 mmHg

Psat ( Agua)= 71, 78 mmHg

Cambiando las unidades:

Psat (1−propanol)= 9,25 Kpa

Psat ( Agua)= 9,57 Kpa

El siguiente paso a realizar es calcular la energía libre de Gibbs, cuando

la solución es binaria, esto se realiza con la siguiente ecuación:

GE
=x 1 ln γ 1 +¿ x2 ln γ 2 ¿
RT
Se realizan los cálculos en Excel Y nos queda que:

Método de Margules

El Método de Margules nos dice que:

GE x1
= A 21 + A12
RT x 1 x22 x2

Analizando la ecuación nos podemos dar cuenta que tiene el mismo

comportamiento de la ecuación de la recta Y=mx+b, Gracias a esto, se
es posible hallar los parámetros, con los datos experimentales hallados
anteriormente, nos queda que:
Al momento de graficar, nos queda que:

Se implementa el método de mínimos cuadrados:

Fo=∑ (Yiexp−Ycal)2 =0

GE
Si se sabe que nuestro Yexp= y que nuestro Ycal= mx +b,
RT x 1 x22
además que los parámetros son igual a: b= A12 y m= A21, se deriva la
ecuación con respecto a cada parámetro:

Dfo
=−2 ∑ (Yiexp−mx−b)2=0
Db
Dfo
=−2 ∑ (Yiexp−mx−b)2=0
Dm
Si se reacomoda las ecuaciones, nos queda que:

∑ Yiexp=m ∑ xi+b∗n
∑ Yiexp xi=m ∑ xi2+b ∑ xi
Remplazando de los datos nos quedó que:

Metodo de Margules
12

10
f(x) = 1.3 x + 2.13
R² = 0.99
8
G^E/RTX1X2^2

0
0 1 2 3 4 5 6 7
x1/X2
 METODO DE VAN LAAR

Este metodo se encuentra dado por:

RTx 1 x 2 A 12− A 21 1
= x 1+
G E
A 12∗A 21 A 12

La ecuacion anterior se asemeja como la de una recta y podemos hacer

que:

A 12−A 21 1
Ycal=mx +b, donde x= x1; m= ; b=
A 12∗A 21 A 12

Remplasando los valores nos queda que:

Van Laar
0.8
f(x) = 0.34 x + 0.43
0.7 R² = 0.97

0.6

0.5
G^E/RTX1X2

0.4

0.3

0.2

0.1

0
0 0.2 0.4 0.6 0.8 1 1.2
RTX1X2/GE
Comparacion de ambos metodos

  Margules Van Laar

Método 1 2 1 2

12 2,1304 2,2725 2,3337 2,3339  

21 1,3002 1,1631 1,2977 1,2977  

R^2 0,9941 - 0,9671 -

Comparación Teórico vs Exp

Exp 1 2 1 2

GE*/RTX1X GE/RTX1X2
2

2,3727 2,0920 2,2213 2,2507 2,2508

2,2029 2,0509 2,1664 2,1681 2,1682

2,0338 1,9850 2,0783 2,0475 2,0476

1,8600 1,8967 1,9602 1,9055 1,9056

1,6219 1,7337 1,7425 1,6893 1,6894

1,5270 1,6285 1,6018 1,5740 1,5740

1,4444 1,4916 1,4189 1,4457 1,4457

1,4462 1,4173 1,3197 1,3844 1,3844

Grafica resultante de la comparación:

Consolidado G^E/RTX1X2
vs X1
2.5
G^E/RTX1X2

1
0 0.2 0.4 0.6 0.8 1
X1

G^E/RTX1X2 Margules G^E/RTX1X2 Van Laar G^E/RTX1X2 Experimental

 2.

Para la columna de destilación que se muestra en la figura 10.7, realice

una iteración del método de BP hasta incluir el cálculo de un nuevo
conjunto de valores de T, a partir de (10-20). Asumir los valores k
independiente de la composición.

Se conoce que F3=100lbmol /h, L5=50 lbmol /h y que L1/U1= 2

Balance General

Liquido destilado =U 1=F 3−L5 = 100-50

U 1=50 lbmol /h
Como L1/U1= 2, entonces

L1 ( 2 )∗50 lbmol 100 lbmol

L 1= ∗U 1 → L 1= =
U1 h h
Balance total de la materia alrededor del condensador es:

V 2=L 1+ U 1
V 2=100 lbmol /h+50 lbmol /h=150 lbmol /h

Las variables iniciales son:

A 100 psia, los valores estimados de k para la temperatura asumida en
cada etapa son:

La ecuación para obtener los valores de la matriz tridiagonal utilizando

los valores del componente C3, con V1=0 y por tanto W=0 es
j−1
A j=V j+ ∑ (F m−U m)
m=1

Ósea,

A5 =V 5+( F 3−U 1 )=150+ ( 100−50 )=200

Lo mismo con las siguientes,

A 4=200

A3 =100

A2=100

Para los valores de B,

j
B j =−(V j+1 + ∑ ( F m−U m ) + U j +V j k ij )=−( F3 −U 1+ V 5 k 1,5)
m= j

B5=−( 100−50+150 x 3.33 ) =−549.5

Mismo procedimiento en las siguientes,

B4 =−605

B3=−525.5

B2=−344.5

B1=−150

Para los valores de C,

C j=V j+1 K 1 , j+1

C 1=V 2+1 K 1,2 150 x 1.63=244.5

Mismo procedimiento en las siguientes,

C 2=325.5

C 3=405

C 4=499.5

Para los valores de D,

D j=−FjZ 1 , j

Así,

D3=−100 x 0.3=−30

D1=D2=D4=D5=0

Organizando la matriz tridiagonal según la forma matricial del algoritmo

de Thomas

−150 244.5 0 0 0 ¿ 0

[ 100 −344.5
0
0
0
100
0
0
200
0
325.5 0 0 ¿
−525.5 405

Resuelto en MatLab 2017

0
−605 499.5
200 −549.5
0
¿ −30
¿
¿
0
0
]
function x = TDMAsolver(a,b,c,d)
%a, b, c columnas, d vector derecho
% N número de columnas
N= length(d);% Modifica coeficientes iniciales
c(1) = c(1) / b(1); % Division cero
d(1) = d(1) / b(1);
for n = 2:1:N
temp = b(n) - a(n) * c(n - 1);
if (n<N)
c(n) = c(n) / temp;
end
d(n) = (d(n) - a(n) * d(n - 1)) / temp;
end
% sustitución atrás.
x(N) = d(N);
for n = (N - 1):-1:1
x(n) = d(n) - c(n) * x(n + 1);
end
end

>> clear
>> a=[0 200 200 100 100]; b=[-150 -344.5 -525.5 -605 -549.5]; c=[244.5
325.5 405 499.5 0]; d=[0 0 -30 0 0 ]
d=
0 0 -30 0 0
>> TDMAsolver(a, b, c, d)

ans =

-0.2801 -0.1718 -0.0098 -0.0019 -0.0003

3. Estime el valor de y(0,5) con h=0,25 mediante el método de Euler,

Heun y RK4.

y ‘=5 x / y – xy , y (0)=2 , y(0,5)=?

Métodos

Por separación de variables

dy 5 x
= −xy
dx y

dy 5− y 2
dx
=x ( )
y

( 5−yy ) dy =xdx
2

y (0,5 ) 0,5
y
∫
2
( 5− y 2 )
dy=∫ xdx
0

−1 2 2 0,52
2
( ln ( 5− y 0,5 ) −ln (5− 2 ))= 2
( ) ( )

( ln ( 5− y ( 0,5 )2) −ln(1))=−0,52

2
2
5− y ( 0,5 ) =e−0,5
2

y ( 0.5 )= √5−e−0,5
y=2.0545557

Método Euler
El método de Euler es el más simple de los métodos numéricos resolver
un problema del siguiente tipo:

Consiste en multiplicar los intervalos que va de   a   en   

subintervalos de ancho  ; osea:

de manera que se obtiene un conjunto discreto de   

puntos:   del intervalo de interes  . Para
cualquiera de estos puntos se cumlple que:

  .
 La condición inicial  , representa el
punto   por donde pasa la curva solución de la
ecuación de el planteamiento inicial, la cual se denotará
como  .

Ya teniendo el punto   se puede evaluar la primera derivada

de   en ese punto; por lo tanto:

Con esta información se traza una recta, aquella que pasa por   y de
pendiente . Esta recta aproxima   en una vecinidad de .
Tómese la recta como reemplazo de   y localícese en ella (la recta)
el valor de y correspondiente a  . Entonces, podemos deducir segun la
Gráfica A:

Se resuelve para  :

Es evidente que la ordenada   calculada de esta manera no es

igual a  , pues existe un pequeño error. Sin embargo, el
 valor   sirve para que se aproxime   en el
punto   y repetir el procedimiento anterior a fin de
generar la sucesión de aproximaciones siguiente:

Solución del punto 3 por el método de Euler empleando Matlab

 |2,03125−2.0545557|
Error= ∗100 %=1,13434 %
2.0545557
Solución del punto 3 con el método de Euler empleando Excel

Gráfica y vs y analítica en el método de Euler

Método de Heun - Método Euler Mejorado

El método de Euler mejorado es una modificación del método de

Euler para resolver EDO´s con condiciones iniciales. La solución que
ofrece este método, es una tabla de la función solución, con valores de
“y” correspondientes a valores específicos de “x”.

Es por esto que uno de los requisitos para este método es especificar el
intervalo de x.

También se requiere de:

- Una ecuación diferencial de primer orden.

Y’= f(x,y)

-         La condición inicial, es decir, el valor de y en un punto

conocido x0.

Y(x0) = y0

El método consiste en usar la ecuación de Euler como una

ecuación predictora y usar este resultado en la ecuación correctora de
Euler-Gauss.

Las ecuaciones del método de Euler Mejorado son las siguientes:

                        

Solución del punto 3 por el método de Euler modificado por medio de

Excel.
Solución del punto 3 por el método de Euler modificado por medio de
Matlab
 |2,0551983−2.0545557|
Error= ∗100 %=0,03128 %
2.0545557
Método RK4 -Métodos de Runge Kutta

Los métodos de Taylor tienen la propiedad de un error local de

truncamiento de orden superior, pero la desventaja de requerir el cálculo
y la evaluación de las derivadas de f(t, y). Esto resulta algo lento y
complicado, en la mayoría de los problemas, razón por la cual, en la
práctica casi no se utilizan. El método de Euler, lamentablemente
requiere de un paso muy pequeño para una precisión razonable.

Los métodos de Runge kutta tienen el error local de truncamiento del

mismo orden que los métodos de Taylor, pero prescinden del cálculo y
evaluación de las derivadas de la función f(t, y).

Se presenta de nuevo el problema de valor inicial cuya solución se

intenta aproximar:
 

(1)

Como en los métodos anteriores, se determina primero la malla {t 0, t1, ... ,

tN} de paso h, donde t0 = a y  tN = b. En estos puntos es donde se va a
obtener la aproximación de la solución.

En esencia, los métodos de Runge-Kutta son generalizaciones de la

fórmula básica de Euler yi+1 = yi + h f(ti, yi) en los que el valor de la
función f se reemplaza por un promedio ponderado de valores de f en el
intervalo ti ≤ t ≤ ti+1, es decir,
 (2)

En esta expresión las ponderaciones wi, i = 1, ..., m son constantes para

las que en general se pide que su suma sea igual a 1, es decir, w 1 +
w2 + ... + wm = 1, y cada kj es la función f evaluada en un punto
seleccionado (t, y) para el cual ti ≤ t ≤ ti+1. Se mostrará que los kj se
definen en forma recursiva.

Se define como orden del método al número m, es decir, a la cantidad de

términos que se usan en el promedio ponderado.

Solución del tercer punto usando el método de RK4 empleando Excel

Grafica Yi vs Solucion exacta

Solución del 3 punto usando el método de RK4 empleando Matlab
 |2,0545580−2,0545557|
Error= ∗100 %=0,00011%
2,0545557