8  Conduction thermique stationnaire

8.1 🌡️ Introduction à la conduction thermique

8.1.1 🎯 Objectif du notebook

Ce notebook a pour but de comprendre comment résoudre un problème de conduction thermique stationnaire à l’aide de la méthode des éléments finis (FEM) avec FEniCS.

Nous allons :

  • introduire le cadre théorique,
  • résoudre un problème modèle (équation de Poisson),
  • implémenter une résolution numérique complète,
  • comparer avec une solution analytique.

8.1.2 📜 Origines historiques

L’étude de la conduction thermique repose sur les travaux de Joseph Fourier (1822), qui introduit l’équation de la chaleur décrivant la diffusion thermique.

Le cas stationnaire conduit à l’équation de Poisson, fondamentale en physique :

  • thermique
  • électrostatique
  • mécanique des fluides

8.1.3 ⚙️ Hypothèses du modèle

Nous travaillons dans le cadre suivant :

  • régime stationnaire
  • matériau homogène
  • pas de convection ni rayonnement
  • conditions de Dirichlet (température imposée)

8.1.4 🧮 Équation étudiée

\[ -k \Delta T = f \quad \text{dans } \Omega \]

👉 Objectif : déterminer le champ de température \(T(x,y)\).

8.2 🧪 Prise en main de FEniCS : équation de Poisson

Avant de traiter un problème physique complet, on commence par un exemple simple.

8.2.1 🎯 Objectif

  • comprendre la structure d’un code FEniCS
  • identifier les étapes clés d’une résolution FEM

8.2.2 🧩 Étapes de la méthode

  1. Définition du maillage
  2. Définition de l’espace de fonctions
  3. Écriture de la formulation variationnelle
  4. Application des conditions aux limites
  5. Résolution
  6. Visualisation

8.2.3 🔷 Maillage du domaine

On considère un carré unité discrétisé en éléments finis.

8.2.4 🧮 Formulation variationnelle

On passe de la forme forte à une forme intégrale adaptée aux éléments finis.

👉 Idée clé :

  • approximer la solution dans un espace discret
  • transformer l’EDP en problème algébrique

🔑 Ce que fait le code

  • définit l’espace \(V\)
  • introduit fonctions test et inconnue
  • construit :
  • forme bilinéaire \(a(u,v)\)
  • forme linéaire \(L(v)\)
  • impose les conditions de Dirichlet

8.2.5 📏 Évaluation de l’erreur

On compare la solution numérique avec la solution analytique connue.

Deux mesures :

  • norme \(L^2\) → erreur globale
  • erreur max → erreur ponctuelle

👉 Cela permet de vérifier :

  • la validité du modèle
  • la qualité du maillage
# Interpolation de la solution exacte pour comparaison
V2 = fem.functionspace(domain, ("Lagrange", 2)) # Définition d'un nouvel espace de fonctions V2 avec un degré de Lagrange 2
uex = fem.Function(V2) # Définition d'une fonction uex dans l'espace de fonctions V2
uex.interpolate(lambda x: 1 + x[0]**2 + 2 * x[1]**2) # Interpolation de uex avec une fonction lambda quadratique

# Calcul de l'erreur L2
import numpy
L2_error = fem.form(ufl.inner(uh - uex, uh - uex) * ufl.dx) # Calcul de l'erreur L2
error_local = fem.assemble_scalar(L2_error) # Assemblage de l'erreur locale sur chaque processus MPI
error_L2 = numpy.sqrt(domain.comm.allreduce(error_local, op=MPI.SUM)) # Calcul de l'erreur L2 totale en prenant la racine carrée de la somme des erreurs locales

# Calcul de l'erreur maximale
#error_max = numpy.max(numpy.abs(u_D.x.array - uh.x.array)) # Calcul de l'erreur maximale entre la solution exacte et la solution numérique
error_max = np.max(np.abs(uex.x.array - uh.x.array))

# Affichage des erreurs seulement sur un processus (rank) spécifique (rank 0 dans ce cas)
if domain.comm.rank == 0:
    print(f"Error_L2 : {error_L2:.2e}")
    print(f"Error_max : {error_max:.2e}")
Error_L2 : 5.51e+02
Error_max : 3.25e+02

8.2.6 📊 Visualisation 3D

La solution est représentée en déformant le maillage selon :

  • la valeur de \(u\)
  • effet “surface thermique”

👉 Très utile pour :

  • interpréter physiquement les résultats
  • détecter des anomalies

8.2.7 ✅ Validation du modèle

Dans cet exemple (carré), la solution analytique utilisée est :

\[ u(x,y) = 1 + x^2 + 2y^2 \]

👉 Elle permet de :

  • valider le code
  • mesurer précisément l’erreur numérique

8.3 🔥 Problème de conduction : disque thermique

Nous passons maintenant à un cas physique réaliste.

8.3.1 🧾 Énoncé

On considère un disque métallique de rayon \(R=1\) m.

8.3.2 🔧 Paramètres physiques

  • Conductivité thermique : \(k = 0.92 \, W\,m^{-1}\,K^{-1}\)

Conditions :

  • Source interne uniforme : \(f = 100 \, W/m^2\)

  • Température imposée au bord : \(T_d = 25^\circ C\)

8.3.3 🎯 Objectif

