Runge Kutta 4: Equation différentielle

[Ch00Ch00]

Hello,
J'ai un TP à faire sur la résolution d'équation différentielle avec plusieurs méthodes, c'est-à-dire Euler, Taylor2 et RungeKutta4. 
Voici l'équation différentielle : 

y′+2y=t²

avec une condition initiale : y(0)=1. 
En résolvant cette équation, on obtient 

y(t)=3/4*e^-2+1/2t²−1/2*t+1/4

 

function[y]=yRungeKutta4(t,pas)
    y(1)=1
    for i=1:(length(t)-1)
        k1=
        k2=
        k3=
        k4=
        y(i+1)=y(i)+(pas/6)*(k1+2*k2+2*k3+k4)
    end
endfunction

Je ne vois pas comment déterminer k1. Or, k1=f(ti,yi). Pourriez vous m'expliquer comment déterminer k1, k2 et k3, et k4 s'il vous plait. En regardant sur internet, je ne comprend toujours pas. 

Je vous remercie d'avance pour vos aides. 
Cordialement, 

 

P.S: Ce post a aussi été posté sur un autre site. 

[julesx]

Si ce n'est pas trop tard...

Je ne suis pas compétent pour juger de la syntaxe de ta function[y], je vois simplement que tu écris y(1)=1 alors qu'au départ, la condition initiale est y(0)=1. A noter également que, dans la solution, tu as oublié un t dans l'exponentielle. Cela dit, pour les k_i...

Comme il s'agit d'un calcul en boucle, il y a systématiquement des valeurs initiales et des valeurs finales pour t et y. Je les note ti, yi, tf et yf. Avec ces notations, on a

k1=ti²-2yi

k2=(ti+pas/2)²-2*(yi+pas*k1/2)

k3=(ti+pas/2)²-2*(yi+pas*k2/2)

k4=(ti+pas)²-2*(yi+pas*k3)

tf=ti+pas

yf=yi+pas*(k1+2*k2+2*k3+k4)/6

A noter que, si tu n'as pas besoin à la fin des valeurs initiales, tu peux écrire directement

ti=ti+pas

yi=yi+pas*(k1+2*k2+2*k3+k4)/6

et te passer de tf et yf.

Moi, j'ai gardé les deux pour pouvoir tracer les segments [(ti,yi) ;(tf;yf)] avec Algobox.

[julesx]

Toujours au cas où...

Comme j'avais fait une simulation avec Algobox, qui ne permet que de définir des fonctions d'une seule variable, je n'avait pas approfondi ma réflexion. J'y reviens donc !

Si ton langage de programmation autorise l'utilisation de fonctions de plusieurs variables, tu définis f(t,y)=t²-2y. Ceci te permet de simplifier les écritures des ki comme suit :

Toujours avec mes notations,

k1=f(ti,yi)

k2=f(ti+pas/2,yi+pas*k1/2)

k3=f(ti+pas/2,yi+pas*k2/2)

k4=f(ti+pas,yi+pas*k3)

Le reste est inchangé. Par contre, j'avais oublié de préciser qu'à la fin, si tu as utilisé tf et yf, il faut "réinitialiser" les variables par ti=tf et yi=yf.

 

 

[Ch00Ch00]

Bonsoir, 

Merci pour votre réponse julesx. Vos explications sont génials !!! J'ai enfin compris la méthode de Runge Lutta 4 grace à vous  et je vous en remercie. 

Pourriez-vous vérifier si j'ai bon pour un autre exercice du même genre s'il vous plait ? 

J'ai du mal pour la méthode de Taylor 2. 

On sait que:

pour la méthode Euler, on a: y(i+1)= y(i) + pas*f(t_i,y_i)

pour la méthode de Taylor 2, on a: y(i+1)= y(i)+pas*(-y(i)+t(i)+1)+(pas²/2)*(y(i)-t(i))

Merci,

 

