Ecuaciones diferenciales en Sage

Sage permite resolver ecuaciones diferenciales ordinarias (EDO).

EDO de primer y segundo orden

Veamos cómo resolver simbólicamente la EDO implícita x'+x-1=0:

sage: var('t') # define la variable independiente t
sage: x = function('x',t) # define x como funcion de t
sage: DE = diff(x, t) + x - 1
sage: desolve(DE, x)
(_C + e^t)*e^(-t)

Es posible indicar una condición inicial con el argumento opcional ics:

sage: desolve(DE, x, [10,2]) # condicion x(10)=2
(e^10 + e^t)*e^(-t)

La EDO puede depender de parámetros:

sage: var('a') # parametro a
sage: desolve(diff(x,t)+x-a, x, ivar=t) # t es la variable independiente
(a*e^t + _C)*e^(-t)
sage: s = desolve(diff(x,t)+x-a, x, [0,2*a], ivar=t) # solucion particular 
(a*e^t + a)*e^(-t) 
sage: s(a=3) 
3*(e^t + 1)*e^(-t) 
sage: s(a=3, t=4).n() 
3.05494691666620

El argumento opcional contrib_ode permite resolver ecuaciones de Clairaut, Lagrange, Riccati y otras:

sage: desolve(diff(x,t)^2+t*diff(x,t)-x, x, contrib_ode=True)
[x(t) == _C^2 + _C*t]

El comando desolve también permite resolver simbólicamente ecuaciones de segundo orden:

sage: var('t')
sage: x = function('x',t)
sage: desolve(diff(x,t,2)-x == t, x)
_K2*e^(-t) + _K1*e^t - t
sage: desolve(diff(x,t,2)-x == t, x, [10,2,1]) # x(10)=2, x'(10)=1
-t + 7*e^(t - 10) + 5*e^(-t + 10)

sage: desolve(diff(x,t,2)+x*(diff(x,t,1))^3==0, x)
1/6*x(t)^3 + _K1*x(t) == _K2 + t

La resolución numérica de EDO de primer orden puede realizarse con el comando desolve_odeint, que usa odeint del paquete scipy.integrate. La EDO debe estar escrita en forma explícita, x'=f(x,t), y debe indicarse como argumento el segundo miembro f. Por ejemplo, en el caso de x'+x-1=0, debemos escribir:

sage: var('t')
sage: tm = [0,0.1..2] # tiempos para el calculo de la solucion
sage: s = desolve_odeint(1-x, 2, tm, x) # cond. inicial x(0) = 2
sage: line(zip(tm,s), marker='o')

odeint

El comando desolve_rk4 resuelve numéricamente EDO de primer orden, tanto en forma implícita como explícita, con el algoritmo clásico de Runge-Kutta de cuarto orden:

sage: s = desolve_rk4(1-x, x, [0,2], t)
sage: p = line(s, marker='o', legend_label='Paso por defecto')
# mas opciones para desolve_rk4:
sage: s1 = desolve_rk4(1-x, x, ics=[0,2], ivar=t, end_points=10, step=1)
sage: p+line(s1, marker='d', color='green', legend_label='Paso 1')

rk4

Con propósitos educativos, es posible usar el comando eulers_method para resolver numéricamente una EDO de primer orden con el método de Euler.

Para algoritmos numéricos más avanzados, puede usarse ode_solver, de la librería científica GNU.

Sistemas de EDO de primer orden

Para resolver simbólicamente sistemas de EDO de primer orden, puede usarse el comando desolve_system:

sage: var('t')
sage: x = function('x', t)
sage: y = function('y', t)
sage: de1 = diff(x,t) + y - 1 == 0
sage: de2 = diff(y,t) - x + 1 == 0
sage: desolve_system([de1, de2], [x,y])
[x(t) == (x(0) - 1)*cos(t) - (y(0) - 1)*sin(t) + 1,
y(t) == (y(0) - 1)*cos(t) + (x(0) - 1)*sin(t) + 1]

sage: s = desolve_system([de1, de2], [x,y], [0,1,2]) # x(0)=1, y(0)=2
[x(t) == -sin(t) + 1, y(t) == cos(t) + 1]

El comando desolve_odeint permite resolver numéricamente sistemas de primer orden como, por ejemplo, la ecuación de Lotka-Volterra:
\begin{cases}\displaystyle\frac{\mathrm{d}x}{\mathrm{d}t}=x(1-y), & x\text{ presas,}\\[0.6 em] \displaystyle\frac{\mathrm{d}y}{\mathrm{d}t}=-y(1-x), & y\text{ predadores.}\end{cases}

sage: var('x y')
sage: des = [x*(1-y), -y*(1-x)]
sage: tm = [0,0.1..10] # tiempos de calculo
sage: s = desolve_odeint(des, [0.5,2], tm, [x,y]) # x(0)=0.5, y(0)=2
sage: line(zip(tm, s[:,0]), legend_label='Presa') + \
....: line(zip(tm, s[:,1]),color='green', legend_label='Predador')

sis-ode-intEl comando desolve_system_rk4 permite resolver numéricamente sistemas de primer orden con el algoritmo clásico de Runge-Kutta de cuarto orden:

sage: var('t')
sage: r = desolve_system_rk4(des, [x,y], [0,0.5,2], t)
sage: presa = [(i,j) for i,j,k in r]
sage: predador = [(i,k) for i,j,k in r]
sage: line(presa, legend_label='Presa') + \
....: line(predador, color='green', legend_label='Predador')

Con propósitos educativos, es posible usar el comando eulers_method_2x2 para resolver numéricamente un sistema de EDO de primer orden con el método de Euler.

Bibliografía

  1. Sage Tutorial v6.3. A Guided Tour
  2. Sage Reference Manual. Solving ordinary differential equations
  3. Sage Constructions v6.3. Calculus
  4. Sage Quickstart for Differential Equations
Advertisements
Esta entrada foi publicada en Sage. Ligazón permanente.

Deixar unha resposta

introduce os teu datos ou preme nunha das iconas:

Logotipo de WordPress.com

Estás a comentar desde a túa conta de WordPress.com. Sair / Cambiar )

Twitter picture

Estás a comentar desde a túa conta de Twitter. Sair / Cambiar )

Facebook photo

Estás a comentar desde a túa conta de Facebook. Sair / Cambiar )

Google+ photo

Estás a comentar desde a túa conta de Google+. Sair / Cambiar )

Conectando a %s