introduction

Une équation différentielle normale est une équation qui comprend les variables indépendantes t, la fonction y (t) qui dépend de t, et sa nième dérivée (n = 0,1,2, ..., N). En d'autres termes

g(t,y,y',y'',y^{(3)},...,y^{(N)}) = 0 \hspace{50pt}(1)

C'est une équation qui peut être décrite sous la forme de.

Ces équations peuvent être approchées par une solution numérique appelée méthode Eular. Je vais expliquer en utilisant un exemple réel.

Se préparer à utiliser la méthode Eular

L'équation différentielle ci-dessus $ (1) $ est en fait sous une forme qui ne convient pas à l'utilisation de la méthode Euler. Vous devez en quelque sorte transformer cela en quelque chose comme $ (2) $ ci-dessous.

\frac{dx(t)}{dt} = f(t,x(t)) \hspace{50pt}(2)
\\pourtant, x(t)=\begin{pmatrix}y\\y'\\y''\\...\\y^{(N-1)}\end{pmatrix}

Exemple réel

Considérons un système avec une masse $ m $ à la pointe d'un ressort avec une constante de ressort de $ k $. Lorsqu'une force vibratoire externe est appliquée au système et que la résistance de l'air est prise en compte, l'équation différentielle est la suivante.

Il est possible de le résoudre analytiquement, mais nous envisagerons de le résoudre par une méthode de calcul numérique.

m\frac{d^2y(t)}{dt^2}+2\gamma \frac{dy(t)}{dt}+ky(t)=a\sin(\omega t)\hspace{50pt}(3)

Lorsque $ (3) $ est transformé,

\frac{d^2y(t)}{dt^2}=-\frac{2\gamma}{m} \frac{dy(t)}{dt}-\frac{k}{m}y(t)+\frac{a}{m} \sin(\omega t) \hspace{50pt}(3-1)

Aussi,

\frac{dy(t)}{dt}=\frac{dy(t)}{dt}\hspace{50pt}(3-2)

Où le vecteur $ x (t) $

x(t)=\begin{pmatrix} x_1(t)  \\ x_2(t) \end{pmatrix} = \begin{pmatrix} y(t)  \\ \frac{dy(t)}{dt} \end{pmatrix}

Si vous mettez, à partir de $ (3-1) $, $ (3-2) $,

\frac{dx(t)}{dt}
=\frac{d}{dt} \begin{pmatrix} y(t) \\ \frac{dy(t)}{dt}  \end{pmatrix}
=\begin{pmatrix} \frac{dy(t)}{dt} \\ \frac{d^2y(t)}{dt^2}  \end{pmatrix}
=\begin{pmatrix} \frac{dy(t)}{dt} \\ -\frac{2\gamma}{m} \frac{dy(t)}{dt}-\frac{k}{m}y(t)+\frac{a}{m} \sin(\omega t)  \end{pmatrix}
=\begin{pmatrix} x_2(t) \\ -\frac{2\gamma}{m} x_2(t)-\frac{k}{m}x_1(t)+\frac{a}{m} \sin(\omega t)  \end{pmatrix}

Par conséquent, il peut être exprimé comme suit.

\frac{dx(t)}{dt} = f(t,x(t))\hspace{50pt}(4)
\\where : f(t,x)=\begin{pmatrix} x_2 \\ -\frac{2\gamma}{m} x_2-\frac{k}{m}x_1+\frac{a}{m} \sin(\omega t)  \end{pmatrix}

Méthode Eular

Grâce aux travaux réalisés jusqu'à présent, nous avons pu créer une belle équation différentielle. Maintenant, parlons de la méthode Eular elle-même.

Premièrement, à partir de la définition de la différenciation

\frac{dx(t)}{dt} = f(t,x(t))

⇒\lim_{\Delta t \to 0} \frac{x(t+\Delta t)-x(t)}{\Delta t} = f(t,x(t))

Ici, si vous utilisez $ \ Delta t $ qui est assez petit,

\frac{x(t+\Delta t)-x(t)}{\Delta t} \fallingdotseq f(t,x(t))

⇒x(t+\Delta t) \fallingdotseq x(t)+f(t,x(t))\Delta t\hspace{50pt}(5)

En utilisant cette formule, si vous connaissez le premier $ x (t) $, vous pouvez calculer $ x (t + \ Delta t) $ après $ \ Delta t $ secondes.

Programme (python)

La méthode Eular est implémentée en python comme suit.

`f.py`


import numpy as np

a = 10           #a(Amplitude de la force externe)
(m,gm2,k) = (30,5,1)#m,2Γ,k
omg = np.pi     #ω

def f(t,x):
    f0 = x[1]
    f1 = -(gm2/m)*x[1] -(k/m)*x[0] + (a/m)*np.sin(omg*t)
    return np.array([f0,f1])

`f_D.py`


import numpy as np
from f import f as f

dt = 0.01

def f_D_eular(t,x):
    return (t + dt, x + f(t,x)*dt)

`main.py`


import numpy as np
from f_D import f_D_eular as f_D

def main():
    (t0,x0)=(0,np.array([1,0]))
    t1 = 100
    (t,x)=(t0,x0)
    while(t < t1):
       print(t,x[0])
       (t,x)=f_D(t,x)

main()

Ce sera comme ça. En passant, le graphique du résultat de l'exécution de ce programme est le suivant.