Méthode Euler Explicite

**Andromedae** · 09/11/2015, 18h33

Bonjour tout le monde;

Voilà, j'ai un petit projet à faire avec l'intégration numérique par plusieurs méthodes.
Je vois bien le principe je dirais analytique (crayon et papier), mais quand à programmer cela j'ai bien plus de mal.

Voilà l'équation, la méthode sur papier que j'ai trouvé et un début de code dont j'aurais besoin de la plus grande aide ^^

L'équation correspond à la première des 3.

Merci d'avance !

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
# usr/bin/python 
 
 
 
from pylab import *
 
import numpy as np
 
 
 
#######################################
 
# ORBITES AUTOUR D'UN TROU NOIR DE SCHWARZSCHILD #
 
#######################################
 
 
 
 
# -------------------------------- #   
#    Methode d'Euler explicite   #
# -------------------------------- #
 
 
Rs = 2				# Rayon de Swharzschild
m = 1				# Masse quelconque d'une particule
c = 1				        # Vitesse lumiere
r = 0.8				# rayon r
E = 1				        # Energie particule
 
 
 
F = ( (1-(Rs/r))**(-1) ) * E/(m*c**2)
 
 
 
def Euler(F, y0, t):

**marco056** · 09/11/2015, 20h36

Il y a pléthore de solutions grâce à moteur de recherche basique.

**Andromedae** · 09/11/2015, 20h58

J'ai déjà tapé Euler sur ce forum, 5 sujets, aucuns ne m'aident réellement :/

**marco056** · 09/11/2015, 22h25

Envoyé par Andromedae

J'ai déjà tapé Euler sur ce forum, 5 sujets, aucuns ne m'aident réellement :/

Pas très perspicace :
https://www.google.fr/search?client=...D9Ss8wft_4SQBA
Il suffit de taper "méthode euler python" et il tombe des dizaines de réponses pertinentes (les solutions sont même données !).

**Andromedae** · 10/11/2015, 09h21

Tu crois que je ne l'ai pas fait ?
J'attendais juste d'avoir quelques pistes en rapport avec mes équations ... qui ne sont pas y' = y ...

**marco056** · 10/11/2015, 10h15

C'est sûr : si tu attends que quelqu'un fasse tes devoirs, tu risques d'attendre longtemps.
Propose un bout de code, montre un peu de bonne volonté et tu auras sans doute de l'aide ici.

**Andromedae** · 10/11/2015, 10h50

Je comprends la discorde, mon code ne s'était pas collé en entier ..
Le code est général, je n'ai pas encore défini y0 et t0 car j'ai un peu de mal à les trouver sachant que ma fonction ne dépend pas de tau :

y'(tau) = constante

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
# usr/bin/python
 
 
from pylab import *
import numpy as np
 
 
 
##################################################
 
# ORBITES AUTOUR D'UN TROU NOIR DE SCHWARZSCHILD #
 
##################################################
 
 
 
 
# ------------------------- #
# Methode d'Euler explicite #
# ------------------------- #
 
Rs = 2                        # Rayon de Swharzschild
m = 1                        # Masse quelconque d'une particule
c = 1                        # Vitesse lumiere
r = 0.8                        # Rayon r
E = 1                        # Energie particule
 
 
F = ( (1-(Rs/r))**(-1) ) * E/(m*c**2)        # Fonction (3a) du sujet
 
def Euler (F, t0, y0, tf, N) :
 
     h = (tf - t0)/N                # Definition du pas fixe
 
     Ly, Lt = [y0],[t0]            # Stockage du point initial
 
     for k in range (N) :
 
         y0 = y0 + h*np.array(F(t0,y0))
         t0+ = h
         Ly.append(y0)
         Lt.append(t0)
 
     return np.array(Lt), np.array(Ly)

**marco056** · 10/11/2015, 12h52

Je comprends mieux aussi

**Andromedae** · 10/11/2015, 13h00

Je pense que :

y(0) = (1- (Rs/r))^-1 * E/mc² avec r = r^n donc n=1

y(tau) = (1- (Rs/r^n))^-1 * E/mc²

**marco056** · 10/11/2015, 19h47

Je ne comprends pas trop : ton pb est résolu ou non (il est marqué comme oui) ?

**tyrtamos** · 11/11/2015, 06h18

Bonjour,

Sans discuter de la validité de l'algorithme, il y a des choses bizarres dans le code Python donné:

- Si le F transmis comme argument dans la fonction Euler est le même que le F du début, il y a un problème: F est écrit comme une variable, mais présenté et utilisé comme une fonction: F(t0,y0)

- je ne comprends pas l'intérêt du np.array de la ligne y0 = y0 + h*np.array(F(t0,y0))

- dans la ligne "t0+ = h", il ne doit pas y avoir d'espace entre le "+" et le "=" (erreur de syntaxe)

**Andromedae** · 11/11/2015, 17h31

Erreur de manipulation : sujet non résolu

Si le F transmis comme argument dans la fonction Euler est le même que le F du début, il y a un problème: F est écrit comme une variable, mais présenté et utilisé comme une fonction: F(t0,y0)

