import numpy as np
from .config import est_pedagogique
[docs]
def euler(f, y0, t0, tf, h, pedagogique=None):
"""
Resout une ODE par la methode d'Euler explicite.
Args:
f (callable): La fonction derivee y' = f(t, y).
y0 (float/array): Condition initiale.
t0 (float): Temps initial.
tf (float): Temps final.
h (float): Pas de temps.
pedagogique (bool, optional): Affiche les premieres etapes si True.
Returns:
tuple: (t, y) ou t est le tableau des temps et y le tableau des solutions.
"""
t = np.arange(t0, tf + h, h)
y = np.zeros((len(t), len(np.atleast_1d(y0))))
y[0] = y0
for i in range(len(t) - 1):
y[i+1] = y[i] + h * f(t[i], y[i])
if est_pedagogique(pedagogique) and i < 5:
print(f"[Pedagogique] Euler Etape {i+1}: y({t[i+1]:.3f}) ≈ {y[i+1]}")
return t, y
def euler_modifie(f, y0, t0, tf, h, pedagogique=None):
t = np.arange(t0, tf + h, h)
y = np.zeros((len(t), len(np.atleast_1d(y0))))
y[0] = y0
for i in range(len(t) - 1):
k1 = f(t[i], y[i])
k2 = f(t[i+1], y[i] + h * k1)
y[i+1] = y[i] + (h/2) * (k1 + k2)
return t, y
[docs]
def rk2(f, y0, t0, tf, h, pedagogique=None):
t = np.arange(t0, tf + h, h)
y = np.zeros((len(t), len(np.atleast_1d(y0))))
y[0] = y0
for i in range(len(t) - 1):
k1 = f(t[i], y[i])
k2 = f(t[i] + h/2, y[i] + (h/2) * k1)
y[i+1] = y[i] + h * k2
return t, y
[docs]
def rk4(f, y0, t0, tf, h, pedagogique=None):
"""
Resout une ODE par la methode de Runge-Kutta d'ordre 4.
Args:
f (callable): La fonction derivee y' = f(t, y).
y0 (float/array): Condition initiale.
t0 (float): Temps initial.
tf (float): Temps final.
h (float): Pas de temps.
pedagogique (bool, optional): Non utilise pour l'instant.
Returns:
tuple: (t, y) ou t est le tableau des temps et y le tableau des solutions.
"""
t = np.arange(t0, tf + h, h)
y = np.zeros((len(t), len(np.atleast_1d(y0))))
y[0] = y0
for i in range(len(t) - 1):
k1 = f(t[i], y[i])
k2 = f(t[i] + h/2, y[i] + (h/2) * k1)
k3 = f(t[i] + h/2, y[i] + (h/2) * k2)
k4 = f(t[i+1], y[i] + h * k3)
y[i+1] = y[i] + (h/6) * (k1 + 2*k2 + 2*k3 + k4)
return t, y
[docs]
def adams_bashforth(f, y0, t0, tf, h, ordre=2):
t = np.arange(t0, tf + h, h)
y = np.zeros((len(t), len(np.atleast_1d(y0))))
y[0] = y0
if ordre > 1:
_, y_boot = rk4(f, y0, t0, t0 + (ordre-1)*h, h)
y[:ordre] = y_boot[:ordre]
for i in range(ordre - 1, len(t) - 1):
if ordre == 1: y[i+1] = y[i] + h * f(t[i], y[i])
elif ordre == 2: y[i+1] = y[i] + (h/2) * (3 * f(t[i], y[i]) - f(t[i-1], y[i-1]))
elif ordre == 3: y[i+1] = y[i] + (h/12) * (23*f(t[i], y[i]) - 16*f(t[i-1], y[i-1]) + 5*f(t[i-2], y[i-2]))
elif ordre == 4: y[i+1] = y[i] + (h/24) * (55*f(t[i], y[i]) - 59*f(t[i-1], y[i-1]) + 37*f(t[i-2], y[i-2]) - 9*f(t[i-3], y[i-3]))
return t, y
[docs]
def adams_moulton(f, y0, t0, tf, h, ordre=2):
t = np.arange(t0, tf + h, h)
y = np.zeros((len(t), len(np.atleast_1d(y0))))
y[0] = y0
if ordre > 1:
_, y_boot = rk4(f, y0, t0, t0 + (ordre-1)*h, h)
y[:ordre] = y_boot[:ordre]
for i in range(ordre - 1, len(t) - 1):
if ordre == 2:
y_pred = y[i] + (h/2) * (3 * f(t[i], y[i]) - f(t[i-1], y[i-1]))
y[i+1] = y[i] + (h/2) * (f(t[i+1], y_pred) + f(t[i], y[i]))
elif ordre == 3:
y_pred = y[i] + (h/12) * (23*f(t[i], y[i]) - 16*f(t[i-1], y[i-1]) + 5*f(t[i-2], y[i-2]))
y[i+1] = y[i] + (h/12) * (5*f(t[i+1], y_pred) + 8*f(t[i], y[i]) - f(t[i-1], y[i-1]))
return t, y