Table des matières PDFPython

Équation de Poisson à deux dimensions

1. Introduction

Un exemple d'équation de Poisson est celle vérifiée par le potentiel électrostatique :

2Vx2+2Vy2=-ρε

ρ(x,y) est la densité volumique de charge électrique.

Un autre exemple est l'équation de la chaleur en régime stationnaire vérifiée par la température :

2Tx2+2Ty2=-σλ

σ(x,y) est la densité volumique de puissance générée localement.

La forme générale de l'équation de Poisson est :

2ux2+2uy2=s(x,y)

s(x,y) est une fonction connue, appelée la source.

Pour la discrétisation, on utilise aussi la forme intégrale de cette équation, pour un volume V délimitée par une surface S :

Sgrad(u)nds=V sdv

La première intégrale est le flux sortant. Cette équation correspond au théorème de Gauss de l'électrostatique, à un bilan d'énergie pour un problème thermique.

On se limite ici à l'équation de Poisson bidimensionnelle en coordonnées cartésiennes. Dans ce cas, le gradient est dans le plan (x,y) et l'intégrale de surface se ramène à une intégrale sur une courbe fermée C; l'intégrale de volume se ramène à une intégrale sur la surface Σ délimitée par la courbe :

Cgrad(u)ndl=ΣsdΣ

Voir aussi Discrétisation de l'équation de Poisson en coordonnées polaires et Discrétisation de l'équation de Poisson en géométrie axiale.

2. Maillage

Les nœuds du maillage sont définis par :

xi=iΔxyj=jΔy

L'indice i varie de 0 à Nx-1, l'indice j de 0 à Ny-1. On pose :

Ui,j=u(iΔx,jΔy)Si,j=s(iΔx,jΔy) figureA.svgFigure pleine page

3. Discrétisation de l'équation

La méthode des différences finies consiste à remplacer la dérivée seconde par une expression approchée à l'aide des valeurs au nœud et aux nœuds voisins. Pour la dérivée seconde par rapport à x, la différence finie la plus utilisée est :

Ui+1,j-2Ui,j+Ui-1,jΔx2

En procédant de même pour l'autre dérivée, on obtient l'équation suivante pour chaque nœud du maillage :

Δy2(Ui+1,j+Ui-1,j)+Δx2(Ui,j+1+Ui,j-1)-2(Δx2+Δy2)Ui,j=Δx2Δy2Si,j

Si U désigne la matrice colonne formée des N points du maillage, la relation précédente s'écrit sous forme matricielle :

AU=F

Il s'agit donc de résoudre un système linéaire à N inconnues.

On utilisera également la notation suivante :

Di,jUi+1,j+Gi,jUi-1,j+Hi,jUi,j+1+Bi,jUi,j-1+Ci,jUi,j=Fi,j

Les matrices D (droite), G (gauche), H (haut), B (bas), C (centre) et S (source) ont les mêmes dimensions que le maillage.

Pour un point non situé sur une frontière, on a donc :

Di,j=Δy2Gi,j=Δy2Hi,j=Δx2Bi,j=Δx2Ci,j=-2(Δx2+Δy2)Fi,j=Si,jΔx2Δy2

La discrétisation peut aussi être obtenue avec la méthode des volumes finis ([1]), qui consiste à discrétiser la forme intégrale exprimée pour un volume de contrôle. Pour un point situé dans le domaine (point A), le volume de contrôle est représenté sur la figure par le contour (a,b,c,d). Il s'agit d'un parallélépipède de côtés Δx,Δy,Δz. Le flux sortant du gradient sur la surface (ab) est discrétisé de la manière suivante :

(ab)grad(u)ndsUi,j+1-Ui,jΔyΔxΔz

On procède de même pour les surfaces (bc), (cd) et (da). L'intégrale de volume de la source est remplacé par :

Si,jΔxΔyΔz

4. Conditions limites

4.a. Condition de Dirichlet

Dans ce cas, le potentiel électrique, ou la température, est imposée sur une frontière. On écrira simplement pour tout point de la frontière :

Di,j=0Gi,j=0Hi,j=0Bi,j=0Ci,j=1Fi,j=U0

4.b. Condition de Neumann

Dans certains cas, le champ électrique peut être imposé sur une frontière. Par exemple sur une frontière perpendiculaire à l'axe Oy :

Vy=-Ey

Dans un problème de thermique, le flux thermique peut être imposé sur une frontière :

Ty=-φyλ

On écrira donc la condition limite sous la forme générale :

uy=wy

Pour le point B situé sur le bord j=0, on écrit un développement de Taylor à l'ordre 2:

u(xi,yj+1)u(xi,yj)+Δyuy(xi,yj)+Δy222uy2(xi,yj)

La dérivée seconde par rapport à y est remplacée par une dérivée seconde par rapport à x à l'aide de l'équation de Poisson :

u(xi,yj+1)u(xi,yj)+Δyuy(xi,yj)+Δy22-2ux2(xi,yj)+s(xi,yj)

Après avoir discrétisé la dérivée seconde comme plus haut, on obtient finalement pour le point de la frontière :