Pour la méthode de Runge Kutta4, celà signifie que c'est d'ordre 4, qui combine quatre itérations successives. Sa précision est la meilleure comparé aux autres méthodes, Euler et Taylor 2 ? Lors de la résolution, si on compare la méthode analytique (c'est à dire faire le calcul sur une feuille) et la méthode Runge Kutta 4, le taux d'erreur doit être très proche de 0% 

Est-ce bien ça ? 

Merci,

[julesx]

Oui pour Euler et Runge Kutta 4 dans le cas de la deuxième équation différentielle, encore que, pour RK4, on aurait pu diminuer un peu le nombre de parenthèses, par exemple -(1/(R*C)) s'écrit plus simplement -1/(R*C). Par ailleurs, utilise uc(i) à la place de y(i). Oui aussi pour ta conclusion, RK4 est évidemment plus précis que Euler.

Par contre, pour Taylor 2, je n'ai pas trouvé d'informations simples. Je ne comprends pas ta relation y(i+1)= y(i)+pas*(-y(i)+t(i)+1)+(pas²/2)*(y(i)-t(i)) qui n'est pas homogène. Tu additionnes en particulier la fonction y, le temps t et 1 qui est sans dimension.

Un intervenant plus compétent ?

 

[Ch00Ch00]

Bonsoir,

Merci pour votre réponse julesx,

Je vais donc dominuer le nombre de parenthèse, il est vrai que c'est '' lourd '' de voir autant de parenthèse inutile. 

Pour la méthode de Taylor 2, en cours on nous a expliqué que cette mééthode, on peut définir le développment de Taylor d'une fonction, on obtient: 

f(x+h)=f(x)+hf'(x) + ... + (hⁿ)/(n!)*f(x)ⁿ 

h étant le pas. Cette méthode est appliqué pour approcher le résultat de la résolution directe.  

y(i+1)= y(i)+pas*(-y(i)+t(i)+1)+(pas²/2)*(y(i)-t(i))

Or (-y(i)+t(i)+1) = f'  et (y(i)-t(i)) = f'' 

A moins que j'ai mal recopié. 

 

Merci, 

[julesx]

il y a 13 minutes, Ch00Ch00 a dit :

y(i+1)= y(i)+pas*(-y(i)+t(i)+1)+(pas²/2)*(y(i)-t(i))

Or (-y(i)+t(i)+1) = f'  et (y(i)-t(i)) = f'' 

J'ai vu sur la toile quelque chose de ce style mais avec f' et f" permuté et surtout, qui s'appliquait au cas particulier y'=t-y. Je ne vois donc pas comment l'adapter aux deux équations de ton exercice.

 

[julesx]

La nuit portant conseil ?

Sans garantie à 100% (confirmation ou infirmation par un autre intervenant), voilà le résultat de mes cogitations.

Avec un développement  au deuxième ordre, Taylor donne

y(x+h)=y(x)+h*y'(x)+h²/2*y"(x)

qu'on transcrit sous forme de suite comme pour Euler.

En me limitant ici au cas  y'=t²-2y, on trouve y" en dérivant cette équation :

y"=2t-2y'

soit, avec y'=t²-2y,

y"=2t-2t²+4y.

Ça rejoint ce que tu avais vu en cours, mais avec les restrictions que j'avais faites précédemment (rappel, interversion de f' et f" et autre second membre).

En termes de relation de récurrence, cela donne donc

y(i+1)=yi+pas*(ti²-2yi)+pas²/2*(2ti-2ti²+4yi)

qu'on peut simplifier un tout petit peu

y(i+1)=yi+pas*(ti²-2yi)+pas²*(ti-ti²+2yi).

N.B.: J'ai fait une simulation avec Algobox, ça semble correct.

 

[Ch00Ch00]

Bonjour,

Merci pour votre aide julesx, 

J'ai verifié sur le logicel Scilab, et c'est exactement ça !!! Merci beaucoup. 

Au début, j'ai dérivé la solution générale de l'équation différentielle, j'obtenais des courbes assez bizzare alors qu'il fallait seulement dérivé à partir de l'équation différentielle de départ. 

Si je comprends bien, la méthode de Taylor 2 fait appel à la méthode d'Euler (pour y'=  y(i+1)= y(i) + pas*f(t_i,y_i)), en suivant le développement limitée, il faut déterminer la dérivée seconde. En fonction de la dérivée seconde on réinjecte y' dans cette équation. 

Dans le cas pour:  

