Table des matières

Équation de Poisson : module Python

1. Introduction

Ce document présente le module poisson, qui permet d'effectuer la résolution numérique de l'équation de Poisson 2D (applications en électromagnétisme et en thermodynamique) par la méthode des différences finies :

où u(x,y) est la fonction inconnue et s(x,y) la fonction source, éventuellement nulle (équation de Laplace). L'équation ci-dessus s'applique à un problème à géométrie XY.

Le module permet également de résoudre l'équation de Poisson en géométrie axiale :

avec u(z,r) la fonction inconnue en coordonnées cylindriques, et l'équation de Poisson vectorielle en géométrie axiale :

L'équation de Poisson en coordonnées polaires :

est en cours d'implémentation.

Deux méthodes itératives de résolution sont possibles :

Ces méthodes itératives sont sous forme de fonctions exécutables programmées en C.

Pour la méthode de Gauss-Seidel avec sur-relaxation, le module comporte un programme tournant sur le processeur central (CPU), et un programme utilisant les capacités de calculs parallèles des cartes graphiques (GPU), via l'interface de programmation OpenCL.

2. Installation du module

Les binaires sous forme de paquets wheel sont aussi disponibles sur pypi.python.org/pypi/poisson/0.2. L'installation se fait avec pip. Pour cela, ouvrir une console puis aller dans le dossier Scripts de la distribution, par exemple C:\Anaconda3\Scripts, et exécuter :

pip install poisson
              

La bibliothèque pthreads est utilisée pour le multi-threading sur le processeur central (CPU). Cette bibliothèque est disponible en standard sur les distributions Linux. Pour la compilation sur Windows, il faut installer la version pthreads-win32. La DLL est incluse dans la distribution binaire pour win32.

Pour compiler et installer sous linux : sudo python setup.py install. Pour compiler sans OpenCL, ouvrir le fichier setup.py et enlever opencl_include et define_macros.

3. Documentation

Exemples d'utilisation :

Scripts de test :

Méthodes numériques

4. Description de l'interface

4.a. Définition du maillage

On définit un maillage réduit, qui servira à définir les objets, et un nombre de niveaux de grilles.

La classe principale est Poisson. Elle se trouve dans le module poisson.main. Son constructeur effectue l'allocation de l'espace mémoire pour le maillage :

Le nombre de niveaux de grilles est le nombre maximal de grilles utilisables dans la méthode multigrille. On aura donc intérêt à définir le problème sur une petite grille et à fair varier le nombre de niveaux de grilles en fonction de la finesse souhaitée pour le calcul.

La figure suivante montre un exemple de maillage avec pw_min=ph_min=2 et levels=3. Il y a trois grilles disponibles pour l'algorithme multigrilles. Les objets seront définis sur la grille 4x4.

../../../../figures/numerique/elliptique/pypoisson/maillage.svg
Figure pleine page

La position du point origine ne joue aucun rôle dans le calcul numérique. Elle est utilisée pour la présentation graphique des résultats. Pour le calcul numérique, seul importe le rapport de la largeur sur la hauteur. Par exemple, si la largeur est deux fois plus grande que la hauteur, et si le nombre de mailles est identique sur les deux dimensions, alors le pas d'espace est deux fois plus grand dans la direction X. En revanche, ces dimensions interviennent dans le calcul des dérivées.

Une seule instance de cette classe est utilisable, car le programme C sous-jacent ne comporte qu'une instance des données.

Pour un problème en coordonnées cylindriques et à symétrie axiale, on utilise la classe PoissonAxial.

La première dimension du maillage (axe Ox en coordonnées cartésiennes) correspond à l'axe de symétrie Oz. Les indices des points sur cet axe sont donc (i,0).

Les fonctions de la classe Poisson décrites ci-dessous sont disponibles aussi pour la classe PoissonAxial.

Pour l'équation de Poisson vectorielle en coordonnées cylindriques et à symétrie axiale, on utilise la classe PoissonAxialVect.

4.d. Polygone

On peut définir un polygone dont les côtés sont parallèles aux axes x et y (ou z et r en géométrie axiale).

La figure suivante montre le polygone obtenu par dirichlet_polygon(1,2,[2,0],[0,1],1).

../../../../figures/numerique/elliptique/pypoisson/polygone.svg
Figure pleine page

On voit que le polygone peut être ouvert. Pour le fermer, il faut se déplacer jusqu'au premier point. On peut aussi définir un polygone sur les bords du domaine, pour modifier la condition limite ou même pour la définir entièrement. Tous les points des grilles jusqu'à la grille la plus fine (ici 32x32) sont bien sûr affectés par cette condition limite.

