Punto Fijo Multivariable

APELLIDOS Y NOMBRES: CAJAHUAMAN RODRIGUEZ MARÍA EUGENIA ESCUELA: INGENIERÍA QUÍMICA PUNTO FIJO MULTIVARIABLE EJERCICIO

Views 96 Downloads 4 File size 680KB

Report DMCA / Copyright

DOWNLOAD FILE

Recommend stories

Citation preview

APELLIDOS Y NOMBRES: CAJAHUAMAN RODRIGUEZ MARÍA EUGENIA ESCUELA: INGENIERÍA QUÍMICA

PUNTO FIJO MULTIVARIABLE

EJERCICIO: DESARROLLAR EL SIGUIENTE SISTEMA DE ECUACIONES POR LOS DIFERENTES MÉTODOS ESTUDIADOS. 𝒇𝟏(𝒙, 𝒚) = 𝒚𝟐 𝒙 − 𝟖𝒙 + 𝟐𝒚 + 𝟓 𝒇𝟏(𝒙, 𝒚) = 𝒙𝟐 𝒚 + 𝒙 − 𝟏𝟎𝒚 + 𝟑 i

x 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16

y 0 0.625 0.70703125 0.73093124 0.73620734 0.73775276 0.73812881 0.73823387 0.73826047 0.7382677 0.73826957 0.73827007 0.7382702 0.73827023 0.73827024 0.73827025 0.73827025

d 0 0.3 0.37421875 0.38941007 0.39389776 0.39497004 0.39527267 0.39534869 0.39536947 0.39537481 0.39537625 0.39537662 0.39537672 0.39537675 0.39537675 0.39537676 0.39537676

0.69327123 0.11062346 0.02831935 0.00692652 0.00188099 0.0004827 0.00012967 3.3755E-05 8.992E-06 2.3563E-06 6.2483E-07 1.6428E-07 4.3459E-08 1.1446E-08 3.0243E-09 7.9724E-10

i 0 1 2 3 4 5 6 7 8 9 10 11 12

x

y

d

0 0.625 0.72589111 0.73684334 0.73810483 0.73825106 0.73826802 0.73826999 0.73827022 0.73827024 0.73827025 0.73827025 0.73827025

0 0.3625 0.39168989 0.39495067 0.39532735 0.39537102 0.39537609 0.39537668 0.39537675 0.39537676 0.39537676 0.39537676 0.39537676

0.7225173 0.10502888 0.01142733 0.00131653 0.00015261 1.7704E-05 2.054E-06 2.3829E-07 2.7646E-08 3.2074E-09 3.7211E-10 4.3172E-11

MÉTODO DE NEWTON RAPHSON %metodo de newton raphson %teniendo en cuenta que Tte=x Tti=y Ts=z clc,clear %la primera ecuacion %8.7459*10^-9*x^4+9.5504*10^-6*x^3+3.9110*10^-3*x^2+4.2556*x-7.7359*10^9*y^4+8.4474*10^-6*y^3-3.4592*10^-3*y^2+0.62958*y-361.58 %la segunda ecuacion %1.0100*10^-9*x^4+1.1029*10^-6*x^3+4.5164*10^-4*x^2+3.626*x+9.896*y4.948*z-440.75 %la tercera ecuacion %1.0100*10^-9*x^4+1.1029*10^-6*x^3+4.5164*10^-4*x^2+3.626*x+3.5326*10^8*z^3-9.3743*10^-5*z^2+0.7002*z-372.76 syms x syms y syms z cf1=input('ingrese la primera ecuacion=','s'); f1=inline(cf1); cf2=input('ingrese la segunda ecuacion=','s'); f2=inline(cf2); cf3=input('ingrese la tercera ecuacion=','s'); f3=inline(cf3); dcf1=diff(cf1,x,1); a=inline(dcf1); dcf2=diff(cf1,y,1); b=inline(dcf2); dcf3=diff(cf1,z,1); %constante c=inline(dcf3); dcf4=diff(cf2,x,1); d=inline(dcf4); dcf5=diff(cf2,y,1); %constante e=inline(dcf5); dcf6=diff(cf2,z,1); %constante f=inline(dcf6); dcf7=diff(cf3,x,1); g=inline(dcf7); dcf8=diff(cf3,y,1); %constante h=inline(dcf8); dcf9=diff(cf3,z,1); j=inline(dcf9); x=input('ingrese el valor de x='); y=input('ingrese el valor de y='); z=input('ingrese el valor de z='); dist=1; k=0; disp('---------------------------------------------------------'); disp('------------------TABLA DE RESPUESTAS--------------------'); disp('---------------------------------------------------------'); while dist>0.0001 m=[a(x) b(y) dcf3 -f1(x,y);d(x) dcf5 dcf6 -f2(x,y,z);g(x) dcf8 j(z) f3(x,z)]; n=rref(m); r=n(1,4); s=n(2,4); t=n(3,4); disp('---------------------------------------------------------');

