8.2 Runge-Kutta 2do Orden dy/dx

Referencia: Burden 5.4 p272 pdf282, Chapra 25.3 p740 pdf164, Rodriguez 9.1.7 p354

Video Revisar: 4’39» + 2 minutos

¿ Que trayectoria siguen los proyectiles?
¿ Que trayectoria siguen los aviones?
¿ Cuál es la relación entre las trayectorias?


Los métodos de Runge-Kutta  mejoran la aproximación a la respuesta sin requerir determinar las expresiones de las derivadas de orden superior, usan una corrección a la derivada tomando valores de puntos alrededor.

Por ejemplo, Runge-Kutta de 2do Orden usa el promedio entre los incrementos xi y xi+h, calculados como términos K1 y K2.

# EDO. Método de RungeKutta 2do Orden 
# estima la solucion para muestras espaciadas h en eje x
# valores iniciales x0,y0
# entrega arreglo [[x,y]]

def rungekutta2(d1y,x0,y0,h,muestras):
    tamano = muestras + 1
    estimado = np.zeros(shape=(tamano,2),dtype=float)
    # incluye el punto [x0,y0]
    estimado[0] = [x0,y0]
    xi = x0
    yi = y0
    for i in range(1,tamano,1):
        K1 = h * d1y(xi,yi)
        K2 = h * d1y(xi+h, yi + K1)

        yi = yi + (K1+K2)/2
        xi = xi + h
        
        estimado[i] = [xi,yi]
    return(estimado)

Ejercicio: Para probar el algoritmo se usa la misma ecuación del problema presentado en  ‘EDO con Taylor‘ :

y'-y -x +x^2 -1 = 0

que aplicado con Runge Kutta, se obtiene:

estimado[xi,yi]
[[ 0.          1.        ]
 [ 0.1         1.2145    ]
 [ 0.2         1.4599725 ]
 [ 0.3         1.73756961]
 [ 0.4         2.04856442]
 [ 0.5         2.39436369]]
Error máximo estimado:  0.00435758459732
entre puntos: 
[ 0.          0.00067092  0.00143026  0.0022892   0.00326028  0.00435758]
>>> 

Compare los resultados con Taylor de 2 y 3 términos.

Los resultados se muestran también en la gráfica:

Se adjunta el programa de prueba que usa la función rungekutta2(d1y,x0,y0,h,muestras)  :

# PROGRAMA PRUEBA
# Ref Rodriguez 9.1.1 p335 ejemplo.
# prueba y'-y-x+(x**2)-1 =0, y(0)=1

# INGRESO.
# d1y = y' = f, d2y = y'' = f'
d1y = lambda x,y: y -x**2 + x + 1
x0 = 0
y0 = 1
h = 0.1
muestras = 5

# PROCEDIMIENTO
puntosRK2 = rungekutta2(d1y,x0,y0,h,muestras)
xi = puntosRK2[:,0]
yiRK2 = puntosRK2[:,1]

# SALIDA
print('estimado[xi,yi]')
print(puntosRK2)

# ERROR vs solución conocida
y_sol = lambda x: ((np.e)**x) + x + x**2

yi_psol = y_sol(xi)
errores = yi_psol - yiRK2
errormax = np.max(np.abs(errores))

# SALIDA
print('Error máximo estimado: ',errormax)
print('entre puntos: ')
print(errores)

# GRAFICA [a,b+2*h]
a = x0
b = h*muestras+2*h
muestreo = 10*muestras+2
xis = np.linspace(a,b,muestreo)
yis = y_sol(xis)

# Gráfica
import matplotlib.pyplot as plt
plt.plot(xis,yis, label='y conocida')
plt.plot(xi[0],yiRK2[0],'o', color='r', label ='[x0,y0]')
plt.plot(xi[1:],yiRK2[1:],'o', color='m', label ='y Runge-Kutta 2 Orden')
plt.title('EDO: Solución con Runge-Kutta 2do Orden')
plt.xlabel('x')
plt.ylabel('y')
plt.legend()
plt.grid()
plt.show()