Les dérivées par rapport à X n'affectent que les côtés parallèles à Y. En principe, un polygone avec une condition de Neumann doit être fermé : il faut donc se déplacer jusqu'au point initial. Les angles recoivent un traitement spécial, qui combine les dérivées selon X et selon Y. Il est important de savoir où se trouve le domaine de calcul par rapport au polygone. Celui-ci est défini par rapport au sens du parcours : si le sens de parcours est le sens horaire alors le domaine se trouve à l'extérieur du polygone.

Ce polygone peut être ouvert. Lorsqu'on parcourt un côté, le milieu de conductivité a1 se trouve à gauche.

4.e. Itérations de Gauss-Seidel sur CPU

Le système linéaire obtenu par discrétisation des équations (équation de Poisson et conditions limites) est résolu par la méthode d'itération de Gauss-Seidel sur le maillage le plus fin. Ce paragraphe présente les fonctions de la classe Poisson qui permettent d'effectuer ce calcul sur le processeur central (CPU).

Pour contrôler la convergence, la fonctions suivante effectue des blocs d'itérations et calcule la norme de la matrice U à chaque bloc.

Sur les processeurs comportant au moins 4 cœurs, l'exécution avec 4 threads est plus rapide.

La fonction suivante fournit le paramètre de sur-relaxation optimale pour le problème de Poisson avec des conditions limites de Dirichlet sur les bords.

Pour les autres types de conditions limites, on peut partir de cette valeur optimale mais il faut en général faire quelques essais.

4.f. Itérations de Gauss-Seidel sur GPU

Ces fonctions ne sont disponibles que sur la version avec OpenCL.

Les itérations peuvent être effectuées sur processeur graphique (GPU). Une plateforme openCL doit être installée.

Sur une carte graphique haut de gamme (type jeu), l'exécution sur GPU est beaucoup plus rapide que sur CPU.

La première fonction permet d'afficher les plateformes openCL présentes sur le système, et pour chaque plateforme les périphériques associés. Dans les cas courants (une seule carte graphique), il y aura une seule plateforme avec un seul périphérique. Les plateformes et les périphériques sont numérotés à partir de 0.

Par défaut, la plateforme 0 et le périphérique 0 sont sélectionnés. La fonction suivante permet de sélectionner une plateforme openCL et un périphérique :

Les deux fonctions suivantes effectuent les itérations, et sont analogues aux fonctions définies plus haut pour le CPU :

4.g. Méthodes multigrilles

Les méthodes multigrilles sont implémentées pour le processeur central (CPU).

Au sujet de ces paramètres, voir Méthode itérative multigrilles

Les nombres d'itérations sont les nombres équivalents pour la grille fine, c'est-à-dire les nombres calculés en tenant compte de la réduction du nombre de points des grilles (4,16, etc). On peut donc comparer la courbe obtenue avec celle donnée par la méthode de Gauss-Seidel. Néanmoins, les transferts intergrilles (prolongation et restriction) ne sont pas pris en compte dans le calcul. On peut estimer que la vitesse de convergence donnée par cette courbe est sur-estimée d'environ 10 à 20 pour cents.

Pour le multigrille complet, la liste des normes ne contient que deux valeurs, la valeur initiale et la valeur finale.

Si les paramètres sont bien choisis, une méthode multigrille nécessite un nombre d'itérations proportionnel au nombre de points du maillage, ce qui en fait une méthode optimale.

4.h. Récupération des données

La matrice U contient les valeurs de u(x,y) aux points de la maille. Elle est fournie sous forme d'une matrice numpy par la fonction suivante :

Les dérivées par rapport à x et y sont obtenues avec les deux fonctions suivantes :

En géométrie axiale, les fonction PoissonAxial.get_derivZ, PoissonAxial.get_derivR, PoissonAxialVect.get_derivZ et PoissonAxialVect.get_derivZ sont analogues au précédentes. Dans ce cas, il peut être nécessaire de calculer la dérivée suivante (par exemple pour calculer le rotationnel du potentiel vecteur en magnétostatique) :

La fonction suivante effectue ce calcul :

Pour faciliter le tracé des résultats, il peut être utile de calculer les valeurs des abscisses et ordonnées des points du maillage. Les fonctions suivantes permettent d'obtenir ces listes. Le point origine du domaine est celui spécifié par les arguments x0 et y0 du constructeur.

Creative Commons LicenseTextes et figures sont mis à disposition sous contrat Creative Commons.