k=k+1 x1=x+r y1=y+s z1=z+t dist=[(x1-x)^2+(y1-y)^2+(z1-z)^2]^0.5 x=x1; y=y1; z=z1; end

Ejecución del programa Ingrese la primera ecuacion=8.7459*10^-9*x^4+9.5504*10^-6*x^3+3.9110*10^3*x^2+4.2556*x-7.7359*10^-9*y^4+8.4474*10^-6*y^3-3.4592*10^-3*y^2+0.62958*y361.58 Ingrese la segunda ecuacion=1.0100*10^-9*x^4+1.1029*10^-6*x^3+4.5164*10^4*x^2+3.626*x+9.896*y-4.948*z-440.75 Ingrese la tercera ecuacion=1.0100*10^-9*x^4+1.1029*10^-6*x^3+4.5164*10^4*x^2+3.626*x+3.5326*10^-8*z^3-9.3743*10^-5*z^2+0.7002*z-372.76 Ingrese el valor de x=0 Ingrese el valor de y=0 Ingrese el valor de z=0 -----------------------------------------------------TABLA DE RESPUESTAS-----------------------------------k=

1

x1 =69.299948940698789243369865307344 y1 =105.89144713612604036963563200557 z1 =173.49098134965179977655079748011 dist =214.74310724545191712126119981639 --------------------------------------------------------k=

2

x1 =70.99038411610885020723037789018 y1 =100.38878772585653283008209830226 z1 =164.26897183442935993090751671722 dist =10.871167884203813915818305471981 ---------------------------------------------------------

k=

3

x1 =70.995306342281621794200821337734 y1 =100.37643500001544245418854084841 z1 =164.2483667228208059206431053477 dist =0.024523229204206096797477430450494 --------------------------------------------------------k=

4

Las raíces de los valores:

x1 =70.995306355501084283436758469309

Tte=X=70.995306355501

y1 =100.37643496802969123484664797628

Tti=Y=100.376434968029

z1 =164.24836667226776621915030426295

Ts=Z=164.248366672267

dist=0.000000061265424936342804582344332710 --------------------------------------------------------

MÉTODO DE NEWTON RAPHSON MODIFICADO %metodo de newton raphson modificado %teniendo en cuenta que Tte=x Tti=y Ts=z %la primera ecuacion %8.7459*10^-9*x^4+9.5504*10^-6*x^3+3.9110*10^-3*x^2+4.2556*x-7.7359*10^9*y^4+8.4474*10^-6*y^3-3.4592*10^-3*y^2+0.62958*y-361.58 %la segunda ecuacion %1.0100*10^-9*x^4+1.1029*10^-6*x^3+4.5164*10^-4*x^2+3.626*x+9.896*y4.948*z-440.75 %la tercera ecuacion %1.0100*10^-9*x^4+1.1029*10^-6*x^3+4.5164*10^-4*x^2+3.626*x+3.5326*10^8*z^3-9.3743*10^-5*z^2+0.7002*z-372.76 clc,clear syms x syms y syms z cf1=input('ingrese la primera ecuacion=','s'); f1=inline(cf1); cf2=input('ingrese la segunda ecuacion=','s'); f2=inline(cf2); cf3=input('ingrese la tercera ecuacion=','s'); f3=inline(cf3); cdf1=diff(cf1,x,1); a=inline(cdf1);