Di,j=Gi,j=12Δy2Hi,j=Δx2Bi,j=0Ci,j=-Δx2-Δy2Fi,j=wy,iΔx2Δy+Si,j12Δx2Δy2

La même équation est obtenue plus simplement par la méthode des volumes finis, en appliquant la forme intégrale sur la surface (e,f,g,h) :

-wy,iΔx+Ui+1,j-Ui,jΔxΔy2+Ui-1,j-Ui,jΔxΔy2+Ui,j+1-Ui,jΔyΔx=12ΔxΔySi,j

La condition limite de Neumann sur les points C et D s'obtient de manière analogue.

Pour le point F situé sur un coin, il y a une dérivée wy imposée sur le côté parallèle à x, une dérivée wy imposée sur le côté parallèle à y. La forme intégrale sur le volume de contrôle (x,y,z,w) conduit à :

-wx,j12Δy-wy,i12Δx+Ui+1,j-Ui,jΔxΔy2+Ui,j+1-Ui,jΔyΔx2=Si,j14ΔxΔy

On obtient donc pour le point F :

Di,j=12Δy2Gi,j=0Hi,j=12Δx2Bi,j=0Ci,j=-12Δx2-12Δy2Fi,j=wx,j12ΔxΔy2+wy,j12Δx2Δy+Si,j14Δx2Δy2

La condition limite de Neumann sur le coin E s'obtient de la même manière. Les coins E et F sont similaires car le domaine de calcul se trouve à l'intérieur du coin. Il n'en est plus de même pour le coin G, sommet d'un objet situé dans le domaine. Le volume de contrôle pour ce point est (a,b,c,d,e,f). On obtient :

Di,j=12Δy2Gi,j=Δy2Hi,j=12Δx2Bi,j=Δx2Ci,j=-32Δx2-32Δy2Fi,j=-wx,j12ΔxΔy2-12wy,iΔx2Δy+Si,j34Δx2Δy2

4.c. Relation de passage entre deux milieux différents

Sur une frontière entre deux milieux différents, il faut tenir compte des différences de conductivité dans un problème de thermique, ou des différences de constantes diélectriques ou magnétiques dans un problème d'électromagnétisme. Cela revient à écrire la forme intégrale de la manière suivante :

Cagradundl=ΣsdΣ

où le coefficient a représente la propriété relative du milieu (conductivité, permittivité, etc), égale à a1 ou a2. On introduit ici deux constantes pour obtenir un traitement symétrique des deux milieux, mais une des deux constantes sera égale à 1.

En un point sur une frontière parallèle à y, la discrétisation de la forme intégrale s'écrit alors :

a1Ui-1,j-Ui,jΔxΔy+a2Ui+1,j-Ui,jΔxΔy+a1+a22Ui,j-1-Ui,jΔyΔyx+a1+a22Ui,j+1-Ui,jΔyΔx=ΔxΔySi,j

On obtient finalement :

Di,j=a2Δy2Gi,j=a1Δy2Hi,j=a1+a22Δx2Bi,j=a1+a22Δx2Ci,j=-(a1+a2)(Δx2+Δy2)Fi,j=Si,jΔx2Δy2

Dans le cas a1=a2=1 , on retrouve la discrétisation de l'équation de Poisson.

5. Méthode itérative de Gauss-Seidel

Dans une méthode itérative (ou méthode de relaxation), le système AU=F est résolu par itérations successives, jusqu'à convergence vers la solution. Au départ, les valeurs données aux Ui,j sont quelconques (généralement nulles). À chaque itération, les N points du maillage sont parcourus dans un ordre déterminé, par exemple par rangées successives. Dans la méthode de Gauss-Seidel, les valeurs Ui,j sont stockées dans un seul tableau. Le noeud courant est modifié par la relation :

Ui,j*=1Ci,jFi,j-Di,jUi+1,j-Gi,jUi-1,j-Hi,jUi,j+1-Bi,jUi,j-1

Les noeuds du membre de droite sont soit des noeuds modifiés dans la précédente itération, soit des noeuds modifiés dans la même itération. Par exemple, dans le cas d'un balayage par lignes successives de gauche à droite et de bas en haut, les termes G et B ont été actualisés dans la même itération.

L'inconvénient de cette méthode est sa lenteur de convergence. Elle peut être accélérée par la méthode de sur-relaxation, qui consiste à adopter une moyenne pondérée de la valeur calculée par le schéma des différences finies et de l'ancienne :

Ui,jnew=ωUi,j*+(1-ω)Ui,jold

Le paramètre de sur-relaxation doit vérifier :

1<ω<2

Pour un problème donné, il faut chercher le paramètre optimal.

Pour plus de détails voir Méthode de Gauss Seidel .

Pour les très grands maillages, la convergence de la méthode de Gauss-Seidel peut être accélérée de manière spectaculaire par les méthodes multigrilles, qui consistent à utiliser successivement des maillages de différentes tailles lors du processus d'itération. Pour plus de détails, voir Méthode multigrilles.

Références
[1]  R.H. Pletcher, J.C. Tannehill, D.A. Anderson,  Computational Fluid Mechanics and Heat Transfer,  (CRC Press, 2013)
Creative Commons LicenseTextes et figures sont mis à disposition sous contrat Creative Commons.