MÉTODO DE HEUN INTRODUCCION Aplica Aplicarr método métodos s numéri numéricos cos para para aprox aproxima imarr soluci solucione ones s de alguna algunas s ecuaciones diferenciales, viendo así la importancia de los métodos numéricos que radica en la aparición de ecuaciones diferenciales que no pued pueden en reso esolver lvers se por mét método odos trad tradic icio iona nale les s, y de ahí ahí la necesidad de implementar algún método de aproximación.
DEFINICION Un método para mejorar la estimación de la pendiente emplea la deter determin minaci ación ón de dos derivada derivadas s en el interv intervalo alo una una en el punto punto inicial y otra en el !nal". #as dos derivadas se promedian después con la !nalidad de o$tener una mejor estimación de la pendiente en todo el intervalo. %ste procedimiento, conocido como método de &eun, se presenta en forma gr'!ca en la siguiente !gura
Primera iteración del método Euler modifcado (ecordar (ecordaremos emos que, en el método de %uler, la pendiente al inicio. yi′ = f ( xi , yi )
...( 1)
se utili)a para extrapolar linealmente a yi +1 yi +1
= y0 + f ( xi , yi ) h
..... ( 2 )
%n el método est'ndar de %uler de$ería parar aquí. *in em$argo, en el métod método o de &eun &eun la yi +1 calculada en la ecuación +" no es la
respuesta !nal, sino una predicción intermedia. or consiguiente, la distinguimos con un superíndice -. #a ecuación +" se llama ecuación predictora o simplemente predictor. a una estimación de yi +1 que permite el c'lculo de una estimación de la pendiente al !nal del intervalo/ yi′+1
=
(
)
f xi +1 , yi +1 .... ( 3 )
Así, se com$inan las dos pendientes 0ecuaciones 1" y 2"3 para o$tener una pendiente promedio en el intervalo/
yî′ + yi′+1 f ( xi , yi ) + f ( xi +1 , yi +1 ) ′ y = = 2
2
%st' pendiente promedio se utili)a después para linealmente desde yi hasta yi +1 con el método de %uler/ yi +1
= yi +
f ( xi , yi )
+ f ( xi +1 , yi +1 ) 2
extrapolar
h .... ( 4 )
A esta ecuación 4" se conoce como ecuación correctora o simplemente corrector %l método de &eun es el único método predictor5corrector de un solo paso.
redictor
→
6orrector
→
yi +1
= y0 + f ( xi , yi ) h
yi +1
= yi +
f ( xi , yi )
+ f ( xi +1 , yi +1 ) 2
h
7$serve que de$ido a que en la ecuación 4" aparece yi +1 a am$os lados del signo igual, entonces se puede aplicar en una forma iterativa. %s decir, una estimación anterior se utili)ar' de manera repetida para proporcionar una estimación mejorada de yi +1
ERROR RE!TI"O PORCENTU! ε t
=
yi j+1 − yij+−11 j i +1
y
E#ERCICIO$ RE$UETO$ E#ERCICIO %
×100%
yi j+−11
= iteracion anterior
yi j+1
= actual del corrector
6on el método de &eun integre y ′ = 4e 0.8 x
x = 0 hasta x = 4 con un tama8o de paso igual a 1. #a condición inicial es en x = 0 y y = 2
→
redictor
yi0+1
→
6orrector
yi +1
− 0.5 y
= yi + f ( xi , yi ) h = yi +
f ( xi , yi )
+ f ( xi +1 , yi0+1 ) 2
h
Antes de resolver el pro$lema numéricamente, se utili)a el c'lculo para determinar la siguiente solución analítica/ y = y = y = y = y = y =
4 1.3 4 1.3 4 1.3 4 1.3 4 1.3 4 1.3
(e
0.8 x
− e −0.5x ) + 2e −0.5x
(e 0.8(0) − e −0.5(0) ) + 2e −0.5(0) (e 0.8(1) − e −0.5(1) ) + 2e −0.5(1)
=2 = 6.1946314
(e 0.8(2) − e −0.5(2) ) + 2e −0.5(2)
= 14.8439219
(e 0.8(3) − e −0.5(3) ) + 2e −0.5(3)
= 33.6771718
(e
0.8(4)
− e−0.5(4) ) + 2e −0.5(4) = 75.3389626
ara la siguiente solución primero calcularemos el predictor yi0+1
y
seguidamente llevaremos este resultado a la ecuación del corrector yi +1 que vendría a ser la solución de &eun i = 0
= 0 ; x1 = 1 ; y0 = 2 ; h = 1 0 yi +1 = yi + f ( xi , yi ) h y10 = 2 + (4e 0,8(0) − 0.5(2))1 0 y1 = 5 f ( xi , yi ) + f ( xi +1 , yi0+1 ) yi +1 = yi + h
x0
2
4e0,8(0) − 0.5(2)! + 4e0,8(1) − 0.5(5)! y1 = 2 + ÷1 2 y1 = 6.7010819 6uando i = 1
= 1 ; x2 = 2 ; y1 = 6.7010819 ; h = 1 yi0+1 = yi + f ( xi , yi ) h y20 = 6, 7010819 + (4 e 0,8(1) − 0.5(6.7010819))1 0 y2 = 12.2527 f ( xi , yi ) + f ( xi +1 , yi0+1 ) yi +1 = yi + h
x1
2
4e0,8(1) − 0.5(6.7010819)! + 4 e0,8(2) − 0.5(12.2527)! y2 = 6.7010819 + ÷1 2 y2 = 16.3197819 6uando i = 2
= 2 ; x3 = 3 ; y2 = 16.3197819 ; h = 1 yi0+1 = y1 + f ( xi , yi ) h y30 = 16.3197819 + (4e 0,8(2) − 0.5(16.31 97819))1 y30 = 27.9720 f ( xi , yi ) + f ( xi +1 , yi0+1 ) yi +1 = yi + h
x2
2
4e0,8(2) − 0.5(16.3197819)! + 4e 0,8(3) − 0.5(27.9720)! y3 = 16.3197819 + ÷1 2 y3 = 37.1992489 6uando i = 3
= 3 ; x4 = 4 ; y3 = 37.1992489 ; h = 1 yi0+1 = yi + f ( xi , yi ) h y40 = 37.1992489 + (4e 0,8(3) − 0.5(37.1992489))1 y40 = 62.6923 f ( xi , yi ) + f ( xi +1 , yi0+1 ) yi +1 = yi + h
x3
2
4e0,8(3) − 0.5(37.19924 89)! + 4e 0,8(4) − 0.5(62.6923)! y4 = 37.1992489 + ÷1 2 y4 = 83.3377674
x 1
y verdadero +.-----;.1<4;21
y&eun +.-----;.=-1->1
9et9:" -.->.1>
+ 2 4
4 14.>42<+ 1< 22.;==1= 1> =?.22><; +;
< 1;.21<=> 1< 2=.1<<+4 >< >2.22==; =4
<.<4 1-.4; 1-.;+
E#ERCICIO & (esolver el pro$lema de valor inicial dy
=
dx
x− y
en 0,3!
2
y(0) = 1
6on n@2 6alculamos h/ h=
b−a n
=
3− 0 3
=1
rimera interacción/ &allamos el predictor yi′+1 = y0 + f ( xi , yi ) h y1
0 −1 = 1 + ÷ ∗1 = 0.5 2
&allamos el corrector yi +1
y1
= yi +
f ( xi , yi )
+ f ( xi +1, yi0+1 ) 2
h
1 ( 0 − 1) + ( 1 − 0.5 ) = 1+ ∗ ÷ ∗1 = 0.875 2 2
&allamos el predictor yi′+1 y2
= y0 + f ( xi , yi ) h
1 − 0.875 = 0.875 + ÷ ∗1 = 0.9375 2
&allamos el corrector
*egunda interacción/
yi +1 y2
= yi +
f ( xi , yi )
+ f ( xi +1 , yi0+1 ) 2
h
1 ( 1 − 0.875 ) + ( 2 − 0.9375 ) = 0.875 + ∗ ÷ ∗1 = 1.171875 2 2
ercera interacción/ &allamos el predictor yi′+1
= y0 + f ( xi , yi ) h 2 − 1.173325 ∗1 = 1.5859375 y3 = 1.173325 + ÷ 2 &allamos el corrector yi +1 y3
= yi +
f ( xi , yi )
+ f ( xi +1 , yi0+1 ) 2
h
1 ( 2 − 1.173325 ) + ( 3 −1.5859375 ) = 1.173325 + ∗ ÷ ∗1 = 1.73242188 2 2
i
Bi
-
1
1
+
+
2
yxi"exact Ci o -.>1=? > 1.1-2;2>2 1.1=1>=? + 1.=2+4+1 1.;;<2<-4 >> >
E#ERCICIO ' Un tanque cilíndrico de fondo plano con un di'metro de 1.?m contiene un líquido de densidad ρ = 1.5 Kg " L a una altura de 2m.*e desea sa$er la altura del líquido dentro del tanque tres minutos después de que se a$re completamente la v'lvula de salida, la cual da un gasto de 0.6 A 2 ga #3 " s donde A es el 'rea seccional del tu$o de salida y es 78.5 ×10 −4 # 2 y $9.81 #"s 2 con ; iteraciones.
$OUCION
G = 0.6 A 2 ga
Acumulación@entrada/
−0.6 A V
dV ρ dt
= 0
*alida/
2 ga
π
π
4
4
= ( 1.5) 2 a
( 1.5)
2
da dt
= −0.6 A
2 ga
e donde da dt
=−
2.4 A 2 ga π ( 1.5)
2
= −0.0026653
2 ga al considerar como tiempo cero el
momento de a$rir la v'lvula y adem's la altura $uscada a un tiempo de 1>da
= −0.0026653
dt s, se llega a a 0 ( )
2 ga
= 3m a ( 180 ) = &
&allamos por el método de &eun A. − yi +1 B. −
1 2
= yi + hf ( xi , yi ) → caso 'redictor
× f ( xi , yi ) + f ( xi +1 , yi +1 )
C. − yi +1
1 = yi + h × f ( xi , yi ) + f ( xi +1 , yi +1 ) → caso corrector 2
i ( % & ' ) * E#ERCICIO )
t i
A
B.
C
2;<1+1?1>-
+.2>;?? 1.>;>=; 1.412=21.-+124 -.;<1?>> -.4+44+<
5-.-1<242 5-.-1=+4< 5-.-1?1?5-.-12-;;? 5-.-1-<= 5-.-->>=?;
+.41<= 1.<-++1 1.44=4? 1.-??4? -.=+;+;
(+)*,,,
Una solución de salmuera de salmuera ra)ón constante de ;#Dmin, hacia el interior de un depósito que inicialmente contiene ?-# de solución de salmuera en la cual se disolvieron :Eg de sal. #a solución contenida en el depósito se mantiene $ien agitada y Fuye hacia el exterior con la misma rapide). *i la concentración de sal presente en el depósito es de -.? EgD#, terminar la cantidad de sal presente en el depósito al ca$o de 1 minutos. G6uenta concentración alcan)ara de sal en el depósito en un tiempo de ? minH &allamos por el método de &eun con h@1.+? y ? iteraciones
= kg de sal dentro del de'osito en el instante t 3 dy y = 6 × 0.5 − 6 × = − 3 y ÷ 25 dt 50 y ( 0 ) = 5 y ( 5 ) = &
x ( t )
A. − yi +1 B. −
1 2
= yi + hf ( xi , yi ) → caso 'redictor
× f ( xi , yi ) + f ( xi +1 , yi +1 )
C. − yi +1
1 = yi + h × f ( xi , yi ) + f 2
i ( % & ' )
( x + , y + ) → caso corrector i 1
i 1
t i
A
B.
C
1.+? +.? 2.=? ?
> 1-.22=? 1+.2=1>< 14.1+2< 1?.;22
+.++ 1.<14=? 1.;4<-> 1.4+-+>4 1.++2++
=.==? 1-.1424 1+.+-4= 12.<>-1?.?-<-
E#ERCICIO$ PROPUE$TO$ 1.5%n un tanque perfectamente agitado se tiene 4--# de una salmuera en la cual est'n disueltos de sal común Ia6#", en cierto momento se hace llegar al tanque un gasto de un >- salmuera que contiene -.? Jg de sal común por litro. *i tiene un gasto de salida de >-#Dmin determine/ a" GKué cantidad de sal hay en el tanque transcurrido 1- minutosH $" GKué cantidad de sal hay en el tanque transcurrido un tiempo muy grandeH
Re.pue.ta
a)
dy
= 40 − 0.2 y
y ( 0)
dx h1
= 25
y ( 10 )
=&
r'ta y ( 10 )
= 175.9466
b)a solucion se o*tiene hasta la cantidad de sal en el tan+ue no ca#*ie con el tie#'o y ( 10 )
= 199.9914
+.5*e hacen reaccionar isotérmicamente +;- g de acetato de etilo CH3COOC2 H 5 con 1=?g de hidróxido de sodio NaOH en solución acuosa ajustando el volumen total a ? litros" para dar acetato de sodio CH 3COONa y el alcohol etílico C2 H 5OH de acuerdo con la siguiente ecuación estequiometria. CH 3COOC2 H 5 + NaOH
→ CH 3COONa + C 2H 5OH
Re.pue.ta dy
= 1.44 × 0.01 ( 0.59 − y ) ( 0.875 − y ) dx h1
y (0)
=0
y ( 30 )
=&
r'ta y ( 30 )
= 0.1692
2.56alcular el tiempo necesario para que el nivel del líquido dentro del tanque esférica con r@?m mostrando en la !gura pase de 4m a 2m.#a velocidad de salida por el ori!cio del fondo es de v@4.>
Re.pue.ta da dt
=−
h10
0.012375 a
( 10a − a 2 )
a ( 0)
=4
a ( &)
=3
r'ta y ( 100 )
= 2.9796
DI!/R!M! DE FU#O
P$EUDOCODI/O
clear all disp('METODO DE HEUN') clc syms x syms y f=inline(input('ingrese la derivada' !'s'))" x=input('ingrese el val#r de x' )" y=input('ingrese el val#r de y' )" $=input('ingrese el val#r de $' )" n=input('ingrese numer# de iteraci#nes' )" clc disp('xi ! yi ! y(i%&) ! (i%&) ')" y(i%&)cas# predict#r (i%&)cas# c#rrect#r f#r i=&n s=$%x" y&=feval(f!x!y)" $y&=$y&" y*=y%$y&" y+=feval(f!s!y*)" $y*=y+$" yn=y%(($y&%$y*),*)" fprintf('-n./&f ./0f ./0f ./0f' !x !y y=yn" x=x%$" end
!y*
!yn
)"
Leamos un ejemplo con el ejercicio 2 para asi compro$ar que las respuestas son iguales
Ahora compro$amos con un ejercicio propuesto la numero 2,ingresamos los datos
C luego nos dar' los resultados