cdf2=diff(cf2,y,1); b=inline(cdf2); cdf3=diff(cf3,z,1); c=inline(cdf3); x=input('ingrese el valor inicial de x='); y=input('ingrese el valor inicial de y='); z=input('ingrese el valor inicial de z='); n=2; d=1; fprintf('n xk yk zk d \n') fprintf('%i %6.6f %6.6f %6.6f ------ \n',n-1,x,y,z) while d>0.0001 x1=x-f1(x,y)/a(x); y1=y-f2(x,y,z)/b(y); z1=z-f3(x,z)/c(z); d=[(x1-x)^2+(y1-y)^2+(z1-z)^2]^0.5; fprintf('%i %6.6f %6.6f %6.6f %6.6f \n',n,x1,y1,z1,d) x=x1; y=y1; z=z1; n=n+1; end

Ejecución del problema ingrese la primera ecuacion=8.7459*10^-9*x^4+9.5504*10^6*x^3+3.9110*10^-3*x^2+4.2556*x-7.7359*10^-9*y^4+8.4474*10^6*y^3-3.4592*10^-3*y^2+0.62958*y-361.58 ingrese la segunda ecuacion=1.0100*10^-9*x^4+1.1029*10^6*x^3+4.5164*10^-4*x^2+3.626*x+9.896*y-4.948*z-440.75 ingrese la tercera ecuacion=1.0100*10^-9*x^4+1.1029*10^6*x^3+4.5164*10^-4*x^2+3.626*x+3.5326*10^-8*z^3-9.3743*10^5*z^2+0.7002*z-372.76 ingrese el valor inicial de x=0 ingrese el valor inicial de y=0 ingrese el valor inicial de z=0 n xk yk zk d 1 0.000000 0.000000 0.000000 ---------2 84.965692 44.538197 532.362182 540.936515 3 74.002805 279.183799 71.026512 517.696189 4 69.634305 52.637871 146.966734 238.975016 5 73.294539 92.245482 171.689745 46.833614 6 71.279841 103.235166 151.532379 23.046726 7 70.906125 93.911814 162.664805 14.525683 8 71.214426 99.618072 164.739826 6.079650 9 71.019812 100.540052 163.039476 1.943994 10 70.990072 99.762807 164.113061 1.325737 11 71.015073 100.310743 164.277234 0.572549 12 70.997413 100.383462 164.139332 0.156898 13 70.995081 100.321128 164.236747 0.115674 14 70.997079 100.370709 164.249608 0.051261

15 16 17 18 19 20 21

70.995490 70.995308 70.995465 70.995323 70.995308 70.995321 70.995308

100.376392 100.371476 100.375928 100.376372 100.375991 100.376389 100.376424

164.238587 164.247355 164.248359 164.247490 164.248277 164.248355 164.248288

0.012501 0.010053 0.004567 0.000986 0.000874 0.000407 0.000077

Las raíces de los valores: Tte=X=70.995308 Tti=Y=100.376424 Ts=Z=164.248288

MÉTODO DE BROYDEN %INTEGRACION POR EL METODO DE BROYDEN clc,clear all disp(' Método de integración de BROYDEN ') syms x syms y for i=1:5 x0=0; y0=0; f1=4*x0+2*y0+x0^2-6; f2=x0*y0^2+x0-10*y0-3; n=[x0; y0]; J=[2*x0+4, 2;y0^2+1, 2*x0*y0-10]; f0=[4*x0+2*y0+x0^2-6; x0*y0^2+x0-10*y0-3]; IJ=inv(J); x=(n-(IJ*f0)); dist=abs(sqrt(n-x)); disp(':::::::::::::::::::::::::::::::::::::::') disp([ x dist ]) n=x; end

1.5714 1.2536 -0.1429 0.3780

1.2341 0.5808 -0.1727 0.1727

1.2171 0.1305

-0.1746 0.0434

1.2170 0.0093 -0.1746 0.0033

1.2170 0.0037 -0.1746 0.0033

POLINOMIO DE LAGRANGE Elaboración del programa en Matlab para la primera pregunta:

