Skip to content
Snippets Groups Projects
Commit caaae495 authored by gergaud's avatar gergaud
Browse files

ajout TP EDP

parent b5033968
Branches
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
EDP/TP_chaleur_python_eleves_2021_2022.ipynb 0 → 100644
+ 584
− 0
View file @ caaae495
%% Cell type:markdown id: tags:
<center>
<h1> Résolution de l'équation de la chaleur par discrétisation différences finies en temps et en espace</h1>
<h1> Année 2021-2022 - IENM2 </h1>
<h1> Nom: </h1>
<h1> Prénom: </h1>
</center>
%% Cell type:markdown id: tags:
## Objectifs
%% Cell type:markdown id: tags:
<div style="background:MistyRose">
<br>
Nous souhaitons résoudre le problème de la chaleur muni de conditions limites de type Dirichlet homogènes suivant:
<br>
$$ (P)~ \left \{
\begin{array}{11}
\ \displaystyle \frac{\partial u}{\partial t}-\Delta u = f \, \, \, \, \, \, (\Omega \times [0, T])\\
\ u(0,t)=0 \, \, \, \, \, \, (\partial \Omega) \\
\ u(1,t)=0 \, \, \, \, \, \, (\partial \Omega) \\
\ u(x,0)=u_{0}(x) \, \, \, \, \, \, (t=0)
\end{array}
\right. $$
<br>
sur un domaine monodimensionnel $\Omega$ par une méthode de discrétisation de type différences finies. Nous vous proposons de détailler pas à pas les différentes étapes en réalisant
systématiquement des tests unitaires afin de valider vos implantations.<br>
Les étapes proposées sont les suivantes:<br>
- Choisir une condition initiale $ u_{0}(x)$ au problème $(P)$.
- Imposer une solution exacte $ u_{ex}(x,t)$ au problème $(P)$.
- En déduire $f$ le terme source associé à $ u_{ex}(x,t)$.
- Résoudre numériquement $(P)$ en utilisant les schémas d'Euler explicite et implicite.
- Déterminer numériquement l'ordre de convergence des schémas de discrétisation en temps.
%% Cell type:markdown id: tags:
### Import des librairies Python nécessaires
%% Cell type:code id: tags:
``` python
import scipy as sc
import scipy.sparse as sparse
import scipy.sparse.linalg
import numpy as np
import matplotlib.pyplot as plt
#import matplotlib.animation as animation
import math
```
%% Cell type:markdown id: tags:
### Définition des paramètres du problème (P)
%% Cell type:code id: tags:
``` python
# Nombre de pas de discretisation en temps
N_t =
# Nombre de pas de discretisation en espace
M =
# Nombre de points de discretisation en espace
N = M + 1
# Temps final [s]
T =
# Pas de discretisation temporel sur [0 T]
delta_t = T/N_t
# Pas de discretisation spatial sur le domaine [0,1]
h = 1/(N-1)
# Nombre CFL
c =
```
%% Output
File "<ipython-input-3-e30d73daff6c>", line 2
N_t =
^
SyntaxError: invalid syntax
%% Cell type:markdown id: tags:
### Définition de la fonction solution exacte choisie, du second membre de l'équation de la chaleur et déduction du schéma de construction associé à la solution exacte
%% Cell type:markdown id: tags:
<div style="background:LightGrey">
Indiquer ici le choix de la solution exacte. En déduire le second membre de l'équation définissant le problème (P).
</div>
%% Cell type:code id: tags:
``` python
def solution_exacte(x,t):
"""
Fonction representant la solution exacte du probleme de la chaleur.
Cette fonction doit verifier les conditions initiales et limites sur
la frontiere du domaine.
Entrees:
x: abscisse du point de maillage [float]
t: temps discret considere [float]
Sortie:
Valeur de la solution exacte evaluee en x et au temps t [float]
"""
# A COMPLETER
def second_membre(x, t):
"""
Fonction representant le second membre associe a la solution exacte du probleme
de la chaleur.
Entrees:
x: abscisse du point de maillage [float]
t: temps discret considere [float]
Sortie:
Valeur du second membre evalue en x et au temps t [float]
"""
# A COMPLETER
def schema_exact(N, h, c, T, delta_t, u0):
"""
Fonction representant le tableau bidimensionnel [espace, temps] correspondant a la
solution en chaque point du maillage et a tout temps discret.
Entrees:
N: Nombre de points de discretisation en espace [int]
h: Pas de maillage spatial [float]
c: Nombre CFL [float]
T: Temps final [float]
delta_t: Pas de discretisation temporel [float]
u0: condition initiale [array]
Sortie:
U_exacte: Tableau bidimensionnel de la solution evaluee en tout x et a tout temps t [float]
"""
# Initialisation de la solution
# Prise en compte de la solution initiale
# Deduction de la solution aux autres instants
```
%% Cell type:markdown id: tags:
### Définition de la condition initiale et affichage de la solution exacte pour différents instants
%% Cell type:code id: tags:
``` python
#
# Definition de la condition initiale comme un tableau monodimensionnel en espace
#
u0 = np.zeros(N, dtype=float)
# A COMPLETER
#
# Deduction de la solution exacte par appel de la fonction schema_exact
#
U = schema_exact(N, h, c, T, delta_t, u0)
#
# Representation graphique a differents instants a choisir
#
x = np.linspace(0, 1, N)
plt.plot(x, U[:,0], label="à T=0")
plt.plot(x, U[:, int(N_t/5)], label="à T/5")
plt.plot(x, U[:, int(N_t/2)], label="à T/2")
plt.plot(x, U[:,int(T/delta_t)-1], label="à T" )
plt.legend()
plt.grid(True)
plt.title("Evolution de la température sur la barre pour T=" + str(T))
plt.show()
```
%% Cell type:markdown id: tags:
### Construction du schéma explicite en temps
%% Cell type:code id: tags:
``` python
def schema_explicite(N, h, c, T, delta_t, U_explicite):
"""
Fonction representant le tableau bidimensionnel [espace, temps] correspondant a la
solution en chaque point du maillage et a tout temps discret obtenue avec le schema
explicite.
Entrees:
N: Nombre de points de discretisation en espace [int]
h: Pas de maillage spatial [float]
c: Nombre CFL [float]
T: Temps final [float]
delta_t: Pas de discretisation temporel [float]
u0: condition initiale [array]
Sortie:
U_explicite: Tableau bidimensionnel de la solution du schema explicite evaluee en tout x et a tout temps t [float]
"""
#
# Création de la matrice A via une matrice creuse (interieur du domaine)
#
#
# Application du schéma sur les points interieurs
#
```
%% Cell type:markdown id: tags:
### Affichage de la solution discrète pour le schéma explicite pour différents instants
%% Cell type:code id: tags:
``` python
#
# Definition de la condition initiale comme un tableau monodimensionnel en espace
#
u0 = np.zeros(N, dtype=float)
# A COMPLETER
#
# Deduction de la solution du schema explicite en temps
#
U_explicite = np.zeros((N, int(T/delta_t)), dtype=float)
U_explicite[:,0] = u0
Ue = schema_explicite(N,h,c, T, delta_t, U_explicite)
#
# Representation graphique a differents instants a choisir
#
x = np.linspace(0, 1, N)
plt.plot(x, Ue[:,0], label="à T=0")
plt.plot(x, Ue[:, int(N_t/5)], label="à T/5")
plt.plot(x, Ue[:, int(N_t/2)], label="à T/2")
plt.plot(x, Ue[:,int(T/delta_t)-1], label="à T" )
plt.legend()
plt.grid(True)
plt.title("Evolution de la température sur la barre pour T=" + str(T))
plt.show()
```
%% Cell type:markdown id: tags:
### Construction du schéma implicite en temps
%% Cell type:code id: tags:
``` python
def schema_implicite(N,h,c, T, delta_t, U_implicite):
"""
Fonction representant le tableau bidimensionnel [espace, temps] correspondant a la
solution en chaque point du maillage et a tout temps discret obtenue avec le schema
explicite.
Entrees:
N: Nombre de points de discretisation en espace [int]
h: Pas de maillage spatial [float]
c: Nombre CFL [float]
T: Temps final [float]
delta_t: Pas de discretisation temporel [float]
u0: condition initiale [array]
Sortie:
U_implicite: Tableau bidimensionnel de la solution du schema implicite evaluee en tout x et a tout temps t [float]
"""
#
# Création de la matrice A via une matrice creuse
#
#
# Application du schéma sur les points interieurs
#
```
%% Cell type:markdown id: tags:
### Affichage de la solution discrète pour le schéma implicite pour différents instants
%% Cell type:code id: tags:
``` python
#
# Definition de la condition initiale comme un tableau monodimensionnel en espace
#
u0 = np.zeros(N, dtype=float)
# A COMPLETER
#
# Deduction de la solution du schema implicite en temps
#
U_implicite = np.zeros((N, int(T/delta_t)), dtype=float)
U_implicite[:,0] = u0
Ui = schema_implicite(N,h,c, T, delta_t, U_implicite)
#
# Representation graphique a differents instants a choisir
#
x = np.linspace(0, 1, N)
plt.plot(x, Ui[:,0], label="à T=0")
plt.plot(x, Ui[:, int(N_t/5)], label="à T/5")
plt.plot(x, Ui[:, int(N_t/2)], label="à T/2")
plt.plot(x, Ui[:,int(T/delta_t)-1], label="à T" )
plt.legend()
plt.grid(True)
plt.title("Evolution de la température sur la barre pour T=" + str(T))
plt.show()
```
%% Cell type:markdown id: tags:
### Détermination de l'ordre de convergence des schémas de discrétisation en temps
%% Cell type:markdown id: tags:
<div style="background:LightGrey">
Indiquer ici comment vous calculez cet ordre de convergence. Quelle est la valeur attendue pour les schémas d'Euler explicite ou implicite ? Vérifier ce résultat graphiquement.
</div>
%% Cell type:code id: tags:
``` python
def erreur(U_exacte, U, h):
"""
Calcul de l'erreur discrete
Entrees:
U_exacte: solution exacte du probleme considere [float]
U : solution discrete pour un schema donne [float]
h : pas de discretisation spatial [float]
Sortie:
norme de l'erreur discrete [float]
"""
# A COMPLETER
```
%% Cell type:code id: tags:
``` python
def ordre_discretisation_temps(h, N, T):
"""
Calcul de l'ordre de discretisation en temps
Entrées :
h: pas de discretisation spatial [float]
N: nombre de points de discretisation en espace [int]
T: temps final [float]
Sortie :
time_step: pas de temps en loi log [float, array]
erreur : erreur en loi log[float, array]
"""
time_step = np.zeros(10)
error = np.zeros(10)
count = 0
for loop in range(20, 120, 10):
delta = 1./loop
c = delta/(h**2)
#
# Initialisation pour la solution exacte
#
u0 = np.zeros(N, dtype=float)
# A COMPLETER
U_exacte = schema_exact(N, h, c, T, delta, u0)
#
# Initialisation pour la solution numérique
#
U_implicite = np.zeros((N, int(T/delta)), dtype=float)
u0 = np.zeros(N, dtype=float)
# A COMPLETER
U_implicite[:,0] = u0
#
#
#
U_implicite = schema_implicite(N, h, c, T, delta, U_implicite)
error[count] = math.log(erreur(U_exacte, U_implicite, h))
time_step[count] = math.log(delta)
count += 1
return time_step, error
```
%% Cell type:code id: tags:
``` python
#
# Construction des elements pour l'analyse du schema
#
time_step, error = ordre_discretisation_temps(h, N, T)
#
# Creation de la figure
#
theta = 1.0
fig = plt.figure()
plt.plot(time_step,error,marker='*',color='blue',label='Euler implicite (theta='"%3.2f"%theta+')')
plt.plot(time_step,1.*time_step,color='red',label='Ordre 1 en temps')
plt.plot(time_step,2.*time_step,color='green',label='Ordre 2 en temps')
plt.xlabel("Log(Delta t)")
plt.ylabel("Log(Norme de l'erreur)")
plt.title('Analyse de convergence: Equation de la chaleur')
plt.legend()
plt.grid(True)
plt.show()
```
%% Cell type:markdown id: tags:
<div style="background:LightGrey">
Indiquer ici l'ordre de convergence obtenu effectivement. Obtenez-vous le comportement obtenu ?
</div>
%% Cell type:code id: tags:
``` python
#
# Affichage des donnees pour validation et inter-comparaison
#
```
%% Cell type:markdown id: tags:
### Schéma de discrétisation en temps du second ordre
%% Cell type:markdown id: tags:
<div style="background:LightGrey">
Résoudre (P) cette fois avec le schéma temporel d'ordre deux implicite Crank-Nicolson. Commenter votre code pour indiquer les changements effectués. Vérifier par une méthodologie identique l'ordre de convergence de ce schéma.
</div>
%% Cell type:markdown id: tags:
### Bonus: Généralisation pour des conditions limites spatiales non homogènes
%% Cell type:markdown id: tags:
<div style="background:MistyRose">
<br>
Nous souhaitons résoudre désormais le problème de la chaleur muni de conditions limites de
type Dirichlet non homogènes suivant:
<br>
$$ (P_{nh}) ~ \left \{
\begin{array}{11}
\ \displaystyle \frac{\partial u}{\partial t}-\Delta u = f \, \, \, \, \, \, (\Omega \times [0, T]),\\
\ u(0,t)= u_g{(t)} \, \, \, \, \, \, (\partial \Omega), \\
\ u(1,t)= u_d{(t)} \, \, \, \, \, \, (\partial \Omega), \\
\ u(x,0)= u_{0}(x) \, \, \, \, \, \, (t=0),
\end{array}
\right. $$
<br>
sur un domaine monodimensionnel $\Omega$ par une méthode de discrétisation de type différences finies.<br>
Reprendre les étapes proposées précédemment:
- Imposer une condition initiale et une solution exacte $ u_{ex}(x,t)$ au problème $(P_{nh})$.
- En déduire $f$ le terme source associé à $ u_{ex}$.
- Résoudre numériquement $(P_{nh})$ en utilisant un schéma de votre choix.
%% Cell type:code id: tags:
``` python
```
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment