Bonjour,
Je me permet de vous contacter aujourd'hui car je rencontre des difficultés à résoudre cet exercice. Voici l'énoncé :
On considère une corde tendue vibrante dont la longueur à l'équilibre est 𝐿 et dont les deux extrémités sont fixées : 𝑢(0,𝑡) = 𝑢(𝐿,𝑡) = 0. La masse de la corde par unité de longueur est 𝜌 et la tension à l'équilibre est donnée par 𝑇0. Le déplacement
transversal de la corde par rapport à sa position d'équilibre est donné par 𝑢(𝑥,𝑡) et satisfait le problème aux conditions limites :
{𝜌𝜕²𝑢/𝜕𝑡²(𝑥,𝑡) = 𝑇0𝜕²𝑢/𝜕𝑥²(𝑥,𝑡)
𝑢(0,𝑡) = 𝑢(𝐿,𝑡) = 0 𝑢(𝑥, 0) = 𝑢0(𝑥) 𝜕𝑢/𝜕𝑡 (𝑥, 0) = 0
Nous allons utiliser les paramètres matériels et la condition initiale donnés ci-dessous : 𝐿 = 1 [𝑚]; 𝜌 = 0.01 [𝐾𝑔⁄𝑚]; 𝑇0 = 0.1 [𝑁]
𝑢0(𝑥) = sin (𝜋𝑥/𝐿) [𝑚𝑚]
La solution du problème est parfaitement périodique dans le temps avec une période
𝑇 = 2𝐿/𝑐. La vitesse de propagation de l'onde est donnée par 𝑐 = √𝑇0/𝜌.
1) Pour les discrétisations spatiale et temporelle
a) Effectuez une discrétisation spatiale par les différences finies
b) Cherchez le système de n équations différentielles du second ordre
𝒖𝒊(𝒕) = 𝒖(𝑿𝒊, 𝒕) 𝒊 = 𝟏, … , 𝒏
c) Effectuer la discrétisation temporelle par les différences finies
2) Résoudre ce problème avec la méthode de Runge Kutta d'ordre 2 et d'ordre
4 (RK2 et RK4)
3) Tracez la solution sur 2T (RK2, RK4 et odeint)
Voici ce que j'ai commencé à coder :
1)
a) On me demande simplement de réaliser une discrétisation spatiale donc je passe par linspace.
def discretisation_spatiale(L, N):
return np.linspace(0, L, N+1)
J'ajusterai N par la suite.
b) Pour cela, je définis les différentes dérivées partielles à l'aide de la méthode des différences finies, donc ∂²u/∂x²(x,t) = ([u i+1,j] −2*[u i,j]+[u i−1,j])/dx**2 et ∂²u/∂t²(x,t) = ([u i,j+1] −2*[u i,j] + u[ i,j−1])/dt**2. J'injecte ces expressions dans l'équation vérifiée par l'énoncé pour obtenir le systèmes de n équations demandé :
[u i,j+1] =2*[u i,j]−u[ i,j−1] + α*([u i+1,j] −2*[u i,j] + [u i−1,j]) en posant α= T0(Δt)²/ρ(Δx)²
c) Je procède comme à la question a) mais en passant par la période cette fois-ci. Je définis donc les paramètres dont j'ai besoin :
T = 2*L / np.sqrt(T0/rho)
ht = T/N
def discretisation_temporelle(T, N):
return np.linspace(0, T, N+1)
2) C'est ici que j'ai du mal. Dans un premier temps, je vais me concentrer sur la méthode de Runge Kutta d'ordre 2 uniquement. Voici les fonctions que j'ai codé :
def derivation(u,v):
v_prime = np.zeros_like(v)
v_prime[1:-1] = (T0 / rho) * np.diff(u, 2) / dx**2
return v,v_prime
def runge_kutta_2(u, v, dt):
k1, k1_prime = derivation(u, v)
k2, k2_prime = derivation(u + dt/2 * k1, v + dt/2 * k1_prime)
u_prochain = u + dt * k2
v_prochain = v + dt * k2_prime
return u_prochain, v_prochain
3) Là je définis les paramètres du problème, mon pas de temps et j'initialise les valeurs.
N = 100
dx = L / N
x = np.linspace(0, L, N+1)
T = 2 * L / np.sqrt(T0 / rho)
dt = 0.1
t = np.arange(0, 2*T, dt)
u_init = u0(x)
v_init = np.zeros_like(u_init)
u = u_init
v = v_init
plt.figure(figsize=(10, 5))
for i in range(len(t)):
if i % 50 == 0:
plt.plot(x,u,label=f"t={t[i]:.2f}s")
u,v = runge_kutta_2(u,v,dt)
plt.legend()
plt.grid(True)
plt.show()
Le hic c'est que je n'obtient absolument pas ce qui m'est demandé, à savoir la solution sur 2T. J'ai du mal à comprendre ce que l'on attend de moi (je ne suis absolument pas avancé en informatique) puisque je ne comprends pas vraiment la méthode de résolution de Runge Kutta mise à part que l'on procède par itération en approchant la valeur de la solution.
Pouvez-vous me dire si je suis sur la bonne voie/ ce qu'il faudrait que j'essaye s'il vous plaît ?
Merci d'avance,
~adddraw
Bonjour,
Je veux bien regarder et vous aider, mais dites moi quel logiciel vous utilisez, je vais le faire pour voir où vous êtes en difficulté.
Bonjour,
C'est très gentil de votre part. Je travaille sur Python. Vous pouvez également le faire sur Spyder .
pour la question 1a) ok
pour la question 1b) j'ai du mal à vous relire : que représente j?,j'ai l'impression que vous répondez à la question c)
je ne suis que ce bon, j'écrirai plutôt
pour la discrétisation spatiale on a :
pour le système de discrétisation temporelle
Oui vous avez raison mais j'ai du mal à l'écrire sur l'ordinateur en fait. C'est aussi bête que ça... Ensuite je combine dans l'équation vérifiée et j'obtiens le système de n équations cherché.
Annule et remplace le précédent ( j'ai répondu trop vite et Latex ne m'a pas corrigé automatiquement bien sûr)
pour la question 1a) Ok
pour la question 1b) j'ai du mal à vous relire :
que représente l'indice "j", j'ai l'impression que vous répondez à la question c.
j'écrirai plutôt
pour la question 1b) la discrétisation spatiale,Système de n équations différentielles ordinaires du second ordre on a :
pour la question 1c) le système de discrétisation temporelle,Système de 2n équations différentielles ordinaires du premier ordre
Ah je comprends mieux. Oui il fallait "juste" appliquer la méthode des différences finies, comme demandé. Merci de me l'avoir fait remarquer !
Il y a des choses qui ne vont pas dans le code aussi. Par exemple, c'est une très mauvaise idée d'appeler np.arange avec un dtype non entier, ça peut provoquer des erreurs d'arrondi subtiles ou des divisions par 0, et tu auras des NaN partout sans comprendre pourquoi.
Méfie toi aussi des expressions du type L / N quand les deux sont entiers, car la division de python est un peu piégeuse. Si ton code est par malheur exécuté par Python 2 ou une implémentation de Python non standard, ça peut très bien te renvoyer L // N à la place, c'est à dire la partie entière de L/N. 2/3 deviendrait alors 2//3, qui vaut 0, et là aussi tu t'exposes à des bugs en tout genre.
C'est toujours un peu mieux dans ce cas de poser L = 1.0 et pas L = 1. Ou bien d'écrire L / float(N).
Pour le reste, Runge-Kutta ne sert pas théoriquement à résoudre les EDP. C'est pour les systèmes autonomes d'ordre 1 uniquement. Ici il faut ruser pour se ramener à quelque chose qui y ressemble. Ca s'appelle "Runge-Kutta lines method" en anglais si tu veux te renseigner là dessus
Je n'étais pas du tout au courant pour la division entre deux entiers. Je ferai attention à l'avenir. Je vais creuser cette méthode alors. Merci bien .
j'ai regardé votre programme, voici des éléments à regarder :
x = np.linspace(0, L, N+1) donc len(x)=101 est-ce bien ce que vous voulez?
Que vaut u0(x)? tel que rédigé il vaut 0 partout, ne faudrait-il pas mettre u0=np.sin(np.pi*x/L)?
A l'intérieur de la fonction dérivation vous utilisez dx,mais vous ne le passez pas en paramètre?
Vous devez être membre accéder à ce service...
Pas encore inscrit ?
1 compte par personne, multi-compte interdit !
Ou identifiez-vous :