function polinomio_de_Lagrange clc,clear % programa para el uso de cinco varialbes % t0=input('ingrese el valor de t0 = '); t1=input('ingrese el valor de t1 = '); t2=input('ingrese el valor de t2 = '); t3=input('ingrese el valor de t3 = '); t4=input('ingrese el valor de t4 = '); c0=input('ingrese el valor de c0 = '); c1=input('ingrese el valor de c1 = '); c2=input('ingrese el valor de c2 = '); c3=input('ingrese el valor de c3 = '); c4=input('ingrese el valor de c4 = '); syms t l0=((t-t1)*(t-t2)*(t-t3)*(t-t4))/((t0-t1)*(t0-t2)*(t0-t3)*(t0-t4)); l1=((t-t0)*(t-t2)*(t-t3)*(t-t4))/((t1-t0)*(t1-t2)*(t1-t3)*(t1-t4)); l2=((t-t0)*(t-t1)*(t-t3)*(t-t4))/((t2-t0)*(t2-t1)*(t2-t3)*(t2-t4)); l3=((t-t0)*(t-t1)*(t-t2)*(t-t4))/((t3-t0)*(t3-t1)*(t3-t2)*(t3-t4)); l4=((t-t0)*(t-t1)*(t-t2)*(t-t3))/((t4-t0)*(t4-t1)*(t4-t2)*(t4-t3)); y=l0*c0+l1*c1+l2*c2+l3*c3+l4*c4; disp('..................................................................' ) disp('........... INTERPOLACIÓN CON EL POLINOMMO DE LAGRANGE............') disp('..................................................................' ) %la ecución es de cuarto orden% %donde t: es el tiempo en segundos% %donde c: esla velocidad de cambio de concentración% ecuacion=simplify(y) t=input('ingrese el valor de t = '); disp('-----------------------------------------------------------------') c=eval(y)

CORRIDA DEL PROGRAMA PARA LA PRIMERA PREGUNTA: ingrese el valor de t0 = 40 ingrese el valor de t1 = 60 ingrese el valor de t2 = 120 ingrese el valor de t3 = 180 ingrese el valor de t4 = 300 ingrese el valor de c0 = 6 ingrese el valor de c1 = 5 ingrese el valor de c2 = 3 ingrese el valor de c3 = 2 ingrese el valor de c4 = 1 ............................................................................................................. ........... INTERPOLACIÓN CON EL POLINOMMO DE LAGRANGE............ .............................................................................................................

ecuacion =

t^4/2358720000 - (157*t^3)/235872000 + (223*t^2)/655200 - (5189*t)/65520 + 1577/182

ingrese el valor de t = 240 ------------------------------------------------------------------

c=

1.4670

DIFERENCIAS DIVIDIDAS function program_diffdivididas clc clear x=[56.5 78.6 113.0 144.5 181.0 205.0 214.5]; y=[1 2 5 10 20 30 40]; xa=x;ya=y; d=zeros(length(y)); d(:,1)=y'; for k=2:length(x) for j=1:length(x)+1-k d(j,k)=(d(j+1,k-1)-d(j,k-1))/(x(j+k-1)-x(j)); end end for w=1:length(x) ds=num2str(abs(d(1,w))); if w>1 if x(w-1)

DIFERENCIAS DIVIDIDAS PARA ADELANTE POLINOMIOS DE NEWTON EN DIFRENCIAS FINITAS puntos T(°C) P(atm)

0 50 40.94

1 60 49.11

2 70 56.05

3 80 62.84

4 90 68.57

5 100 72.3

puntos

x 0

fx 50

1er

2do

3ro

4to

5to

40.94 8.17

1

60

49.11

-1.23 6.94

2

70

56.05

3

80

62.84

1.08 -0.15

6.79

-0.91 -1.06

5.73 4

90

-1.99 1.96 -0.03 -0.94

68.57

-2 3.73

5

100

72.3

hacia adelante h x xo s

10 92 50 4.2

p1 p2 p3 p4 p5

75.254 66.9884 72.31064 71.79904 71.9149466

