Práctica 3-Problemas de valores iniciales en ecuaciones diferenciales ordinarias PDF

Preview

Summary

Download Práctica 3-Problemas de valores iniciales en ecuaciones diferenciales ordinarias PDF

Description

Práctica 3. Problemas de valores iniciales en ecuaciones diferenciales ordinarias

La solución analítica de la ecuación diferencial ordinaria (EDO) de primer orden la podemos hallar con los siguientes pasos: 1. Primero habrá que expresar y´(t) como la derivada de y respecto a t (dy/dt) y con posterioridad aplicaremos el método de variables separables dejando a un lado de la ecuación “t” y al otro lado “y”: 𝑦´(𝑡) =

𝑃 ∗ 𝑦(𝑡); 𝑦(0) = 𝑋 100 𝑃 𝑑𝑦 = ∗𝑦 𝑑𝑡 100 𝑑𝑦 𝑃 ∗𝑡 = 100 𝑦

2. Cuando tenemos las variables separadas en cada lado procederemos a realizar el despeje aplicando la integral a cada uno de los lados: ∫

𝑑𝑦 𝑃 ∗ 𝑑𝑡 =∫ 100 𝑦

𝑙𝑛(𝑦) =

𝑃

100

∗𝑡+𝐶

3. Ahora aplicaremos la exponencial en ambos lados de la ecuación, lo cual nos ayudará a despejar 𝑒 𝑙𝑛 dejando “y” a un solo lado de la ecuación: 𝑃

𝑒 ln(𝑦) = 𝑒 100 ∗𝑡 ∗ 𝑒 C 𝑃

𝑦 = 𝑒100 ∗𝑡 ∗ 𝑒 C

(1)

4. Teniendo que “𝑒 C " se puede entender como una constante sustituiremos dicha expresión por "𝐶1 " y también sustituiremos el valor inicial en dicha ecuación y(0)=X: 𝑃

𝑦 = 𝑒 100 ∗𝑡 ∗ 𝐶1 𝑃

y(0)=X; X= 𝐶1 * 𝑒100

𝑋 = 𝐶1

∗0

1

5. Así pues, sustituyendo nuestro valor de 𝐶 1 en la expresión (1) tendremos la solución analítica a la EDO: 𝑃

𝑦 = 𝐶 ∗ 𝑒 100 ∗𝑡

(2)

Expresión en R de solución analítica y ejemplo: y=function (X,P,t) X*exp(P/100*t) y(5000,2.25,3) Para la representación gráfica implementaremos el siguiente código: y=rep(1,times=t) y[0]=X for (i in 1:t) { y[i+0]=y[i]*exp(P/100) } return(matrix(c(y),nrow=t+1)) } y (5000,2.25,50) windows() plot(y,type="l",xlab='x',ylab='y',col='red') lines(y,col="red") legend('topleft', c("Solución Analítica"),bty='n', col = c('red'),lty=c(1))

Método de Euler Partiendo de t0 e y(0)=5000€ con un número de saltos de longitud h=1 (o lo que es lo mismo diferencial de t) los cálculos se realizarán en bucle partiendo de la ecuación inicial hasta llegar a N que será nuestra última t. Sabiendo N, h se puede calcular como h=

𝑡𝑓−𝑡0 𝑁

Haremos un primer ejemplo del primer paso a mano y después indicaremos el código a usar en R: 𝑃

Y1=Y0+h*100 * Y0=5000+1*0.0225*5000=5112,5€ Código en R: euler=function(fnc,x0,xf,N,y0){ # y'=fnc(x,y) en (x0,xf], y(x0)=y0 utilizando n+1 puntos h=(xf-x0)/N x=seq(x0,xf,by=h) y=rep(y0,times=N+1) for (i in 1:N) { y[i+1]=y[i]+h*fnc(x[i],y[i]) } return(matrix(c(x,y),nrow=2,byrow=TRUE))

t1=t0 +h=0+1=2

} A mayor N como veremos en los resultados finales Euler ofrece un mejor resultado, ya que menor será h con lo cual de esta manera podremos disminuir el error de aproximación a la solución exacta. En la tabla final, se puede observar el resultado.

Método de Euler mejorado Da mejores resultados que Euler. El código a implementar en R es el siguiente: 

Código en R:

eulermejorado=function(fnc,x0,xf,N,y0){ h=(xf-x0)/N x=seq(x0,xf,by=h) y=rep(y0,times=N+1) for (i in 1:N) { y[i+1]=y[i]+(h/2)*(fnc(x[i],y[i])+fnc(x[i]+h,y[i])+h*fnc(x[i],y[i])) } return(matrix(c(x,y),nrow=2,byrow=TRUE)) } Método de Runge-kutta



Para la aplicación de dicho método se programará en R el siguiente código:

RK4=function(fnc,x0,xf,N,y0){ h=(xf-x0)/N x=seq(x0,xf,by=h) y=rep(y0,times=N+1) j=h/2 for (i in 1:N) { k1=fnc(x[i],y[i]) k2=fnc(x[i]+j,y+j*k1) k3=fnc(x[i]+j,y[i]+j*k2)

k4=fnc(x[i]+h,y[i]+h*k3) y[i+1]=y[i]+(h/6)*(k1+2*k2+2*k3+k4) } return(matrix(c(x,y),nrow=2,byrow=TRUE)) } Por último, implementaremos para N=10 y N=100 los métodos anteriores y analizaremos cuáles ofrecen mejores resultados. El código que vamos a usar para ello es el siguiente: source("eulermejorado.R") source("euler.R") source("RungeKutta4.R")

fun=function(x,y){return (2.25/100*y)} sole=euler(fun,0,3,10,5000) sole2=euler(fun,0,3,100,5000)

solem=eulermejorado(fun,0,3,10,5000) solem2=eulermejorado(fun,0,3,100,5000)

solrk4=RK4(fun,0,3,10,5000) solrk42=RK4(fun,0,3,100,5000)

Resumen de los resultados para h=1

Métodos EULER EULER MEJORADO Runge-Kutta

h=1 5345,150703€ 5349,121338€ 5349,147823€

100, % 0,0000747645% 0,0000053582% 0,000000649636%

El resultado analítico es 5349,151298€ por lo que el método que mejor se aproxima al resultado es el método de R unge-Kutta. De igual manera, evaluando en un N superior podemos observar

como al realizar un mayor número de iteraciones los resultados van convergiendo en todos los métodos, sobre todo para el método de Euler....