Déterminer :

  • le champ de température \(T(x,y)\)
  • la distribution thermique dans le disque

8.3.4 ⚖️ Formulation forte

8.3.4.1 Équation :

\[ \nabla \cdot \mathbf{q} = f \]

8.3.4.2 Loi de Fourier :

\[ \mathbf{q} = -k \nabla T \]

👉 On obtient :

\[ -k \Delta T = f \]

8.3.4.3 Condition limite :

\[ T = T_d \text{ sur } \Gamma \]

👉 Problème de type Dirichlet pur

8.3.5 🧮 Résolution analytique

Nous cherchons ici une solution exacte du problème afin de :

  • valider la solution numérique
  • mieux comprendre le comportement physique

8.3.5.1 🔁 Passage à l’équation de Poisson

En combinant :

  • l’équation d’équilibre
  • la loi de Fourier

on obtient :

\[ -k\,\Delta T(x,y)=f \]

8.3.5.2 🔄 Passage en coordonnées polaires

Compte tenu de la symétrie du disque, le problème ne dépend que de \(r\) :

\[ \frac{1}{r} \frac{\partial}{\partial r} \left(r \frac{\partial T(r)}{\partial r}\right)=-\frac{f}{k} \]

👉 Cela simplifie fortement le problème.

8.3.5.3 ✏️ Intégration de l’équation

Après double intégration :

\[ T(r)=-\dfrac{fr^2}{4k}+A\ln(r)+B \]

8.3.5.4 📌 Détermination des constantes

  • En \(r=0\) : la température doit rester finie ⇒ \(A=0\)
  • En \(r=R\) : \(T(R)=T_d\)

👉 On obtient la solution finale :

\[ T(r)=\dfrac{f}{4k}(R^2-r^2)+T_d \]

💡 À retenir

  • solution parabolique radiale
  • maximum au centre
  • dépend de :
  • la source \(f\)
  • la conductivité \(k\)

8.3.6 🧩 Formulation faible du problème

La formulation faible permet de transformer l’EDP en un problème adapté aux éléments finis.

8.3.6.1 🔁 Principe

On multiplie l’équation par une fonction test \(T^{\ast}\) et on intègre sur \(\Omega\).

8.3.6.2 📐 Espaces fonctionnels

  • Espace admissible :

\[ \Theta = \{T^{\ast} \in H^1(\Omega) \, | \, T^{\ast}=T_d \text{ sur } \Gamma_{ext} \} \]

  • Espace test :

\[ \hat{\Theta} = \{T^{\ast} \in H^1_0(\Omega) \, | \, T^{\ast}=0 \text{ sur } \Gamma_{ext} \} \]

8.3.6.3 🧮 Formulation variationnelle

Trouver \(T \in \Theta\) tel que :

\[ a(T,T^{\ast}) = L(T^{\ast}) \quad \forall T^{\ast} \in \hat{\Theta} \]

avec :

  • forme bilinéaire :

\[ a(T,T^{\ast})=\int_{\Omega}k\,\nabla T \cdot \nabla T^{\ast}\, dV \]

  • forme linéaire :

\[ L(T^{\ast})= \int_{\Omega}f\,T^{\ast}\, dV \]

⚠️ Remarque importante

👉 Cette solution analytique n’existe que pour des cas simples (géométrie régulière).

Dans les cas réels :

  • géométries complexes
  • matériaux hétérogènes

➡️ on utilise la méthode des éléments finis pour obtenir une solution approchée.

8.4 💻 Résolution numérique

Nous allons maintenant :

  1. créer un maillage avec GMSH
  2. résoudre avec FEniCS

8.4.1 🧱 Maillage du disque

Discrétisation en triangles → adaptée aux géométries courbes.

Info    : Meshing 1D...
Info    : Meshing curve 1 (Ellipse)
Info    : Done meshing 1D (Wall 0.000111379s, CPU 0.002472s)
Info    : Meshing 2D...
Info    : Meshing surface 1 (Plane, Frontal-Delaunay)
Info    : Done meshing 2D (Wall 0.106236s, CPU 0.126485s)
Info    : 1550 nodes 3099 elements

8.4.2 ⚙️ Implémentation

Le code réalise :

  • définition de l’espace \(V\)
  • écriture de la formulation
  • application des conditions limites
  • résolution du système linéaire

👉 Résultat : champ de température discret

Valeurs de la condition aux limites: [298. 298. 298. ... 298. 298. 298.]
Solution obtenue: [298.         301.15427399 298.         ... 299.89482035 299.38706021
 299.38414955]

8.4.3 💾 Export des résultats

Les résultats sont exportés au format XDMF pour visualisation sous ParaView.

👉 Permet :

  • visualisation avancée
  • post-traitement

8.5 ✅ Conclusion

Ce notebook a permis de :

  • comprendre la conduction thermique stationnaire
  • manipuler la méthode des éléments finis
  • utiliser FEniCS + GMSH
  • comparer solution analytique vs numérique

8.5.1 🚀 À retenir

  • FEM = outil puissant pour EDP
  • précision dépend du maillage
  • indispensable en ingénierie

8.5.2 🔍 Perspectives

  • raffinement de maillage
  • conditions de Neumann
  • problèmes 3D
  • couplage thermique / mécanique