DIFERENCIAS FINITAS PARA ATRÁS %metodo de diferencias finitas hacia atras %para obtener una ecuacion cuadratica y el valor de Cp clc,clear %la temperatura en x esta en grados celsiuss x=[-17.78 93.33 204.44 3155.56 426.67 537.78 648.89 760 871.11 982.22 1093.33]; fx=[0.2498 0.2529 0.2552 0.2595 0.2628 0.2669 0.27 0.2718 0.2765 0.28 0.285]; n=11; h=111.11; for i=1:n-1 t(i,1)=fx(i+1)-fx(i); end for j=2:n-1 for i=j:n-1 t(i,j)=t(i,j-1)-t(i-1,j-1); end end % el valor de xint se debe ser 980 para comprobar la ecuacion formada xint=input('ingrese el valor de xint='); s=(xint-x(11))/h; Cp=fx(11)+s*t(10,1)+s*(s+1)/2*t(9,2); fprintf('t(%6.6f)=%6.6f \n',xint,Cp) disp('------------------------------------------------------------------------------')

disp('------------------------------------------------------------------------------') disp('---------------------LA ECUACION QUE SE HA DE FORMAR ES-----------------------') disp(' Cp=0.285+s*0.005+s*(s+1)/2*0.0015.........unidades en btu/lb °C') disp(' Cp=6.0751*10^-8*xint^2-8.1092*10^-5*xint+0.30104.....unidades en btu/lb °c') m=6.0751*10^-8*xint^2-8.1092*10^-5*xint+0.30104 disp('----------------------------LA ECUACION CORRECTA ES---------------------------') disp('Cp=1.413*10^-7*xint^2-1.8862*10^4*xint+0.70021..............unidades en Ws/g°C') %SABIENDO QUE EL VALOR DE Ts=z disp('---------------------------sabiendo que xint=Tinf2=(16+z)/2-----------------') disp(' Cp=1.413*10^-7*xint^2-1.8862*10^-4*xint+0.70021') disp('----------LA ECUACION QUE SE HA DE USAR EN LOS DEMAS PROGRAMAS----------------') disp(' Cp=1.413*10^-7*((16+z)/2)^2-1.8862*10^4*(16+z)/2+0.70021')

Ejecucion del programa t = t(10,1)=0.0050 0.0031 0 0.0023 -0.0008 t(9,2)=0.0015 0.0043 0.0020 0.0033 -0.0010 xint=T∞2 0.0041 0.0008 0.0031 -0.0010 0.0018 -0.0013 0.0047 0.0029 0.0035 -0.0012 0.0050 0.0015 ingrese el valor de xint=980 t(980.000000)=0.279888 -------------------------------------------------------------------------------------------------------------------------------------------------------------------------------LA ECUACION QUE SE HA DE FORMAR ES----------------------Cp=0.285+s*0.005+s*(s+1)/2*0.0015.........unidades en btu/lb °C Cp=6.0751*10^-8*xint^2-8.1092*10^-5*xint+0.30104.....unidades en btu/lb °c m = 0.2799 ----------------------------LA ECUACION CORRECTA ES--------------------------Cp=1.413*10^-7*xint^2-1.8862*10^-4*xint+0.70021..............unidades en Ws/g°C

---------------------------sabiendo que xint=Tinf2=(16+z)/2----------------Cp=1.413*10^-7*xint^2-1.8862*10^-4*xint+0.70021 ----------LA ECUACION QUE SE HA DE USAR EN LOS DEMAS PROGRAMAS---------------Cp=1.413*10^-7*((16+z)/2)^2-1.8862*10^-4*(16+z)/2+0.70021

INTEGRACIÓN NUMERICA 5

∫ (cos 𝑥 + 𝑥 2 + 3)𝑑𝑥 2

n= a= b= h=

3 2 5 1

x0= x1= x2= x3=

simpson 3/8

2 3 4 5

TRAPEZOIDAL COMPUESTO

I= 46.1012276

n= a= b= h=

f(x) 6.58385316 11.0100075 18.3463564 28.2836622

4 2 5 0.75

simpson 3/8 I= 46.132755

I= 46.7901216

x0= x1= x2= x3= x4=

2 2.75 3.5 4.25 5

f(x) 6.58385316 9.63819762 14.3135433 20.6164125 28.2836622

SIMPSON COMPUESTO I= 46.1282606

n= a= b= h=

5 2 5 0.6

x0= x1= x2= x3= x4= x5=

simpson 3/8

