Shared3.sagewsOpen in CoCalc

Thème 3 : Solutions approchées d'une équation différentielle

Base de travail : le TD3 de 2015-16 (Y. Demay)

Egalement ces notes de cours : tableau 1, tableau 2

1. Méthode solve existante dans Sagemath ?

Lire le livre Sagemath, le tutoriel, faire une recherche sur internet "sagemath equations différentielles" ou "sagemath solve differential equations"

Exemple avec l'équation différentielle y'(x)=-(x-1)y(x) et la condition initiale y(0)=1/\sqrt{e} .

#voir par ex. http://doc.sagemath.org/html/fr/tutorial/tour_algebra.html
#desolve?
y = function('y')(x)
f(x,z)=-(x-1)*z;x0=0;y0=1/sqrt(e) #exemple. définir f(x,y) crée un conflit de variables
h(x)=desolve(diff(y,x)==f(x,y),y,ics=[x0,y0]);h
plot1=plot(h,-5,5,color='black');plot1
x |--> e^(-1/2*x^2 + x - 1/2)

2. Implémentatation de la méthode d'Euler

Ecrire une fonction prenant comme argument une fonction x,y\mapsto f(x,y) , les nombres x_0 , y_0 , h , N , rendant une liste des N+1 coordonnées (x_k,y_k)_{0\leq k\leq N} x_k=x_0+kh et (y_k) vérifie la relation de récurrence

y_{k+1}=y_k+hf(x_k,y_k)

Dessiner la réunion des segments reliant deux points (x_k,y_k) consécutifs

Essayer avec f(x,y)=y , x_0=0 , y_0=1 , différentes valeurs de h , N=2/h

#listes de points (x_k,y_k) par méthode d'Euler
def euler(f,x0,y0,h,n):
    l=[[x0,y0]]
    for i in range(n):
        x=l[i][0];y=l[i][1]
        l.append([x+h,y+h*f(x,y)]) #l=l+[[x+h,y+h*f(x,y)]]
        #print(l)
    return(l)
#test1
#f(x,y)=-(x-1)*y;x0=0;y0=1/sqrt(e)
print f,',',x0,',',y0
plot2=line(euler(f,x0,y0,0.1,50))+line(euler(f,x0,y0,0.05,100),color='green');plot2
#comparaison avec solution formelle
plot2+plot1
(x, z) |--> -(x - 1)*z , 0 , 1/sqrt(e)

Devoir

Une équation différentielle d'ordre 2 est une équation de la forme

y''=f(x,y,y')
On peut la convertir en une équation diff. vectorielle d'ordre 1 : Y'=F(x,Y) en posant

Y=\left(\begin{array}{l} y\\ y' \end{array}\right), \qquad F\big(x,\left(\begin{array}{l} y\\ z \end{array}\right)\big)=\left(\begin{array}{l} z\\ f(x,y,z) \end{array}\right)
L'équation différentielle vectorielle d'ordre 1 admet une solution approchée par la méthode d'Euler exposée plus haut.

1. Ecrire une fonction prenant comme argument une fonction (x,y,z)\mapsto f(x,y,z) , les nombres x_0,y_0,y'_0,h,N , rendant une liste de N+1 triplets (x_k,y_k,y'_k) correspondant par la méthode d'Euler à une approximation de la solution de l'équation y''=f(x,y,y') vérifiant les conditions initiales y(x_0)=y_0 et y'(x_0)=y'_0 .

Ecrire une fonction dessinant la réunion des segments reliant deux points (x_k,y_k) consécutifs.

2. Tester les deux fonctions avec l'équation

y''+y'+2y=cos(4x),\qquad x_0=0,\ y(0)=1,\ y'(0)=2

Qu'obtient t-on si on varie h ?

La fonction sagemath desolve donne t-elle une solution formelle de cette équation ?

Comment se compare le graphe de la solution formelle avec le graphe de la solution approchée de l'équation pour h=0.1, N=100 puis h=0.01, N=1000 ?