F correspond à mon équation différentielle. Il faudrait donc que je supprime de l'argument Euler alors ?

je ne comprends pas l'intérêt du np.array de la ligne y0 = y0 + h*np.array(F(t0,y0))

Je pensais lister les valeurs dans un tableau pour après en faire un graphique (avec matplotlib)

dans la ligne "t0+ = h", il ne doit pas y avoir d'espace entre le "+" et le "=" (erreur de syntaxe)

Autant pour moi, erreur de syntaxe en effet

Le code donnerait cela alors avec déjà ces quelques modifications :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
# usr/bin/python
 
 
from pylab import *
import numpy as np
 
 
 
##################################################
 
# ORBITES AUTOUR D'UN TROU NOIR DE SCHWARZSCHILD #
 
##################################################
 
 
 
 
# ------------------------- #
# Methode d'Euler explicite #
# ------------------------- #
 
Rs = 2                        # Rayon de Swharzschild
m = 1                        # Masse quelconque d'une particule
c = 1                        # Vitesse lumiere
r = 0.8                        # Rayon r
E = 1                        # Energie particule
 
 
F = ( (1-(Rs/r))**(-1) ) * E/(m*c**2)        # Fonction (3a) du sujet
 
def Euler (t0, y0, tf, N) :
 
     h = (tf - t0)/N                # Definition du pas fixe
 
     Ly, Lt = [y0],[t0]            # Stockage du point initial
 
     for k in range (N) :
 
         y0 = y0 + h*np.array(F(t0,y0))
         t0 += h
         Ly.append(y0)
         Lt.append(t0)
 
     return np.array(Lt), np.array(Ly)

**marco056** · 11/11/2015, 17h53

J'ai un gros doute : r ne varie pas ?

**Andromedae** · 11/11/2015, 17h55

Ce que je pense, c'est que pour y0 ==> r^n = r donc n = 1
Mais à chaque itération, n varie donc r varie

Mais je ne vois pas comment introduire cette fonction r(tau)^n

La méthode d'Euler devant donner je pense : Fonction (r^(n+h)) - Fonction(r^(n)) / h

**marco056** · 11/11/2015, 19h22

Le souci, à mon avis, est en physique, avant d'en effectuer la programmation.
Que cherchez-vous ? A déterminer r ?

**Andromedae** · 20/11/2015, 15h09

Bonjour,

Je reviens vers vous car j'ai continué à chercher, et peut-être à légèrement progresser.

Voici le code amélioré par rapport à la dernière fois :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
#!usr/bin/python 
 
from pylab import *
from scipy import *
from scipy .integrate import odeint		# Pour resoudre les equations differentielles
import numpy as np
import math as m
import matplotlib .pyplot as plt		# Pour tracer des graphiques
 
 
										##################################################
 
										# ORBITES AUTOUR D'UN TROU NOIR DE SCHWARZSCHILD #
 
										##################################################
 
 
 
# Definition des variables :
#---------------------------
 
Rs = 1.0						# Rayon de Swharzschild
r = 20.0*Rs						# Rayon r
E = 10.0						# Energie particule = mc**2
L = 0.0							# Moment cinetique de la particule
 
# Implicitement on a :
m = 1.0							# Masse quelconque d'une particule
c = 1.0							# Vitesse lumiere
 
 
V0 = 0.0
 
 
# Definition du systeme d'equations :
#------------------------------------
 
t_point = 1/(1-(Rs/r)) * E						# Equation differentielle (3a) du sujet / t_point = dt/dtau
 
phi_point = L/(m*r**2)							# Equation differentielle (3b) du sujet / phi_point = dphi/dtau
 
V0**2 == E**2 - (1-(Rs/r))*(1+(L**2/r**2))				# Equation differentielle (3c) du sujet
 
A0 = (Rs/(2*r**2)) + (3/2)*((Rs*L**2 - r*L**2)/(r**4))			# Equation differentielle (3c) du sujet
 
 
print ' '
print 't_point vaut :' ,t_point
print 'phi_point vaut :' ,phi_point
print 'V0**2 vaut :' ,V0**2
print 'A0 vaut :' ,A0
print ' '
 
 
# ----------------------------------------------------------------- # 
# Methode d'Euler explicite pour le mouvement de chute libe radiale #
# ----------------------------------------------------------------- # 
 
t0 = 0.0
tf = 100.0
y0 = 0.0
n = 10000.0
 
def Euler (t_point, phi_point, V0, A0, t0, y0, tf, n) :
 
	h = (tf - t0)/n				# Definition du pas
 
	t = np.arange(t0,tf,h)
 
	y = np.zeros(n)
 
	y[0] = y0
 
	for i in range (n-1) :
 
		x[i+1] = x[i] + h
		y[i+1] = y[i] + h*t_point(y[i],x[i])
 
		print x[i+1]	
	return(x,y)

==> Mes fonctions sont bien définies pour le problème que je pose.
==> Dans mon programme :

