N.B. : le candidat attachera la plus grande importance à la clarté, à la précision et à la concision de la rédaction. Si un candidat est amené à repérer ce qui peut lui sembler être une erreur d'énoncé, il le signalera sur sa copie et devra poursuivre sa composition en expliquant les raisons des initiatives qu'il a été amené à prendre.
Les calculatrices sont autorisées
Le sujet est composé de deux parties, largement indépendantes.
Autour de l'équation de Poisson
Ce problème s'intéresse à la résolution numérique de quelques problèmes d'électrostatique. Il se compose de deux parties.
I. Étude de l'équation de Poisson et de différentes méthodes de résolution numérique.
II. Deux études de cas : fil infini chargé et mouvement d'un électron entre les plaques d'un condensateur.
Les différentes parties sont largement indépendantes.
Un aide-mémoire numpy/matplotlib/pyplot est présent à la fin du sujet.
Partie I - Équation de Poisson
I. 1 - Établissement de l’équation
Q1. Rappeler l'équation de Maxwell-Gauss ainsi que la relation entre le champ et le potentiel électrostatique . En déduire l'équation de Poisson :
Préciser les noms et les unités usuelles de et .
Q2. Citer plusieurs situations physiques en dehors de l'électrostatique pour lesquelles il existe une équation analogue.
I. 2 - Équation adimensionnée pour un problème plan
On veut résoudre l'équation de Poisson dans une portion de plan carrée de côté . On pose :
Q3. Montrer qu'on peut écrire l'équation sous la forme suivante :
où sera exprimé en fonction de et .
I. 3 - Discrétisation
Afin de résoudre numériquement l'équation de Poisson, on va utiliser un maillage de , de pas , et on va transformer les dérivées partielles par des différences entre les valeurs de aux différents points du maillage (on parle aussi des nœuds du maillage). La figure 1 (page suivante) représente le maillage de pour .
Figure 1 - Maillage de pour
Q4. En faisant un développement limité à l'ordre 2 autour du point de coordonnées ( ), montrer qu'on peut exprimer la valeur de en ce point sous la forme suivante : .
Q5. Comme et , on note désormais le potentiel en un point du maillage. Montrer alors qu'on peut écrire l'équation de Poisson sous la forme suivante :
étant une fonction à définir en fonction de et .
I. 4 - Résolution
La fonction étant connue, on montre en mathématiques que la solution de l'équation de Poisson est unique si on fixe les conditions aux limites sur la frontière du domaine . Ces conditions sont essentiellement de deux types:
on impose le potentiel en tout point de (conditions de Dirichlet),
on impose une condition sur les dérivées partielles de en tout point de (conditions de Neumann).
Dans ce problème, on ne va considérer que des conditions de Dirichlet.
La frontière contient naturellement les points du bord de (donc appartenant aux quatre côtés du carré), mais elle peut aussi contenir certains points à l'intérieur de où le potentiel est fixé en raison de la présence d'électrodes.
L'ensemble des points de coordonnées est donc composé de deux sous-ensembles :
ceux dont le potentiel est connu, appartenant à la frontière ,
ceux dont le potentiel est inconnu, appartenant à mais pas à (donc dans ).
Méthode de Jacobi
À partir de l'équation (1), on peut exprimer :
La résolution s'effectue alors en deux étapes.
Initialisation
a) On fixe le potentiel des points de à la valeur imposée physiquement (bords et électrodes).
b) On donne aux points de potentiel inconnu, donc appartenant à , une valeur arbitraire , en général nulle.
Itérations
On calcule une nouvelle valeur des potentiels en appliquant l'équation (2) pour tous les points de , tandis que pour les points de .
Le processus est répété jusqu'à obtenir des valeurs du potentiel quasiment stables. En notant le nombre d'itérations, on a donc pour le point de coordonnées n'appartenant pas à la frontière :
La convergence de la méthode est vérifiée à l'aide du critère de convergence , défini par :
Le calcul sera stoppé au bout de itérations, quand deviendra inférieur à un seuil de convergence fixé arbitrairement.
Implémentation informatique
On va utiliser la bibliothèque numpy permettant une utilisation simple des tableaux de flottants à deux dimensions ; un aide-mémoire est disponible en fin de sujet.
Le chargement des bibliothèques classiques est assuré par les lignes suivantes :
# importation des bibliothèques
import numpy as np
import matplotlib.pyplot as plt
import math
On supposera que les tableaux numpy suivants, utilisés comme arguments dans les fonctions à définir dans les questions qui suivent, ont pour signification :
V : tableau courant du potentiel en un point de ,
rhos : tableau contenant la densité de charge en un point de ,
frontiere : tableau de booléens indiquant si le point de coordonnées ( ) appartient ou non à . En particulier, tous les points du bord du domaine seront tels que frontiere True.
Q6. Écrire la fonction nouveau_potentiel(V, rhos, frontiere, i, j) retournant la nouvelle valeur du potentiel au point selon l'équation (3).
Q7. Montrer que pour modifier toutes les valeurs contenues dans pendant une itération, il est nécessaire de disposer d'une copie de ce tableau.
On rappelle que l'attribut shape permet de récupérer les dimensions d'un tableau numpy.
Q8. Écrire la fonction itere_J (V, rhos, frontiere) modifiant la totalité du tableau V[i,j] lors d'une seule itération et retournant l'erreur calculée conformément à l'équation (4).
Q9. Écrire la fonction poisson(f_iter, V, rhos, frontiere, eps) ayant pour premier argument une fonction du même type que celle définie à la question précédente, pour dernier argument eps le seuil arbitraire de convergence et dont le rôle est de modifier le tableau des potentiels jusqu'à convergence.
I. 5 - Améliorations
Méthode de Gauss-Seidel
C'est une modification de la méthode de Jacobi, pour laquelle on montre que la convergence est légèrement plus rapide. Supposons que l'on balaye le tableau des potentiels selon les indices et croissants : dans ces conditions, les points situés à gauche et en dessous du point courant ont déjà été calculés. On va utiliser ces nouvelles valeurs, probablement plus proches de la solution, dans la formule permettant le calcul de . Ceci donne l'algorithme de Gauss-Seidel :
Q10. Montrer qu'il n'est plus nécessaire de copier le tableau pour la mise à jour lors d'une itération en utilisant l'équation (5). Faut-il modifier la fonction nouveau_potentiel pour passer de la méthode de Jacobi à celle de Gauss-Seidel ?
Q11. Écrire la fonction itere_GS(V, rhos, frontiere) modifiant la totalité du tableau V[i,j] lors d'une seule itération et retournant l'erreur calculée conformément à l'équation (4).
Méthode de Gauss-Seidel adaptative
Les méthodes de Jacobi et de Gauss-Seidel n'utilisent pas la valeur de pour calculer . La méthode de sur-relaxation (Successive Over Relaxation method) consiste à calculer la nouvelle valeur d'un nœud comme une combinaison linéaire de la valeur courante et de celle donnée par le schéma de Gauss-Seidel. En introduisant le paramètre de relaxation , on a alors :
L'étude mathématique de cette relation permet de montrer les résultats suivants :
la méthode converge uniquement si et elle converge plus rapidement que la méthode de Gauss-Seidel si ,
il existe une valeur optimale de qui permet la convergence avec un nombre d'itérations en pour une valeur de fixée.
Pour la résolution de l'équation de Poisson envisagée dans ce problème (conditions de Dirichlet sur un maillage carré), on montre que la valeur optimale est :
Q12. Écrire la fonction nouveau_potentiel_SOR(V, rhos, frontiere, i, j, omega) retournant la nouvelle valeur du potentiel au point ( ) selon l'équation (6).
Q13. Écrire la fonction itere_SOR(V, rhos, frontiere) optimale modifiant la totalité du tableau V [i,j] lors d'une seule itération et retournant l'erreur calculée conformément à l'équation (4).
La résolution du problème peut alors se faire par un appel de la forme
poisson(itere_SOR, V, rhos, frontiere, eps),
les tableaux carrés V, rhos, frontiere étant de dimensions convenables pour représenter un maillage comportant nœuds.
Q14. Quelle est la complexité temporelle de l'appel précédent quand ?
La figure 2 représente, pour , la durée d'exécution (en secondes) en fonction de . Cette courbe est-elle en accord avec la complexité temporelle attendue ?
Quelle serait la durée d'exécution pour ? Commenter.
Figure 2 - Durée d'exécution (en secondes) de poisson(itere_SOR, V, rhos, frontiere, eps) en fonction de
I. 6 - Détermination du champ électrique
Connaissant le potentiel V , il est souvent nécessaire de calculer numériquement les composantes et du champ électrique au niveau des nœuds du maillage, qui sont alors conservées dans deux tableaux Ex et Ey dimensionnés correctement (ce sont donc des tableaux carrés de éléments).
Q15. Expliquer rapidement comment il serait possible de définir la fonction calc_ExEy(Ex, Ey, V , h ), permettant, à partir du tableau V et du pas du maillage h , le remplissage des deux tableaux Ex et Ey passés en arguments.
Remarque : on ne demande pas d'écrire le code de la fonction, juste de décrire précisément les étapes de calcul, ainsi que les différents cas à considérer.
Partie II - Deux études de cas
II. 1 - Fil cylindrique chargé uniformément
Étude théorique
On considère dans le vide un fil cylindrique infini d'axe et de rayon , portant une charge volumique constante .
Q16. En se plaçant en coordonnées cylindriques d'axe , montrer par des considérations de symétrie et d'invariance que le champ . En déduire la forme des surfaces équipotentielles.
Q17. En appliquant le théorème de Gauss, calculer le champ dans tout l'espace. Tracer rapidement l'allure de en fonction de .
Q18. On donne : . Calculer la valeur maximale de la norme du champ électrique, ainsi que la valeur pour .
Étude numérique
Pour pouvoir utiliser la méthode de Gauss-Seidel adaptative, on place le fil infini au centre d'une enceinte de longueur infinie et de section carrée ( ), portée au potentiel nul (figure 3).
Dans la suite, on prendra .
Figure 3 - Fil infini dans une enceinte de section carrée, portée au potentiel nul
Le programme permettant la résolution de ce problème commence ainsi :
# initialisations
eps0 = 8.85e-12 # epsilon_0
L = 20.0e-2 # 20 cm
N = 100 ; h = L/N # définition du maillage
rho = 1.00e-5 # densité vol. de charge rho
# les tableaux globaux numpy pour le cylindre chargé
rhos_cyl = np.zeros((N+1,N+1)) # tableau des valeurs de rho''
V_cyl = np.zeros((N+1,N+1)) # le potentiel
Ex_cyl = np.zeros((N+1,N+1)) # la composante Ex
Ey_cyl = np.zeros((N+1,N+1)) # la composante Ey
# le tableau définissant la frontière est initialement
# rempli entièrement par la valeur False
frontiere_cyl = np.zeros((N+1,N+1), bool)
Q19. Écrire la fonction dans_cylindre( ) retournant un résultat booléen indiquant si le point de coordonnées est à l'intérieur ou sur le bord du cercle de centre et de rayon .
Q20. Écrire la fonction initialise_rhos_cylindre(tab_rhos), initialisant le tableau rhos_cyl contenant les valeurs pour les nœuds du maillage.
Q21. Écrire la fonction initialise_frontiere_cylindre(tab_f), mettant à True les points appartenant à la frontière, donc de potentiel fixé.
La résolution numérique avec la méthode de Gauss-Seidel adaptative, utilisant les valeurs numériques précédentes et un seuil de convergence , mène à la figure 4 , où on a tracé un réseau de courbes équipotentielles, le potentiel et la composante du champ le long de l'axe de symétrie défini par . En outre, la valeur calculée de Ex_cyl[50, 50] est égale à -0.0023321214257521206.
Q22. Commenter le plus complètement possible ces résultats; on veillera, en particulier, à les comparer au modèle théorique (allure des courbes, valeurs numériques ...).
Figure 4 - Équipotentielles, et en fonction de pour
Une autre résolution est effectuée, avec une répartition de charges dans le cylindre différente de la précédente, utilisant la même valeur de la densité volumique . Elle mène aux courbes de la figure 5 (page suivante).
Figure 5 - Équipotentielles, et en fonction de pour la nouvelle répartition de charge
Q23. Déduire de ces courbes la répartition de charges dans le cylindre dans cette deuxième situation. Calculer le champ électrique en tout point pour cette répartition de charges dans le cas d'un cylindre infini seul dans l'espace.
On vérifiera que la valeur maximale du champ électrique calculée à l'aide de cette modélisation est compatible avec celle déduite de la figure 5.
II. 2 - Mouvement d'un électron dans un tube d'oscilloscope
La figure 6 montre un tube d'oscilloscope de petite dimension, dans lequel des électrons émis par la cathode sont accélérés et déviés vers un écran luminescent. La déviation est assurée par le passage des électrons entre les plaques de deux condensateurs plans : un pour la déviation horizontale, l'autre pour la déviation verticale. L'étude qui suit ne concernera que le condensateur responsable de la déviation verticale.
Figure 6 - Petit tube d'oscilloscope, de longueur d'environ 20 cm
On modélise la trajectoire d'un électron de la façon suivante (figure 7 page suivante, où la zone de déviation est grisée) :
on négligera l'effet de la pesanteur,
émis à vitesse nulle par effet thermo-électronique au niveau de la cathode portée au potentiel nul, l'électron est accéléré à l'aide d'une tension afin d'acquérir à l'entrée de la zone de déviation une vitesse ,
pendant son trajet dans la zone de déviation, il est soumis à un champ électrique lié aux potentiels des plaques du condensateur, de longueur et séparées par une distance ,
poursuivant son mouvement, il arrive sur la surface de l'écran à une distance de l'axe , l'écran étant situé à la distance du centre du condensateur.
Figure 7 - Schéma de la zone de déviation
Valeurs numériques :
masse d'un électron ; charge élémentaire .
Étude physique
Q24. En appliquant la conservation de l'énergie, calculer la vitesse de l'électron à l'entrée de la zone de déviation. Faire l'application numérique. Commenter.
Q25. On modélise les plaques de déviation comme un condensateur sans effets de bord : le champ électrique est donc considéré comme nul si et uniforme si , ses lignes de champ étant parallèles à l'axe . Exprimer le champ entre les plaques en fonction de et d.
Q26. On suppose que la vitesse d'entrée de l'électron dans la zone de déviation est . En appliquant les lois de la mécanique, établir l'équation de la trajectoire de l'électron entre les plaques (pour ).
Q27. Montrer que l'équation de la trajectoire pour est donnée par : .
En déduire que l'ordonnée du spot sur l'écran est : .
Faire l'application numérique pour .
Étude numérique
Pour savoir si la modélisation précédente est pertinente, on va envisager une détermination numérique de la trajectoire de l'électron. Pour cela, on place le condensateur de déviation dans une enceinte carrée au potentiel nul de côté , le centre du condensateur étant à du bord gauche de l'enceinte (figure 8 page suivante). Les autres caractéristiques électriques et géométriques sont les mêmes que précédemment.
La résolution numérique se déroule alors en deux étapes :
calcul du potentiel et du champ électrique par la méthode de Gauss-Seidel adaptative dans l'enceinte,
calcul de la trajectoire de l'électron à l'aide de la méthode d'Euler.
Figure 8 - Modélisation des plaques de déviation dans l'enceinte
Le programme permettant cette résolution numérique commence ainsi :
L = 0.100; N = 100; h = L/N
V0 = 950 ; Vp = 180
m = 9.11e-31 ; e = 1.60e-19
rhos_osc = np.zeros((N+1,N+1))
V_osc = np.zeros((N+1,N+1))
Ex_osc = np.zeros((N+1,N+1))
Ey_osc = np.zeros((N+1,N+1))
frontiere_osc = np.zeros((N+1,N+1), dtype=bool)
Calcul du potentiel et du champ électrique dans l'enceinte
Q28. Quelles sont les valeurs qui doivent être contenues dans le tableau rhos pour le problème considéré?
Q29. Écrire la fonction initialise_frontiere_condensateur(tab_V, tab_f), permettant l'initialisation des tableaux V_osc et frontiere_osc à l'aide de la ligne de code suivante :
Pour pouvoir utiliser la méthode d'Euler, il est nécessaire de pouvoir calculer les composantes et du champ pour et . Cependant, la méthode de résolution (associée à la fonction calc_ExEy définie dans la question 15) ne permet de calculer les composantes Ex_osc et Ey_osc qu'aux nœuds du maillage.
Soit un point de coordonnées . Ce point est dans la cellule , où et . Posons et .
Q30. Montrer alors que :
Écrire de même la formule permettant de calculer . En déduire qu'il est possible de calculer les composantes du champ en tout point de .
On supposera dans la suite que les fonctions val_Ex(Ex,Ey,x,y,h) et val_Ey(Ex,Ey,x,y,h) sont définies et retournent les valeurs des composantes du champ électrique pour le point de coordonnées calculées à l'aide des formules précédentes.
Calcul de la trajectoire par la méthode d'Euler
Q31. Montrer que les équations permettant de décrire le mouvement de l'électron par la méthode d'Euler sont les suivantes, avec comme petit incrément temporel et les variations pendant des grandeurs :
Q32. Compte tenu du changement de la position d'origine du repère (imposée par la résolution numérique de l'équation de Poisson, figure 8), quelles sont les conditions initiales du mouvement de l'électron?
Déterminer pour calculer environ 200 points successifs le long de la trajectoire. Faire l'application numérique.
Q33. En tenant compte des réponses aux questions précédentes, compléter le code d'initialisation des variables de la simulation (****** dans le code suivant) :
Npts = 200 # nombre de points pour le tracé de la trajectoire
v0 = ****** # vitesse initiale de l'électron
dt = ****** # incrément temporel
# tableaux des coordonnées x et y de l'électron
lx = np.zeros(Npts) ; ly = np.zeros(Npts)
# tableaux des vitesses en x et en y
lvx = np.zeros(Npts) ; lvy = np.zeros(Npts)
#conditions initiales
lx[0] = ******; ly[0] = ******
lvx[0] = *******; lvy[0] = *******
Q34. Écrire les lignes de code implémentant la boucle de remplissage des tableaux , lvy selon la méthode d'Euler.
Comparaison théorie/simulation
La figure 9 montre le résultat de la simulation précédente. On y voit le réseau de courbes équipotentielles, ainsi que deux trajectoires 1 et 2 , l'une étant associée au calcul théorique, l'autre à la simulation numérique. Chaque trajectoire est constituée de 200 points de calcul séparés d'une durée .
Q35. Reproduire sommairement sur la copie la figure 9, y ajouter le tracé de quelques lignes de champ orientées dans les différentes parties de la zone de déviation.
Identifier, en le justifiant, chaque trajectoire. Expliquer pourquoi la trajectoire 1 est plus courte que la trajectoire 2.
À votre avis, peut-on se contenter de l'étude théorique pour prévoir le point d'impact de l'électron sur l'écran? (On attend une réponse chiffrée.)
Figure 9 - Tracé dans la zone de déviation de quelques courbes équipotentielles, des trajectoires de l'électron (théorique et simulée numériquement)
Aide-mémoire numpy/matplotlib/pyplot
Importation des bibliothèques
Les bibliothèques sont importées de la façon suivante :
import math
import numpy as np
import matplotlib.pyplot as plt
Manipulation des tableaux numpy
La création d'un tableau numpy à deux dimensions dont toutes les valeurs sont initialisées à 0 est faite par l'instruction np.zeros(format), format étant un doublet de la forme (n_lignes , n_colonnes) :