f(x) 6.58385316 8.90311125 12.2417052 16.6490323 22.0526671 28.2836622

2 2.6 3.2 3.8 4.4 5

TRAPEZOIDAL COMPUESTO

I= 46.1323234

I= 46.3681641

MÉTODO SIMPSON 3/8 Hallando la ecuación donde K esta e función de T:

T 0 100 200 300 400 600 800 1000

K 73 67 62 55 48 40 36 35

80 70 60 50 40 30 20 10 0

y = 4E-05x2 - 0.079x + 74.304 R² = 0.9943

0

200

400

600

800

1000

1200

𝐾 = 4 ∗ 10−5 ∗ 𝑇 2 − 0.079 ∗ 𝑇 + 74.304 … … … … … … 𝑒𝑐𝑢𝑎𝑐𝑖ó𝑛 1

Del enunciado se tiene: 𝐴 = 2∗𝜋∗𝑟∗𝐿 ,

𝐿 = 2𝑚

𝑅1 = 20𝑐𝑚 , 𝑅2 = 120𝑐𝑚 (En metros) 𝑇1 = 900°𝐶 , 𝑇2 = 250°𝐶

𝑞 = −𝐾𝐴

𝑑𝑇 𝑑𝑟

… … … … . 𝑒𝑐𝑢𝑎𝑐𝑖ó𝑛 2

Reemplazando los valores:

𝑞 = −(4 ∗ 10−5 ∗ 𝑇2 − 0.079 ∗ 𝑇 + 74.304) ∗ (2 ∗ 𝜋 ∗ 𝑟 ∗ 2)

𝑑𝑇 𝑑𝑟

𝑅2 =120𝑐𝑚

∫𝑅

1 =20𝑐𝑚

𝑞

𝑑𝑟 𝑟

q[ln(𝑟)]1.2 0.2 q=

𝑇 =250°𝐶

2 = ∫𝑇 =900°𝐶 −(4 ∗ 10−5 ∗ 𝑇2 − 0.079 ∗ 𝑇 + 74.304) ∗ (4 ∗ 𝜋)𝑑𝑇 1

= (4 ∗ 𝜋) ∫𝑇

𝑇2 =250°𝐶

1 =900°𝐶

(4∗𝜋)

−(4 ∗ 10−5 ∗ 𝑇 2 − 0.079 ∗ 𝑇 + 74.304)𝑑𝑇

𝑇 =250°𝐶

2 −(4 ∗ 10−5 ∗ 𝑇 2 − 0.079 ∗ 𝑇 + 74.304)𝑑𝑇 ∫ 𝑇 =900°𝐶 ln(1.2)−ln(0.2) 1

Se utiliza el algoritmo de SIMPSON 3/8. Elaboración del programa en Matlab para la

segunda pregunta: function metodo_simpson_tres_octavos clc, clear all syms x ax=input('ingrese el valor del limite inferior: bx=input('ingrese el valor del limite superior: n=input('ingrese el valor de n: '); hx=(bx-ax)/n; x0=input('ingrese el valor de x0: '); x1=x0+hx; x2=x1+hx; x3=bx; f0=input('ingrese la funcion a integrar f0(x)=' %f0=4*10^(-5)*x^2-0.079*x+74.304% x=x0; fx0=eval(f0) x=x1; fx1=eval(f0) x=x2; fx2=eval(f0) x=x3; fx3=eval(f0) disp('------------------------------------') disp('--------Hallando la integral------- ') disp('------------------------------------') I1=(3*hx/8)*(fx0+3*fx1+3*fx2+fx3) q=((-4*pi)/(1.791759469))*I1

'); ');

,'s');

CORRIDA DEL PROGRAMA PARA LA SEGUNDA PREGUNTA: ingrese el valor del limite inferior: 900 ingrese el valor del limite superior: 250 ingrese el valor de n: 3 ingrese el valor de x0: 900 ingrese la funcion a integrar f0(x)=4*10^(-5)*x^2-0.079*x+74.304

fx0 =

35.6040

fx1 =

38.9984

fx2 = 46.1484 fx3 = 57.0540 -------------------------------------------Hallando la integral-----------------------------------------I1 =

-2.8283e+04

q = 1.9836e+05 (=198361.2794)