Sdgdf

Solo disponible en BuenasTareas
  • Páginas : 7 (1681 palabras )
  • Descarga(s) : 0
  • Publicado : 26 de junio de 2010
Leer documento completo
Vista previa del texto
Ecuaciones en Derivadas Parciales y Análisis Numérico Prácticas
Capítulo 2. Ecuaciones diferenciales ordinarias (EDOs). 2.1 Resolución de una ecuación diferencial ordinaria.
Vamos a resolver numéricamente la siguiente ecuación diferencial:

dy =y 1−y  dt
con condición inicial y(0) = 0.2. Para empezar vamos a utilizar el método de Euler explícito. Este algoritmo, que estudiaréis en teoría con más detalle, se basa en la aproximación:

y tdt ≈y t dt ∙ y ' t =y tdt ∙ y t1−y t 
para un valor de dt pequeño. Para empezar, tomamos un conjunto de instantes de tiempo  t n=t 0 n dt y  definimos  y 0=y  0 . A continuación, definimos el valor de los puntos y n de forma recursiva,  mediante el esquema numérico:

y n1=y n dt ∙ y n 1−y n Lo implementamos del siguiente modo tsteps=1000; % Subdivisiones del intervalo t_f=10; % valor final de tiempo dt=t_f/tsteps;  % el paso de tiempo t=0:dt:t_f; % definimos los diferentes tiempos y=zeros(1,tsteps+1);% inicializamos y y(1)=0.2;  % definimos dato inicial for n=1:tsteps y(n+1)=y(n)+dt*y(n)*(1­y(n)); end  %metodo de Euler explicito Dibujamos los resultados plot(t,y);Matlab tiene implementados varios métodos numéricos que podemos aplicar directamente. Crea un fichero  de Matlab (File/New/M­file) que contenga las siguientes líneas: function deriv=edo(t,y) deriv=y*(1­y); Guarda el archivo con el nombre edo.m y establece el directorio donde lo has guardado como directorio de  trabajo. Vuelve a la línea de comandos y teclea:

t_0=0;  t_f=10; y0=0.2; [t,y]=ode23('edo',[t_0,t_f],y0); plot(t,y)La function Matlab ode23 es una implementación de un método Runge­Kutta de orden (2,3). En este caso  no es necesario poner el paso de tiempo porque la misma función ode23 lo va reduciendo sobre la marcha si  es necesario. Otra function Matlab es ode45 basado en un método explícito Runge­Kutta de orden (4,5). En  general, ode45 es la mejor función para aplicar en un primer momento a muchos de los problemas:[t,y]=ode45('edo',[t_0,t_f],y0); plot(t,y) En este caso, es fácil comprobar que la función:

y t=

et 4et

es la solución exacta del problema de valor inicial anterior. Podemos usar esta referencia para comprobar qué  tal precisión hemos obtenido con los distintos métodos. y_exacta=exp(t)./(4+exp(t)); plot(t,abs(y­y_exacta)); %repetir cuando y se ha calculado con % distintos metodos

2.2 Errores.Vamos a estudiar la precisión que obtenemos en función del paso de tiempo utilizado. Para empezar creamos  una función en matlab que devuelva sólo el valor aproximado de y(1) tomando como argumento el número de  pasos temporales: function y=euler1(tsteps) t_f=1; % valor final de tiempo dt=t_f/tsteps;  % el paso de tiempo t=0:dt:t_f; % definimos los diferentes tiempos y=0.2;  % definimos dato inicialfor n=1:tsteps y=y+dt*y*(1­y); end  %metodo de Euler explicito Ahora podemos hacer una gráfica de la precisión obtenida (comparamos con el valor exacto,

e ≈0,404609675 ) en función del paso temporal dt. 4e
yexacto=exp(1)/(4+exp(1)); pasos=[];errores=[]; %vectores con los pasos temporales y los %errores for k=[4 10 30 100 300 1000] pasos=[pasos, 1/k]; %engordamos los vectores sobre la marchaerrores=[errores, abs(euler1(k)­yexacto)]; end plot(pasos, errores)

En la siguiente gráfica a la izquierda puedes ver el resultado, mientras que a la derecha se muestra la gráfica  correspondiente a un método de Runge­Kutta de orden 2 (el método del punto medio). Como puedes ver, la  convergencia del método de Euler es lineal mientras que la del método del punto medio es cuadrática (la gráfica es aproximadamente una parábola).

En general, para estimar el orden de convergencia  de un método, se escribe el error como una  potencia del paso de tiempo: y se toma logaritmos para obtener una ecuación  lineal:

error=C ∙ dt 

log  error= ∙ logdt logC 
Al aplicar este procedimiento al error del método  de Euler, se obtiene una recta (a la derecha, en ...
tracking img