Méthode PN

La méthode P_N permet la résolution de l'équation du transfert radiatif (ou équation de Boltzmann) utilisée pour la propagation des particules telles photons, neutrons, neutrinos, etc. Cette méthode été introduite[1] par James Jeans (1917)[2]. Elle consiste en un développement de la solution sur une base de polynômes orthogonaux.

L'équation du transfert radiatif

Dans le cas d'un milieu stationnaire unidimensionnel en espace comportant émission, diffusion élastique (sans changement de fréquence) et absorption la propagation est décrite par l'équation équation intégro-différentielle linéaire suivante

\mu {\frac {\mathrm {d} L(x,\mu )}{\mathrm {d} x}}+\kappa _{t}(x)L(x,\mu )=S(x,\mu )+2\pi \kappa _{d}\int _{-1}^{1}{\mathcal {P}}(\mu ,\mu ')L(x,\mu )\mathrm {d} \mu

où

$L(x,\mu )$	luminance (spectrale ou intégrée),
$x$	variable d'espace
$\mu =\cos(\theta )$	donne la direction de propagation pour une fonction de phase supposée de révolution,
$\theta$	angle de colatitude pour des coordonnées sphériques,
$\kappa _{t}(x)=\kappa _{a}(x)+\kappa _{d}(x)$	coefficient d'extinction totale,
$\kappa _{a}(x)$	coefficient d'absorption,
$\kappa _{d}(x)$	coefficient d'extinction par diffusion,
${\mathcal {P}}(\mu ,\mu ')$	fonction de phase (supposée entièrement définie par la déviation au cours d'une interaction),
$S(x,\mu )$	fonction source volumique.

Mise en œuvre de la méthode

La base du développement choisie sera :

en bidimentionnel des harmoniques sphériques,
si la solution est à symétrie azimutale (cas traité ici) on choisit généralement des polynômes de Legendre notés P, quoique beaucoup d'autres choix soient possibles.

Le développement

La luminance est supposée de forme suivante, à variables séparées

L(x,\mu )=\sum _{i=0}^{N}{\frac {2i+1}{2}}P_{i}(\mu )L_{i}(x)

L'orthogonalité des polynômes s'écrit

\int _{-1}^{1}P_{i}(\mu )P_{j}(\mu )\mathrm {d} \mu ={\frac {2}{2i+1}}\delta _{ij}

où δ_ij est la fonction de Dirac.

En multipliant l'équation donnant L(x, μ) par P_j et en intégrant sur μ il vient, compte tenu de la relation d'orthogonalité

L_{i}(x)=\int _{-1}^{1}P_{i}(\mu )L(x,\mu )\mathrm {d} \mu

Les premiers polynômes ont une interprétation physique

l'énergie volumique correspond à L₀

E(x)={\frac {2\pi }{c}}\int _{-1}^{1}L(x,\mu )\mathrm {d} \mu ={\frac {2\pi }{c}}L_{0}(x)

où c est la vitesse de la lumière.

l'exitance (flux) correspond à L₁ (loi de Lambert)

M(x)=2\pi \int _{-1}^{1}L(x,\mu )\mu \mathrm {d} \mu =\pi L_{1}(x)

la pression

P(x)={\frac {2\pi }{c}}\int _{-1}^{1}L(x,\mu )\mu ^{2}\mathrm {d} \mu

Les moments de Legendre de la fonction de phase

La déviation au cours d'une interaction est supposée ne dépendre que de du produit scalaire des angles définissant les directions de propagation Ω = (θ, Φ) et Ω' = (θ', Φ') et on définit la déviation par son cosinus α tel que

\alpha =\mathbf {\Omega } \cdot \mathbf {\Omega } '=\mu \mu '+{\sqrt {1-\mu ^{2}}}{\sqrt {1-\mu '^{2}}}\cos {(\phi -\phi ')}

La fonction de phase est développée en polynômes de Legendre

{\mathcal {P}}(\alpha )={\frac {1}{2\pi }}\sum _{i=0}^{\infty }{\frac {2i+1}{2}}g^{i}P_{i}(\alpha )

où les g_i sont les moments de Legendre de la fonction de phase

g_{i}=\int _{-1}^{1}{\mathcal {P}}(\alpha )P_{i}(\alpha )\mathrm {d} \alpha

Le premier moment est souvent utilisé pour caractériser la dissymétrie de la fonction de phase, voire pour construire une fonction standard comme celle de Henyey-Greenstein[3].

Les système P_N

En multipliant l'équation de transfert par chacun des $P_{j}$ et en intégrant sur μ en tenant compte de la propriété d'orthogonalité des polynômes on obtient un système de N+1 équations pour N+2 inconnues $L_{0},L_{1},...L_{N+1}$ . On fera donc une hypothèse pour clore le système. La plus simple consiste à imposer $L_{N+1}=0$ mais il existe des méthodes plus élaborées exprimant L_N+1 en fonction des termes connus[4].

{\begin{array}{rcl}{\frac {\mathrm {d} L_{1}}{\mathrm {d} x}}+(\kappa _{a}+\kappa _{d})L_{0}&=&S\\[0.6em]{\frac {i+1}{2i+1}}{\frac {\mathrm {d} L_{i+1}}{\mathrm {d} x}}+\kappa _{i}L_{i}+{\frac {i}{2i+1}}{\frac {\mathrm {d} L_{i-1}}{\mathrm {d} x}}&=&0\,,\qquad 1<i<N-2\\[0.6em]{\frac {N}{2N+1}}{\frac {\mathrm {d} L_{N-1}}{\mathrm {d} x}}+\kappa _{N}L_{N}&=&0\end{array}}

où

\kappa _{i}(x)=\kappa _{a}(x)+(1-g_{i})\kappa _{d}(x)

On montre[5] que la solution tend vers la solution physique lorsque N → ∞.

Elle est équivalente en tous points (difficulté de mise en œuvre, performances en durée de résolution) à la méthode S_N, à l'exception des pathologies, différentes dans les deux cas : cette méthode est sensible au phénomène de Gibbs[6].

Lien avec l'approximation d'Eddington

L'approximation d'Eddington a été obtenue indépendamment des travaux de James Jeans mais elle constitue en fait la méthode P₁. Cette approximation est due à Arthur Eddington[7]. En utilisant l'expression des deux premiers polynômes de Legendre et les quantités définies plus haut le développement s'écrit

L(\mathbf {x} ,\nu )={\frac {1}{2}}\left[L_{0}(\mathbf {x} )+3\mu L_{1}(\mathbf {x} )\right]

D'où les quantités

E={\frac {2\pi }{c}}L_{0}\,,\qquad M=\pi L_{1}\,,\qquad P={\frac {2\pi }{3c}}L_{0}={\frac {E}{3}}

Cette dernière expression constitue l'approximation d'Eddington.

Méthode SP_N

Cette méthode constitue une simplification de la méthode P_N (S pour simplified) proposée par E. M. Gelbart de manière heuristique[8] et justifié ultérieurement par l'analyse asymptotique[9] ou par une approche variationnelle[10]. Dans cette méthode on omet la dérivée pour les valeurs paires de i, d'où

L_{i}=-{\frac {1}{\kappa _{i}}}\left({\frac {i+1}{2i+1}}{\frac {\mathrm {d} L_{i+1}}{\mathrm {d} x}}+{\frac {i}{2i+1}}{\frac {\mathrm {d} L_{i-1}}{\mathrm {d} x}}\right)

En rapportant cette expression dans le système P_N on obtient le système SP_N : pour i impair

{\begin{array}{rcl}{\frac {\mathrm {d} }{\mathrm {d} x}}\left({\frac {1}{3\kappa _{1}}}{\frac {\mathrm {d} L_{0}}{\mathrm {d} x}}\right)+{\frac {\mathrm {d} }{\mathrm {d} x}}\left({\frac {2}{3\kappa _{1}}}{\frac {\mathrm {d} L_{2}}{\mathrm {d} x}}\right)&=&-\kappa _{a}S\\[0.6em]{\frac {i(i-1)}{(2i+1)(2i-1)}}{\frac {\mathrm {d} }{\mathrm {d} x}}\left({\frac {1}{\kappa _{i-1}}}{\frac {\mathrm {d} L_{i-2}}{\mathrm {d} x}}\right)+{\frac {(i+1)(i+2)}{(2i+1)(2i+3)}}{\frac {\mathrm {d} }{\mathrm {d} x}}\left({\frac {1}{\kappa _{i+1}}}{\frac {\mathrm {d} L_{i+2}}{\mathrm {d} x}}\right)&&\\[0.6em]+{\frac {\mathrm {d} }{\mathrm {d} x}}\left[\left({\frac {(i+1)^{2}}{(2i+1)(2i+3)\kappa _{i+1}}}+{\frac {i^{2}}{(2i+1)(2i-1)\kappa _{i-1}}}\right){\frac {\mathrm {d} L_{i}}{\mathrm {d} x}}\right]-\kappa _{i}L_{i}&=&0\end{array}}

Il s'agit d'un système d'équations de type diffusion donc sans difficulté numérique. Le nombre d'équations a été réduit par un facteur 2.

Un exemple simple

Milieu homogène semi-infini ; exposants en fonction de N.

Supposons un milieu homogène semi-infini à diffusion isotrope. Le système P_N est alors le système linéaire suivant

\left|{\begin{array}{ccccccc}0&{\frac {1}{\kappa _{a}}}&0&0&...&0\\[0.6em]{\frac {1}{3\kappa _{1}}}&0&{\frac {2}{3\kappa _{1}}}&0&...&0\\[0.6em]0&{\frac {2}{5\kappa _{2}}}&0&{\frac {3}{5\kappa _{2}}}&...&0\\[0.6em]0&0&{\frac {3}{7\kappa 4}}&0&...&0\\[0.6em]...&...&...&...&...&{\frac {N}{(2N+1)\kappa _{N}}}\\[0.6em]0&0&0&...&{\frac {N+1}{(2N+3)\kappa _{N}}}&0\end{array}}\right|\left|{\begin{array}{c}{\frac {\mathrm {d} L_{0}}{\mathrm {d} x}}\\[0.6em]{\frac {\mathrm {d} L_{1}}{\mathrm {d} x}}\\[0.6em]{\frac {\mathrm {d} L_{2}}{\mathrm {d} x}}\\[0.6em]...\\[0.6em]{\frac {\mathrm {d} L_{N-1}}{\mathrm {d} x}}\\[0.6em]{\frac {\mathrm {d} L_{N}}{\mathrm {d} x}}\end{array}}\right|=-\left|{\begin{array}{c}L_{0}\\[0.6em]L_{1}\\[0.6em]L_{2}\\[0.6em]...\\[0.6em]L_{N-1}\\[0.6em]L_{N}\end{array}}\right|+\left|{\begin{array}{c}S\\[0.6em]0\\[0.6em]0\\[0.6em]0\\[0.6em]0\\[0.6em]0\\[0.6em]0\end{array}}\right|

On suppose que N est impair : avec cette condition la matrice ci-dessus qui est diagonalisable a N + 1 valeurs propres réelles par paires de signe opposé, ce qui justifie le choix d'une valeur impaire. Sa solution se met sous la forme d'une solution générale du système homogène constitué de (N + 1) / 2 exponentielles négatives (les valeurs positives étant exclues pour obtenir une solution finie) et du terme constant S / κ_a constituant une solution particulière

L_{i}={\frac {S}{\kappa _{a}}}+\sum _{1}^{\frac {N+1}{2}}\beta _{i}e^{-\gamma _{i}x}

La courbe donnant les valeurs de γ_i montre la convergence de la méthode avec l'ordre. Les faibles valeurs correspondantes aux termes variant le plus lentement convergent plus rapidement. D'une façon générale l'erreur varie comme N^-2[11].

Références

(en) Michael M. Modest, Radiative Heat Transfer, Academic Press, 2003, 822 p. (ISBN 0-12-503163-7, lire en ligne)
(en) J. H. Jeans, « The Equations of Radiative Transfer of Energy », Monthly Notices of the Royal Astronomical Society, vol. 78,‎ 1917, p. 28-36 (lire en ligne)
(en) J. Patrick Harrington, « The Henyey-Greenstein phase function », sur University of Maryland
(en) M. Schäfer, M. Frank et C. D. Levermore, « Diffusive Corrections to P_N Approximations », Multiscale Modeling and Simulation, vol. 9, n^o 1,‎ 2011, p. 1-28
(en) A. M. Sultagazin, « Convergence of the Method of Spherical Harmonics for the Non-Stationary Transport Equation », USSR Computational Mathematics and Mathematical Physics, vol. 14, n^o 1,‎ 1974, p. 165-176
(en) Benjamin Seibold, « StaRMAP », sur Université Temple
(en) A. Eddington, The Internal Constitution of the Stars, Dover, 1926, p. 322
(en) E. M. Gelbart, Simplified Spherical Harmonics Equations and their Use in Shielding Problems (Report WAPD-T-1182 (Rev.1)), Bettis AtomicPower Laboratory, 1961
(en) G. C. Pomraning, « Asymptotic and Variational Derivation of the Simplified PN Approximation », Annals of Nuclear Energy, vol. 20, n^o 9,‎ 1993, p. 623-637
(en) R. G. McClaren, « Theoretical Aspects of the Simplified P_n Equations », Transport Theory and Statistical Physics, vol. 39,‎ 2011, p. 73-109
(en) B. Davison, « On the Rate of Convergence of the Spherical Harmonics, Isotropic Scattering », Canadian Journal of Physics, vol. 38, n^o 11,‎ 1960, p. 1536-1545

Voir aussi

Méthode MN

Portail de la physique
Portail de l’optique
Portail de l’astronomie

Cet article est issu de Wikipedia. Le texte est sous licence Creative Commons - Attribution - Partage dans les Mêmes. Des conditions supplémentaires peuvent s'appliquer aux fichiers multimédias.

[Modest-1] (en) Michael M. Modest, Radiative Heat Transfer, Academic Press, 2003, 822 p. (ISBN 0-12-503163-7, lire en ligne)

[2] (en) J. H. Jeans, « The Equations of Radiative Transfer of Energy », Monthly Notices of the Royal Astronomical Society, vol. 78,‎ 1917, p. 28-36 (lire en ligne)

[3] (en) J. Patrick Harrington, « The Henyey-Greenstein phase function », sur University of Maryland

[4] (en) M. Schäfer, M. Frank et C. D. Levermore, « Diffusive Corrections to P_N Approximations », Multiscale Modeling and Simulation, vol. 9, n^o 1,‎ 2011, p. 1-28

[5] (en) A. M. Sultagazin, « Convergence of the Method of Spherical Harmonics for the Non-Stationary Transport Equation », USSR Computational Mathematics and Mathematical Physics, vol. 14, n^o 1,‎ 1974, p. 165-176

[6] (en) Benjamin Seibold, « StaRMAP », sur Université Temple

[7] (en) A. Eddington, The Internal Constitution of the Stars, Dover, 1926, p. 322

[8] (en) E. M. Gelbart, Simplified Spherical Harmonics Equations and their Use in Shielding Problems (Report WAPD-T-1182 (Rev.1)), Bettis AtomicPower Laboratory, 1961

[9] (en) G. C. Pomraning, « Asymptotic and Variational Derivation of the Simplified PN Approximation », Annals of Nuclear Energy, vol. 20, n^o 9,‎ 1993, p. 623-637

[10] (en) R. G. McClaren, « Theoretical Aspects of the Simplified P_n Equations », Transport Theory and Statistical Physics, vol. 39,‎ 2011, p. 73-109

[11] (en) B. Davison, « On the Rate of Convergence of the Spherical Harmonics, Isotropic Scattering », Canadian Journal of Physics, vol. 38, n^o 11,‎ 1960, p. 1536-1545