* t_point <=> dt/dtau
* phi_point <=> dphi/dtau
* Vo² <=> (dr/dtau)²
* Ao <=> d/dtau (Vo²)

J'ai donc un système de 4 équations du premier ordre

==> Par ailleurs, dans mon programme, m = c = 1 Ce qui permet d'alléger mes équations.

==> Le code fonctionne jusqu'au print de mes fonctions, mais j'ai de gros soucis pour mettre en place ma méthode d'Euler explicite.
==> La quantité à faire varier étant "tau", je bloque sur Euler.

Quelqu'un pourrait-il m'aiguiller un peu ?

Merci d'avance !

**Clodion** · 20/11/2015, 15h58

Envoyé par Andromedae

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
#!usr/bin/python 
 
from pylab import *
from scipy import *
from scipy .integrate import odeint		# Pour resoudre les equations differentielles
import numpy as np
import math as m
import matplotlib .pyplot as plt		# Pour tracer des graphiques
 
 
										##################################################
 
										# ORBITES AUTOUR D'UN TROU NOIR DE SCHWARZSCHILD #
 
										##################################################
 
 
 
# Definition des variables :
#---------------------------
 
Rs = 1.0						# Rayon de Swharzschild
r = 20.0*Rs						# Rayon r
E = 10.0						# Energie particule = mc**2
L = 0.0							# Moment cinetique de la particule
 
# Implicitement on a :
m = 1.0							# Masse quelconque d'une particule
c = 1.0							# Vitesse lumiere
 
 
V0 = 0.0
 
 
# Definition du systeme d'equations :
#------------------------------------
 
t_point = 1/(1-(Rs/r)) * E						# Equation differentielle (3a) du sujet / t_point = dt/dtau
 
phi_point = L/(m*r**2)							# Equation differentielle (3b) du sujet / phi_point = dphi/dtau
 
V0**2 == E**2 - (1-(Rs/r))*(1+(L**2/r**2))				# Equation differentielle (3c) du sujet
 
A0 = (Rs/(2*r**2)) + (3/2)*((Rs*L**2 - r*L**2)/(r**4))			# Equation differentielle (3c) du sujet
 
 
print ' '
print 't_point vaut :' ,t_point
print 'phi_point vaut :' ,phi_point
print 'V0**2 vaut :' ,V0**2
print 'A0 vaut :' ,A0
print ' '
 
 
# ----------------------------------------------------------------- # 
# Methode d'Euler explicite pour le mouvement de chute libe radiale #
# ----------------------------------------------------------------- # 
 
t0 = 0.0
tf = 100.0
y0 = 0.0
n = 10000.0
 
def Euler (t_point, phi_point, V0, A0, t0, y0, tf, n) :
 
	h = (tf - t0)/n				# Definition du pas
 
	t = np.arange(t0,tf,h)
 
	y = np.zeros(n)
 
	y[0] = y0
 
	for i in range (n-1) :
 
		x[i+1] = x[i] + h
		y[i+1] = y[i] + h*t_point(y[i],x[i])
 
		print x[i+1]	
	return(x,y)

==> Mes fonctions sont bien définies pour le problème que je pose.
==> Dans mon programme :

* t_point <=> dt/dtau
* phi_point <=> dphi/dtau
* Vo² <=> (dr/dtau)²
* Ao <=> d/dtau (Vo²)

J'ai donc un système de 4 équations du premier ordre

==> Par ailleurs, dans mon programme, m = c = 1 Ce qui permet d'alléger mes équations.

==> Le code fonctionne jusqu'au print de mes fonctions, mais j'ai de gros soucis pour mettre en place ma méthode d'Euler explicite.
==> La quantité à faire varier étant "tau", je bloque sur Euler.

Quelqu'un pourrait-il m'aiguiller un peu ?

Merci d'avance !

Bonjour,
Je ne comprends pas grand chose, mais:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

V0**2 == E**2 - (1-(Rs/r))*(1+(L**2/r**2))

Signifie tester l'égalité entre "0.0**2" et le reste. C'est le but recherché?
Cela met une erreur en mettant "=" comme pour les autres? Alors on essaye des trucs…
D'autre part:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
t_point = 1/(1-(Rs/r)) * E
phi_point = L/(m*r**2)…

ne sont pas des équations différentielles, mais des constantes.

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
		x[i+1] = x[i] + h
		y[i+1] = y[i] + h*t_point(y[i],x[i])

D'où vient la liste "x"?
Ainsi que dit plus haut, "t_point" est une constante, pas une liste.

L'usage du shell Python est fortement conseillé!

Aiguillage: avant d'envisager pylab, scipy, numpy, matplotlib… il serait approprié de commencer par un tutoriel Python, tel que celui de G.Swinnen. Gratuit.

Clodion

PS:
-il n'y a qu'une fonction de définie et elle n'est pas appelée.
-non, il n'y a pas quatre équations de premier ordre. Il y a quatre constantes.
-le code fonctionne (vraiment?) pour les imports et les déclarations de constantes (et encore!).

Méthode Euler Explicite

Calcul scientifique Python

Discussions similaires

Partager

Partager