Ejemplos de solución numérica de EDO de orden superior

En estos ejemplos mencionaremos los comandos más frecuentemente usados para la solución numérica de EDO de orden superior de los programas de software MaximaPython y Matlab, sin entrar en detalles sobre los mismos.  Para una mayor información, usted deberá consultar los manuales de cada programa.

Si requiere información para solucionar numéricamente EDO de primer orden, visite:

Ejemplos de solución numérica de EDO de 1er orden.

Ejemplo 1: solución de EDO de orden superior.

Las herramientas para el cálculo numérico de sistemas de ecuaciones diferenciales (ED) de primer orden, nos permiten obtener de una manera fácil, la solución para una ED de orden superior.

Solucionar el PVI: \frac{d^2y}{dx^2}+x\frac{dy}{dx}+20\sin(y)=0, con y(0)=0 y y'(0)=1.

Solución:

Si hacemos y_1=y y y_2=\frac{dy_1}{dx}, la ED planteada se puede expresar como un conjunto de ED así:

\frac{dy_1}{dx}=y_2

\frac{dy_2}{dx}=-xy_2-20\sin(y_1)

El PVI se solucionará en el intervalo 0<=x<6.  Observe que el TEU garantiza solución en dicho intervalo.

Python:

Para solucionar el sistema planteado con Python, usamos el comando odeint(df,y0,x) de la librería scipy.integrate. Aquí, df es la expresión al otro lado del igual del sistema de ecuaciones diferenciales escritas como \dfrac{dy}{dx}=..., y0 es el vector de condiciones iniciales y x es el vector de variable independiente(dominio) en donde se calculará la solución.  De esta manera, para el sistema planteado el código será algo como:

from matplotlib import *python
from numpy import *
from scipy.integrate import odeint
# Definicion de E.D a resolver: y''+x*y'+20*sin(y)=0
def df(y,x):
    y1, y2= y[0], y[1] #variables de paso
    dy1=y2  # derivada variable 1
    dy2=-x*y2-20*sin(y1) # derivada variable 2
    return [dy1,dy2]

# Condiciones iniciales 
y0 =[0,1]

# Definicion del rango 
x = linspace(0,6,500)

# Para mayor info: help(odeint)
# solucion numerica
sol = odeint(df, y0, x)
y=sol[:,0] #toma el vector correspondiente a la solucion de y1
plot(x,y)

La gráfica obtenida es:

Sln numerica ED

Máxima:

En este caso, usamos el comando rk(ode,var,initial, domain).  En donde ode, al igual que Python, es la expresión al otro lado del igual del sistema de ecuaciones definidas como \frac{dy}{dx}=..., var son las variables a resolver, initial es el vector de condiciones iniciales y domain el vector del dominio de solución.

El código es:

sol:rk([-x*y2-20*sin(y1),y2],[y2,y1],[1,0],[x,0,6,0.1])$ Maxima
y:makelist([sol[k][1],sol[k][3]],k,1,length(sol))$
plot2d([discrete,y]);

Cuyo resultado gráfico es:

Sln numérica con máxima

En este caso, el comando rk() entrega “tripletas” [xk,y2k,y1k]. De tal forma que para graficar, se debe construir un vector o “lista”, con el comando makelist() con los n elementos de la 1era y 3era “columna”de la solución.

Matlab:

El ejemplo dado se puede resolver mediante la secuencia de comandos:

>>f=@(x,y) [y(2);-x*y(2)-20*sin(y(1))]; Matlab
xspan=[0 6]; y0=[0,1]; 
[x,y]=ode45(f,xspan,y0);
plot(x,y(:,1))

 

Para más información pueden consultar esta presentación en la carpeta de recursos:

Responder

Introduce tus datos o haz clic en un icono para iniciar sesión:

Logo de WordPress.com

Estás comentando usando tu cuenta de WordPress.com. Cerrar sesión / Cambiar )

Imagen de Twitter

Estás comentando usando tu cuenta de Twitter. Cerrar sesión / Cambiar )

Foto de Facebook

Estás comentando usando tu cuenta de Facebook. Cerrar sesión / Cambiar )

Google+ photo

Estás comentando usando tu cuenta de Google+. Cerrar sesión / Cambiar )

Conectando a %s