Si on dérive une deuxième, on obtient u_c''= (-1/(R*C))*(u_c') ? 

On obtient alors:  u_c''= (-1/(R*C))*((-1/(R*C))(u_c))+(u₀/(R*C))

Est-ce bien ça ? 

 

Merci beaucoup !!!!  

 

[julesx]

Oui, mais, encore une fois, simplifie l'écriture en supprimant les parenthèses inutiles. Par ailleurs, ici, vu l'écriture "en ligne", il vaut mieux, lorsque c'est possible, rejeter les dénominateurs à la fin :

u_c''= u_c'/(R*C)

u_c''= (u_c+u₀)/(R*C)²

J'ai également fait une simulation pour cette équation, là encore, c'est cohérent.

[Ch00Ch00]

Merci beaucoup pour votre aide !

 

[julesx]

De rien, ce fut un plaisir de correspondre avec quelqu'un qui réagissait aux réponses, contrairement à beaucoup d'autres demandeurs qui ne donnent plus signe de vie après le premier message.

Meilleurs vœux pour 2018.

[Ch00Ch00]

Bonsoir, 

Il se fait un peu tard mais: Bonne et heureuse année 2018 !!!

[julesx]

Pas de problème, on a droit à tout le mois de janvier pour envoyer ses vœux ! Cela dit, merci beaucoup.

[Ch00Ch00]

Bonsoir,

Grâce à vous, j'ai pu comprendre comment programmer Runge Kutta 4 et Taylor 2. Pourriez-vous jeter un coup d'oeil sur ce que j'ai fait s'il vous plait ?

 (équation différentielle)

D'après mon professeur, je dois obtenir des courbes sinusoïdales. Lorsque je lance, je trouve des courbes assez bizarres. Je pense avoir fait une erreur, mais je ne sais pas où.

 

Merci d'avance, 

 

 

 

[julesx]

Désolé, je n'ai jamais pratiqué ce langage de programmation, je ne pourrais donc pas te dire s'il y a des erreurs de syntaxe ou autres. La seule chose que je constate, c'est qu'a priori, tu n'as donné aucune valeur numérique aux différents paramètres. Pour info, voilà ce que ça donne avec RC=1 s, uo=1 V et w=1 rad/s, soit y'+y=sin(t), avec la condition initiale y(0)=0.

[Ch00Ch00]

Bonsoir, 

Merci pour votre réponse. 

On a: u_c '(t)= (-1/(R*C))*u_c(t)+(u₀.sin(wt))/(R*C)

u''_c(t)= (-1/(R*C))*u'_c(t) 

=  (-1/(R*C))*((-1/(R*C))*u_c(t)+(u₀.sin(wt))/(R*C))

= (1/(R*C)^2).u_c(t)- (u₀.sin(wt))/(R*C)²) 

Est-ce bien ça ? 

[julesx]

Pas tout à fait, car, comme u₀.sin(wt) est aussi une fonction du temps, il faut rajouter sa dérivée. Je ne l'avais pas fait dans mes réponses précédentes car, implicitement, j'avais considéré u₀ comme constant.

u_c '(t)=- u_c(t)/(RC)+u₀.sin(wt)/(RC)

=>

u_c"(t)= -u'_c(t)/RC+wu₀.cos(wt)/(RC)

soit

u_c"(t)= -(-u_c(t)/(RC)+u₀.sin(wt)/(RC) )/RC+wu₀.cos(wt)/(RC)=u_c(t)/(RC)²-u₀.sin(wt)/(RC)²+wu₀.cos(wt)/(RC)

sauf erreur...

 

 

[julesx]

J'ai regardé un peu comment fonctionne Scilab. Vu ton code, t est un tableau. Comme est-ce que tu le définit ? Et pour l'affichage de la courbe, quelle est la syntaxe ?

[Ch00Ch00]

Bonjour julesx,

Au final, j'ai réussi grace à vous, car j'avais mal dérivé ma fonction pour Taylor 2. Et le problème des parenthèses. 

Sur Scilab, pour l'affichage des courbes: 

#Plot : afficher les courbes
plot(t,y,'k-')  #k: code couleur: noir
plot(t,y1,'b-') #b: code couleur: bleu
plot(t,y2,'g-') #g: code couleur: vert
plot(t,y3,'r-') #r: rouge
xgrid(2) #grille 
xtitle('y=exp(-t)+t') #Titre de la figure affichée
xlabel('t') #axe des abscisses
ylabel('y(t)') #axe des ordonnées 
length(t) 
legend('Analytique','Euler','Taylor','RungeKutta4') #légende

avec y,y1,y2,y3 les fonctions définies précédemments. 

P.S: Je suis désolée d'avoir mis du temps à répondre

[julesx]

Pas de problème pour le temps mis à répondre, je sais bien bien qu'en tant qu'étudiante, on ne fait pas toujours ce qu'on veut.

Par contre, tu n'as pas mentionné comment tu définis t.

[Ch00Ch00]

On définie t comme ceci: 

vi=input("donner la valeur initiale:")
vf=input("donner la valeur finale :")
pas=input("donner la valeur du pas:")
t=[vi:pas:vf]
disp(t)

 

[julesx]

OK, je m'attendais d'ailleurs à quelque chose de ce genre car c'est de cette façon que j'avais programmé mes simulations.

Je suppose évidemment que tu introduis également avec des input les valeurs de R, C, w et uo.

[Ch00Ch00]

Oui, c’est bien cela.

[julesx]

OK, bonne continuation.