# Influence de la variation de la charge moyenne dans une simulation numérique d'un arc sous vide en présence d'un AMF

Claire Simonnet, Pierre Freton, Jean-Jacques J Gonzalez, Frank Reichert, Arkadz Petchanka

► To cite this version:

Claire Simonnet, Pierre Freton, Jean-Jacques J Gonzalez, Frank Reichert, Arkadz Petchanka. Influence de la variation de la charge moyenne dans une simulation numérique d'un arc sous vide en présence d'un AMF. Journal International de Technologie, de l'Innovation, de la Physique, de l'Energie et de l'Environnement, 2025, 9, 10.52497/jitipee.v9i1.370. hal-05356888

HAL Id: hal-05356888

https://hal.science/hal-05356888v1

Submitted on 10 Nov 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
JOURNAL INTERNATIONAL DE TECHNOLOGIE, DE L'INNOVATION,
DE LA PHYSIQUE, DE L'ENERGIE ET DE L'ENVIRONNEMENT

Influence de la Variation de la Charge
Moyenne dans une Simulation Numérique
d'un Arc Sous Vide en Présence d'un AMF.

C. Simonnet, P. Freton, J.J. Gonzalez, F. Reichert, A. Petchanka

<figure><img src="image_2.png" /></figure>

ISSN : 2428-8500

DOI : 10.52497/jitipee.v9i1.370

Le sujet de cet article a été présenté lors du colloque sur les arcs électriques (CAE XVII) les 17 et 18 mars 2025 à Limoges.

L'article publié a fait l'objet d'une expertise indépendante par deux spécialistes du domaine.
# Influence de la variation de la charge moyenne dans une simulation numérique d'un arc sous vide en présence d'un AMF.

Claire Simonnet(1), Pierre Freton(1), Jean-Jacques Gonzalez(1),

Frank Reichert(2), Arkadz Petchanka(2)

claire.simonnet@laplace.univ-tlse.fr

**Résumé** – Un moment délaissée, l'étude des disjoncteurs moyenne et haute tension utilisant le vide comme milieu de coupure connait un nouvel essor depuis quelques années. En effet, les normes environnementales tendent à diminuer drastiquement l'utilisation de gaz à effet de serre dans les disjoncteurs et à favoriser les alternatives. De ce fait, l'utilisation de dispositifs dans le vide comme moyen de coupure s'avèrent être une alternative pertinente. Cependant, peu d'études théoriques discutant des hypothèses utilisées existent sur ces dispositifs. De plus, la difficulté d'effectuer des mesures expérimentales dans un environnement à basse pression (environ $10^{-5}$ Pa) et fort courant (plusieurs kA) contraint les chercheurs et les industriels à développer des modèles numériques afin de comprendre les phénomènes physiques et améliorer ces systèmes de coupure. Dans cette communication nous présentons un modèle numérique tridimensionnel (3D) qui décrit un arc sous vide à fort courant (15 kA) en régime subsonique entre des contacts en cuivre sous l'influence d'un champ magnétique axial (AMF). Le modèle est basé sur une approche MHD du plasma inter-électrodes à deux températures (électronique et ionique). L'étude est réalisée sur le logiciel commercial Ansys Fluent 23 R2 utilisant des développements en langage C et des UDF. Nous décrivons le modèle, les conditions aux limites et les principaux résultats de simulation. Les champs de température ionique et électronique, de pression absolue et de vitesse ionique seront donnés. Sur la base de cette modélisation d'arc, nous présentons une étude paramétrique de la charge moyenne du plasma en lien avec les hypothèses faites dans la littérature. Les résultats montrent qu'en considérant une charge moyenne plus élevée et proche de la réalité, la température du milieu augmente de plusieurs eV.

**Mots clés** : arc électrique, plasma, vide, disjoncteur, simulation

DOI : 10.18145/jitipee.v9i1.370

<footnote>(1) Laplace, UMR 5213 CNRS-UPS-INP, Université Paul Sabatier, 118 route de Narbonne, 31062, Toulouse, France.</footnote>

<footnote>(2) Siemens Energy, Paulsternstr. 26, 13629, Berlin, Germany</footnote>
## Introduction

De manière générale les disjoncteurs sont des appareils permettant d'interrompre le courant dans un réseau électrique lors d'une surcharge ou d'un court-circuit [1]. Dans les deux cas, ces phénomènes sont néfastes pour l'utilisateur, les installations et l'environnement. L'objectif du disjoncteur consiste en la protection du système électrique, l'isolement des composants défectueux du reste du réseau et la garantie de fonctionnement des autres parties du réseau de manière sécurisée [2]. Le temps de coupure doit être proche du temps caractéristique du court-circuit allant de la microseconde à la milliseconde [3]. Les disjoncteurs moyenne tension (MT) jusqu'à 52 kV et haute tension (HT) plus de 52 kV, utilisent comme milieux de coupure principalement le SF₆ et le vide [4]. Les sous-produits de décomposition du SF₆ sont néfastes pour l'environnement [5] et les normes environnementales imposent une diminution de leur utilisation dans les systèmes de coupure. Par ailleurs, les nouveaux gaz envisagés n'ont pas à l'heure actuelle de performances équivalentes. Le vide est donc un bon candidat comme milieu de coupure à étudier pour améliorer les disjoncteurs car c'est un bon isolant (le libre parcours moyen des molécules entre deux collisions est de l'ordre de la centaine de mètre dans un environnement à basse pression) [3], leur cout de fonctionnement est faible et ils ont un moindre impact sur l'environnement. Dans un disjoncteur à vide (figure1) se trouve une ampoule à vide où l'arc électrique se forme entre les deux électrodes cylindriques, de quelques centimètres de longueur, encapsulées sous une pression de 10⁻¹ à 10⁻⁵ Pa [3] [6] [7].

<figure><img src="image_3.png" /><figcaption>**Figure 1** : Schéma d'une ampoule à vide de Siemens Energy : vue extérieure à gauche, vue intérieure au centre et un zoom de l'espace inter-électrodes à droite [6].</figcaption></figure>

Idéalement, en position fermée, la résistance de l'ampoule est nulle et le courant circule normalement. Lors d'un courant de défaut dans le circuit, la résistance du disjoncteur doit être infinie pour ne plus laisser passer le courant. Dans un cas réel, la synchronisation entre la vitesse d'ouverture des contacts (environ 1 m/s) [8] et le passage du courant n'est pas parfaite. Ce déphasage engendre un arc électrique durant le processus de coupure qui se maintient jusqu'au passage à zéro du courant alternatif de fréquence 50 Hz [1] [9]. Lors du passage du courant, la densité de courant augmente brutalement, et aux derniers points de contact entre les surfaces des électrodes, on observe une fusion locale des électrodes et un pont métallique liquide qui se forme sous l'influence du chauffage par effet Joule [3]. Les contacts continuent de se séparer
et le pont métallique chauffe et s'étire dans l'espace inter-électrode. Son instabilité conduit à sa rupture et il apparaît un arc électrique concentré dans un plasma de vapeurs métalliques issues de l'explosion du pont de métal fondu [3]. En fonction de l'amplitude de la densité de courant, l'arc adopte un régime de fonctionnement diffus (dont l'écoulement des ions est supersonique) ou concentré représenté respectivement sur la figure 2 a) et b) (sans AMF expliqué plus tard) et se maintient tant que le courant ne passe pas par zéro [9].

<figure><img src="image_4.png" /><figcaption>**Figure 2 :** Schéma du mode de fonctionnement de l'arc sans AMF a) arc diffus supersonique et b) arc concentré ; et avec AMF c) arcs multiple, d) arc diffus subsonique et e) arc diffus colonnaire.</figcaption></figure>

Le régime diffus s'établit pour des valeurs de courant de quelques A à quelques kA [9], voir figure 2 a). La cathode émet, via des spots cathodiques (SC) des points lumineux circulaires, mobiles [10] de quelques micromètre de diamètre [11], des jets de plasma composés de neutres, d'électrons et d'ions. Ces jets se mélangent ensuite dans l'espace inter-électrodes pour former un plasma diffus et homogène réparti dans tout le volume. L'anode est un collecteur passif de charges. Dans ce contexte, l'érosion de la cathode est faible et est comprise entre 35 et 40 40µg/C [12].

Pour de fortes intensités de courant, plusieurs kA, on retrouve l'arc en mode concentré où la colonne de plasma se contracte fortement par effet Hall. L'effet Hall traduit la déviation des charges par un champ magnétique azimutale créé par les lignes de courant plus intenses [13]. Le resserrement des lignes de courant vers l'axe concentre l'arc à l'anode sur une région plus petite que pour le mode diffus. Les SC se regroupent et le plasma forme un canal conducteur plus lumineux. L'anode est considérée active, c'est à dire émettrice de vapeurs métalliques et de goutes (elles peuvent être responsables de l'échec de la coupure) [3]. L'énergie déposée à l'anode augmente fortement sur une zone réduite, le matériau s'échauffe et émet des particules neutres qui sont ionisées par les électrons incidents.

La transition entre le régime d'arc diffus supersonique et l'arc concentré, est encore sujet à réflexion dans la littérature, mais plusieurs hypothèses sont en discussion. Avec l'augmentation du courant de plusieurs kA, lorsque le régime des ions devient transsonique à l'anode, une onde de choc s'y produit et se déplace rapidement vers la cathode [9] [14]. Suite à ce phénomène, les effets magnétohydrodynamiques, dont l'effet des forces de Lorentz affecte la décharge, où
l'interaction entre le courant d'arc et le champ magnétique auto-induit joue sur l'augmentation de la pression durant l'écoulement du plasma de la cathode à l'anode. La diminution de la vitesse d'écoulement des ions, et la diminution du nombre de Mach conduit à la constriction du plasma et à l'augmentation de la densité de courant [14]. La transition entre ces deux régimes dépend fortement du rayon des électrodes et de la distance inter-électrodes en plus de la valeur du courant [9].

Afin de diminuer l'échauffement des contacts, leur dégradation et d'augmenter la fiabilité du disjoncteur, on « joue » sur les champs magnétiques externes pour contrôler l'arc. L'idée est d'augmenter la plage de courant où l'arc est en mode diffus. Le champ magnétique peut être transverse (TMF) [15], radial (RMF) [16] ou axial (AMF) [17] et dépend de l'usinage spécifique des électrodes. Dans certains cas, le champ magnétique externe peut aussi être généré par une bobine de Helmholtz enroulées autour de l'ampoule à vide [1]. Les plasmas d'arc sous vide combinés avec l'utilisation d'un champ magnétique axial externe (AMF : B<sub>AMF</sub>) ont été utilisés avec succès dans les interrupteurs à vide commerciaux pour les réseaux de distribution moyenne tension depuis plusieurs décennies.

Pour notre projet nous nous consacrerons à l'étude de la configuration d'une ampoule à vide avec un champ magnétique axial externe AMF pour lequel trois régimes de fonctionnement de l'arc sont possibles. A faible courant et AMF élevé, ce sont des arcs multiples que l'on observe, schématisés sur la figure 2 c). Les jets de plasma issus des SC ne se mélangent pas et atteignent l'anode séparément [18]. A fort courant et fort AMF, on distingue un arc diffus subsonique dont l'aspect est similaire à l'arc diffus subsonique. Par la force de l'AMF, l'onde de choc ne peut pas se produire et l'arc ne transite pas vers le mode concentré, comme sur la figure 2 d). Puis pour un courant et un AMF fortement élevés nous observons une colonne centrale lumineuse plus dense plongée dans un plasma plus diffus [14], voir sur la figure 2 e). Les SC sont très denses au niveau de la cathode et la densité locale du plasma est très élevée. Dans ce mode, l'anode est toujours fortement émettrice et les électrodes sont de nouveau sujettes à une intense érosion. Il est important de noter que les plages de distinction des modes d'arc varient d'un auteur à l'autre.

<figure><img src="image_5.png" /><figcaption>**Figure 3 :** Illustration d'un arc diffus.</figcaption></figure>

Le plasma diffus subsonique est finalement le régime de fonctionnement qui supportera le mieux la coupure pour des plages de courant MT et HT, et limitera la destruction des électrodes. Sur la figure 3 nous retrouvons un schéma de la représentation des zones composant le plasma
diffus subsonique. L'extinction de ces arcs de commutation à des courants élevés (plusieurs dizaines de kA) et à des tensions élevées (plus de 100 kV) est un défi technologique qui nécessite l'aide de modèles numériques. L'amélioration de la conception et du fonctionnement des interrupteurs à vide dépend en grande partie d'une bonne compréhension de la physique de l'arc électrique. Néanmoins, il est encore difficile de conclure sur l'influence de chaque paramètre et les bases des modèles utilisés ne sont parfois pas toujours bien expliquées dans la littérature. Ceci nous oblige donc à revenir sur des travaux existants [14] [19] [20] [21] pour mettre en œuvre notre propre modèle numérique 3D sur une géométrie simplifiée afin de faire des études paramétriques et étudier l'influence des hypothèses classiquement utilisées. Parmi les hypothèses utilisées nous avons choisi de discuter l'hypothèse d'une charge constante en espace et en temps.

Ainsi dans une première partie nous posons les bases du modèle avec le système d'équations et les hypothèses, puis nous présentons pour une configuration donnée accessible dans la littérature les champs des grandeurs principales avant de discuter l'influence de la charge Z sur les résultats.

# 1. Simulation numérique

## 1.1. Configuration

Nous étudions une géométrie 3D présentée sur la figure 4. Nous considérons deux électrodes cylindriques de cuivre, parallèles, de diamètre D = 50 mm et de hauteur d = 5 mm, et espacées d'une hauteur h = 10 mm. Seule une demi-alternance sera étudiée, ainsi, sur la figure 4 on définit l'anode en haut et la cathode en bas. Seul le domaine inter-électrodes est modélisé, mais des conditions aux limites (C.L) à l'anode (en rouge), à la cathode (en bleu) ainsi qu'aux frontières latérales, nous permettent de fermer le système.

<figure><img src="image_6.png" /><figcaption>**Figure 4 :** Schéma de la configuration du modèle numérique.</figcaption></figure>

Un champ magnétique externe est appliqué dans l'espace inter-électrode noté $\overrightarrow{B_{AMF}}$. Lors du processus de coupure, les flux de particules ont un mouvement hélicoïdal et sont dirigés de la cathode vers l'anode. Ce type de géométrie est classiquement étudié dans la littérature, nous la retrouvons notamment dans les travaux de E. Schade [14] basés sur les équations développées par Braginskii [22]. Les jets de plasmas, les SC, et la gaine anodique ne sont pas modélisés mais pris en compte via des C.L surfaciques. Le modèle 3D stationnaire/transitoire est mis en place sur le logiciel Ansys Fluent 23 R2.
## 1.2. Hypothèses

La physique de l'arc électrique dans le vide dans les disjoncteurs moyenne et haute tension est très peu documentée. Les hypothèses du modèle classiquement utilisées dans la littérature sont rappelées ici :

### La zone fluide (Plasma) :

* Seule la partie plasma inter-électrodes est décrite par les équations. Les électrodes ne sont pas réellement représentées mais sont prises en compte par des conditions aux limites.

* Le plasma est un mélange d'ions et d'électrons respectivement de masse $m_i$ et $m_e$ et de charge moyenne fixe $Z = 1.9$ [23] [24]. Ce milieu est considéré totalement ionisé et quasi neutre avec $n_e = Zn_i$ où $n_i$ est la densité d'ions et $n_e$ est la densité d'électrons. Nous considérons les ions du cuivre suivant : $Cu^+$, $Cu^{2+}$, $Cu^{3+}$. Les particules neutres ne sont pas prises en comptes.

* Le libre parcours moyen concernant les interactions entre particules chargées de la même espèce (électrons-électrons, ions-ions) est toujours plus petit ($10^{-5}$ m pour les ions à $10^{-7}$ m pour les électrons) que les dimensions caractéristiques de l'ordre du cm [4]. Dans ces conditions, le nombre K de Knuden est $K \ll 1$ pour les électrons et les ions, correspondant à un milieu continu. Nous utilisons donc l'approximation fluide pour les électrons et les ions [10]. Les fluides sont considérés comme des gaz parfait.

* Le plasma est considéré hors équilibre thermique : les espèces peuvent être représentées par deux températures, une pour les ions ($T_i$) et une pour les électrons ($T_e$). Les ions et les électrons sont caractérisés par une fonction de distribution Maxwellienne. Nous avons supposé un écoulement laminaire.

* On considère également une conductivité électronique anisotrope due à la présence du champ magnétique axial (Ces développements ne sont pas systématiquement pris en compte dans les travaux de la littérature).

* Les pertes d'énergies par rayonnement ne sont pas prises en comptes dans cette étude.

### Zone cathodique :

* La surface de la cathode est considérée comme une entrée de flux de masse dans le fluide. Un profil de densité de courant Gaussien est utilisé.

* Les particules dans le milieu plasma proviennent de l'érosion de la cathode. Un flux de masse de cuivre est imposé. C. W. Kimblin [12] mentionne pour le cuivre, en se basant sur les travaux de la littérature un taux d'érosion compris entre 35 et 40 µg/C, d'autres auteurs utilisent 50 µg/C [1, 9, 25], dans notre cas nous avons pris cette dernière valeur car elle fait référence à des travaux plus récents.

### Zone anodique :

* L'électrode collecte les flux de plasma – elle est considérée comme une sortie avec une condition de flux. Les vapeurs métalliques qui potentiellement pourraient être émises ne sont pas considérées dans notre modèle.

* Une gaine existe et régule les flux d'électrons et d'ions. Un potentiel spécifique existe au sein de la gaine.
Pour le champ magnétique axial externe :

* Un champ uniforme est considéré et une valeur de 5 mT/kA est généralement utilisée [14]. Dans un cas de référence où l'on imposera un courant de 15 kA, le champ sera de 75 mT.

1.3. Calcul des propriétés

Dans ce modèle, toutes les propriétés des ions et des électrons sont directement calculées durant
la simulation à partir des équations provenant de l'ouvrage de Braginskii [22]. Les coefficients
des équations dépendent de la valeur de la charge moyenne du plasma Z=1.9 [26]. Les temps
de relaxation entre deux collisions électrons-électrons τ<sub>e-e</sub> et des ions-ions τ<sub>i-i</sub> sont
respectivement données par ces relations :

$$
\tau_{e-e} = \tau_e = \frac{12\pi\epsilon_0^2 \sqrt{\frac{\pi}{2}m_e k_B T_e k_B T_e}}{n_e e^4 Z \log(\Lambda)} \quad (1)
$$

$$
\tau_{i-i} = \tau_i = \frac{12\pi\epsilon_0^2 \sqrt{\frac{\pi}{2}m_i k_B T_i k_B T_i}}{n_e e^4 Z^4 \log(\Lambda)} \quad (2)
$$

On retrouve : $\epsilon_0$ la permittivité relative du vide, $k_B$ la constante de Boltzmann et $e$ la charge élémentaire d'un proton. Le terme $\log(\Lambda)$ est le logarithme Coulombien qui est donné en cgs par Braginskii :

$$
\log(\Lambda) = 23,4 - 1,15\log_{10}(n_e \times 10^{-6}) + 3,35\log_{10}(T_e(\text{eV})) \quad (3)
$$

Braginskii donne également une formulation permettant de considérer l'anisotropie de la conductivité électrique en fonction des directions (parallèle ou perpendiculaire au champ B<sub>AMF</sub>) et en fonction de la valeur du B<sub>AMF</sub>. En effet, suivant la valeur du champ magnétique externe, les trajectoires des électrons peuvent être modifiées entraînant des changements des composantes de la conductivité électrique. Ce phénomène est également pris en compte dans nos développements par les relations suivantes.

$$
\sigma_{\parallel} = \frac{\sigma_0}{\alpha_0} \tag{4}
$$

$$
\sigma_{\perp} = \frac{\sigma_0}{a} \tag{5}
$$

$$
\sigma_0 = \frac{e^2 n_e \tau_e}{m_e} \qquad (6)
$$

$$
a = \left( 1 - \frac{\alpha'_{1} \cdot x^2 + \alpha'_{0}}{x^4 + \delta_{1} \cdot x^2 + \delta_{0}} \right) \qquad (7)
$$

$$
x = \omega_e \tau_e \tag{8}
$$

$$
\omega_e = \frac{eB_{AMF}}{c m_e} \qquad (9)
$$

Avec les composantes suivantes : $\sigma_{\parallel}$ la conductivité électrique parallèle à l'axe du champ magnétique, $\sigma_{\perp}$ la conductivité électrique perpendiculaire à l'axe du champ magnétique, $\sigma_0$ est une constante, les coefficients de l'équation a : $\alpha_0, \alpha'_1, \alpha'_0, \delta_1, \delta_0$ dépendent de la charge Z, nous utilisons ceux calculés pour un Z proche de 1.9 et proviennent de l'ouvrage de Braginskii [22] et son rassemblés dans le tableau 1. Le terme $\omega_e$ est la pulsation cyclotronique des électrons et dépend de $B_{AMF}$ avec la constante c la vitesse de la lumière.
<table><caption>**Tableau 1 :** Valeur des coefficients utilisés pour le calcul de la conductivité anisotropique [22].</caption><thead><tr><th>Z</th><th>α<sub>0</sub></th><th>α'<sub>1</sub></th><th>α'<sub>0</sub></th><th>δ<sub>1</sub></th><th>δ<sub>0</sub></th></tr></thead><tbody><tr><td>2</td><td>0.4408</td><td>5.523</td><td>0.5956</td><td>10.80</td><td>1.0465</td></tr></tbody></table>

Le tenseur des contraintes visqueuses des électrons est défini par la relation suivante avec $u_e$ la vitesse des électrons et $L_H$ la longueur caractéristique du problème (distance inter-électrode de 1 cm) :

$$ \bar{\pi}_e = 0.73 n_e k_B T_e \tau_e \frac{u_e}{L_H} \qquad (10) $$

Le rapport du gradient du tenseur des contraintes visqueuses des électrons sur le gradient de pression des électrons $P_e$ est proportionnel au rapport du libre parcours moyen associé aux collisions électrons-électrons sur la longueur caractéristique. Ce libre parcours moyen étant très faible (environ $10^{-5}$ à $10^{-4}$ m), le terme de contraintes visqueuses associé aux électrons peut être négligé devant le gradient de pression des électrons dans l'équation de quantité de mouvement [1].

La conductivité thermique des ions et des électrons s'écrit :

$$ k_e = 4.9 \frac{n_e k_B^2 T_e \tau_e}{m_e} \qquad (11) $$

$$ k_i = 3.9 \frac{n_e k_B^2 T_i \tau_i}{m_i} \qquad (12) $$

La vitesse du son qui dépend de la température des ions et des électrons est donnée par :

$$ v_{son} = \sqrt{\frac{k_B}{m_i} (T_e + \gamma T_i)} \qquad (13) $$

On retrouve dans la relation précédente pour un gaz parfait monoatomique le coefficient $\gamma = \frac{c_P}{c_V} = \frac{5}{3}$. Il est à noter que pour la viscosité et la conductivité thermique des électrons et des ions, nous considérons seulement le calcul des termes de translation en accord avec la théorie de Braginskii. De ce fait, les termes de réaction chimique, en lien avec une charge moyenne fixe du plasma ne sont pas considérés.

## 1.4. Système d'équation

Les équations magnétohydrodynamiques qui décrivent l'arc dans cette configuration sont couplées avec les équations de Maxwell dans le cas du vide. La mise en œuvre d'un modèle numérique 3D décrivant le milieu plasma comme un fluide à deux températures permettra de prédire les champs de température, de pression, de vitesse pour les ions et les électrons. Pour l'instant, le modèle en place ne décrit que le régime subsonique. Pour comprendre l'approche et les équations établies pour obtenir les résultats suivants, il faut s'inspirer des travaux de Braginskii [22] et de E. Schade [14] et B. Tezenas [19] Y. Langlois [20]. Nous rappelons dans cette étude que nous sommes en stationnaire, aucun terme temporel n'apparaîtra dans les équations utilisées.
## Partie ions :

Les caractéristiques du fluide sur Ansys Fluent sont celles des ions. La vitesse calculée par fluent est la vitesse des ions $\vec{u}_i$ et elle permet d'écrire la loi de conservation de la masse suivante :

$$ \vec{\nabla}(\rho_i \vec{\nabla} \cdot \vec{u}_i) = 0 \qquad (14) $$

La densité de masse du fluide est celle des ions, car $m_i \gg m_e$ et s'écrit :

$$ \rho \approx \rho_i = n_i m_i = \frac{m_i P}{k_B (T_i + Z T_e)} \qquad (15) $$

La pression du fluide calculée par fluent est la pression totale $P$, c'est-à-dire celle des ions $P_i$ additionnée à celle des électrons $P_e$ :

$$ P = P_i + P_e = n_i k (T_i + Z T_e) \qquad (16) $$

On établit également le bilan de quantité de mouvement des ions suivant :

$$ \vec{\nabla}(\rho_i \cdot \vec{u}_i \cdot u_{ij}) = -\vec{\nabla}(P) - \vec{\nabla} \cdot \bar{\pi}_i + (\vec{j} \times \vec{B})_j \qquad (17) $$

Les résolutions des équations de l'énergie des ions et des électrons sont traitées en fonction de l'enthalpie. L'enthalpie massique des ions $h_i$ calculée par Fluent est donnée par la relation suivante :

$$ h_i = \frac{5}{2} \frac{n_i k_B T_i}{\rho_i} = \frac{5}{2} k_B \frac{T_i}{m_i} \qquad (18) $$

On définit l'équation de l'énergie des ions :

$$ \vec{\nabla}(\rho_i \vec{u}_i(h_i)) = \vec{u}_i \vec{\nabla}(P_i) + \vec{\nabla}(k_i \vec{\nabla} \cdot T_i) - \bar{\pi}_i : \vec{\nabla} \cdot \vec{u}_i + Q_i \qquad (19) $$

Les termes sources de l'équation de l'énergie sont respectivement le terme convectif, le terme diffusion, le tenseur des forces visqueuses des ions, $Q_i$ le terme d'échange d'énergie par collisions élastiques. On ne prend pas en compte les pertes d'énergie par rayonnement dans cette étude.

$$ Q_i = 3 \frac{m_e}{m_i} \frac{n_e}{\tau_e} k_B (T_e - T_i) \qquad (20) $$

## Partie électrons :

L'enthalpie massique des électrons $h_e$ est donnée par la relation suivante :

$$ h_e = \frac{5}{2} \frac{n_e k T_e}{\rho_i} = \frac{5}{2} k Z \frac{T_e}{m_i} \qquad (21) $$

On définit l'équation de l'énergie des électrons :

$$ \vec{\nabla}(\rho_i \vec{u}_i h_e) = \vec{u}_i \vec{\nabla}(P_e) + \vec{\nabla}(k_e \vec{\nabla} \cdot T_e) + \left( \frac{j_x^2 + j_y^2}{\sigma_{\perp}} \right) + \left( \frac{j_z^2}{\sigma_{\parallel}} \right) - Q_i + F_r \qquad (22) $$

Le terme source $F_r$ s'exprime suivant l'équation (23) où le coefficient $\beta = 2g_T$ et $g_T$ est une constante de Braginskii [22] :

$$ F_r = \frac{2\beta}{e} \vec{j} \cdot \vec{\nabla}(T_e) \qquad (23) $$
Finalement, nous définissons la vitesse des électrons à partir de celle des ions et de la densité de courant :

$$ \vec{u}_e = \vec{u}_i - \frac{\vec{j}}{\epsilon n_e} \qquad (24) $$

Sur la base de la loi d'Ohm généralisée [14], on peut définir la densité de courant par :

$$ \vec{j} = \sigma [-\vec{\nabla}V_{\text{eff}} + \vec{u}_e \times \vec{B}] \qquad (25) $$

Avec :

$$ -\vec{\nabla}V_{\text{eff}} = -\vec{\nabla}V + \frac{\vec{\nabla}P_e}{e n_e} \qquad (26) $$

La loi de conservation du courant impose que :

$$ \vec{\nabla}(\vec{j}) = 0 \qquad (27) $$

On obtient donc la condition suivante :

$$ \vec{\nabla}(\sigma\vec{\nabla}V_{\text{eff}}) - \vec{\nabla}(\sigma(\vec{u}_e \times \vec{B})) = 0 \qquad (28) $$

Le calcul du champ magnétique total $\vec{B}_{tot}$ sera la somme du champ magnétique auto-induit $\vec{B}_{self}$ et du champ magnétique externe $\vec{B}_{AMF}$.

$$ \vec{B}_{tot} = \vec{B}_{self} + \vec{B}_{AMF} \qquad (29) $$

Le $\vec{B}_{self}$ est obtenu en calculant le gradient du potentiel vecteur dans le domaine et en utilisant la loi de Biot & Savart pour fermer le système.

## 1.5. Conditions aux limites

* A la cathode, la zone ablatée représente 75% de la surface. Dans notre cas, en première approche, nous considérons que l'ablation est régie par un profil gaussien imposé à l'électrode entière présenté sur la figure 5. Il est réglé de manière à ce que la valeur sur le bord de la cathode soit proche de 10% de la valeur de l'axe.

<figure><img src="image_7.png" /><figcaption>**Figure 5** : Densité de courant à la cathode J (A/m²) obtenu par un profil de Gauss.</figcaption></figure>
* Comme suggéré par E. Schade [14], nous supposons que les températures des électrons et des ions sont imposées aux niveaux des électrodes.

* Enfin, la pression de sortie à l'anode correspondant à la pression totale est calculée en fonction de la vitesse des ions. Cette stratégie de calcul est issue de la thèse de Y. Langlois [1]. Tout d'abord, nous résolvons l'équation de conservation de l'énergie des ions dans la pré-gaine. Ensuite, nous calculons l'équation de densité ionique dans la pré-gaine par l'équation de Boltzmann. Enfin, nous établissons la loi de conservation de la masse à travers la pré-gaine. En reliant ces relations, nous pouvons déduire la pression de sortie en supposant que la vitesse de Böhm est constante dans la pré-gaine.

## 2. Résultat de la simulation

### 2.1. Logiciel de simulation numérique Ansys Fluent 23 R2

Nous utilisons Ansys Fluent version 23 R2 afin de modéliser le disjoncteur à arc sous vide (VACB en anglais). Ansys Fluent est un logiciel de simulation numérique de type CFD utilisant la méthode des volumes finis afin de modéliser des écoulements de fluides en résolvant les équations de Navier-Stokes. Dans notre étude, nous couplons ce logiciel avec des développements (en langage C) contenant des UDFs (User Defined Functions) et des UDSs (User Defined Scalars) spécifiques pour décrire la physique de l'arc électrique [27].

### 2.2. Maillage et représentation des résultats

Les électrodes (anode en haut et cathode en bas) et le domaine fluide, sont représentés avec un maillage fin comme le montre la figure 6 a), de 728000 cellules, avec un pas de 0.1 mm au voisinage des électrodes et un pas de 0.22 mm au milieu. Les calculs de simulation sont réalisés seulement dans la partie fluide et non dans les électrodes.

<figure><img src="image_8.png" /><figcaption>Figure 6 : Maillage utilisé dans la simulation</figcaption></figure>

a) vue en 3D b) dans le plan xz de Ansys Fluent.

Les résultats obtenus seront présentés sur un plan de coupe pour faciliter leur compréhension comme illustré sur la figure 6 b).

### 2.3. Arc diffus subsonique

Les résultats de simulation de l'arc sous vide en régime subsonique dans l'espace inter-électrodes sont présentés. Le cas utilisé dans l'article de E. Schade et D.L. Shmelev [14] est notre référence pour cette étude. Nous fixons une valeur de courant de I = 15 kA, un B<sub>AMF</sub> de 75 mT, un taux d'ablation de 50 µg/C, une charge fixe de Z = 1.9. Nous observons d'abord le champ de la densité de courant J<sub>z</sub> sur la figure 7.
<figure><img src="image_9.png" /><figcaption>**Figure 7 :** Champ de la densité de courant $J_z$ pour un courant de 15 kA, un champ magnétique appliqué $B_{AMF}$ de 75 mT.</figcaption></figure>

La valeur maximale de la densité de courant $J_z$ est observée au centre du plasma avec une valeur de $1.1 \times 10^7 A/m^2$. Sur la figure 8), on observe respectivement la température ionique $T_i$ et électronique $T_e$ dans l'espace inter-électrode.

<figure><img src="image_10.png" /><figcaption>**Figure 8 :** Champs de la température a) ionique $T_i$ en eV et b) électronique $T_e$ en eV pour un courant de 15 kA, un champ magnétique appliqué $B_{AMF}$ de 75 mT.</figcaption></figure>

Les champs sont symétriques et les ordres de grandeurs des températures sont sensiblement les mêmes. En comparaison avec la littérature et pour les mêmes conditions, E. Schade obtient également des température autour de 3.8eV [14]. La figure 9 représente la pression absolue $P$ dans l'ampoule.
<figure><img src="image_11.png" /><figcaption>**Figure 9** : Champ de la pression absolue P (mbar) pour un courant de 15 kA, un champ magnétique appliqué B<sub>AMF</sub> de 75 mT.</figcaption></figure>

La pression est maximale au centre du plasma là où la densité de courant est maximale. Les interactions entre particules sont fortes, la densité d'ions au centre est d'environ 1. 10<sup>22</sup> m<sup>-3</sup> et le flux de particules chargées est maximal. La figure suivante correspond à la vitesse moyenne des ions v<sub>i</sub> en m/s.

<figure><img src="image_12.png" /><figcaption>**Figure 10** : Champ de la vitesse moyenne des ions pour un courant de 15 kA, un champ magnétique appliqué B<sub>AMF</sub> de 75 mT.</figcaption></figure>

C'est la composante de la vitesse azimutale qui prédomine. En effet la force (j × B) de l'équation 17, provoque la rotation des ions autour de l'axe z. L'augmentation de la pression de la cathode à l'anode entraîne et accélère le flux de plasma vers l'anode. Nous remarquons aussi une augmentation lente de la vitesse des ions dans la direction radiale pour finalement s'accélérer drastiquement aux niveaux des parois, d'où une vitesse élevée de l'ordre de 3100 m/s sur les bords.

**2.4. Discussion**

Les résultats de notre modèle sont intéressants et en accord avec les résultats de la littérature. Ils nous permettent de bien appréhender le comportement du plasma et de réaliser à présent des études paramétriques. Néanmoins il est important de souligner les limites de notre modèle même si elles sont classiquement utilisées par la plupart des auteurs.

Nous considérons à la cathode un taux d'érosion classiquement utilisé dans la littérature fonction de la densité de courant [25]. Nous fixons un profil de densité de courant pour une intensité donnée. La densité de courant peut être représentée par une valeur constante sur le rayon [15] ou bien par un profil radial [1]. Cependant, durant la coupure, le courant varie de plusieurs kA, il faudrait considérer alors un profil de densité de courant évoluant en fonction du temps et du régime d'écoulement des ions (diffus ou concentré par exemple).
Les pertes d'énergies par rayonnement ne sont pas systématiquement prises en compte dans les modèles. La difficulté réside dans la mise en place des hypothèses de calculs permettant d'estimer le rayonnement d'un plasma à deux températures. Généralement les auteurs associent la contribution des raies aux ions et celle du continuum aux électrons. De plus, le calcul du rayonnement est incomplet avec le NIST car la base de données ne prend pas en compte les ions Cu<sup>+3</sup> et Cu<sup>+4</sup>. Ainsi devant ces lacunes nous ne considérerons pas dans notre modélisation la prise en compte du rayonnement. Ces pertes par rayonnement peuvent être cependant conséquentes, par exemple dans le cadre des plasmas thermiques, le rayonnement sur l'axe est du même ordre de grandeur que le chauffage par effet Joule.

Le champ magnétique extérieur dépend de la géométrie, de la nature des électrodes, et du courant qui les traverses. Dans notre simulation, nous considérons un champ B<sub>AMF</sub> constant dans tout l'espace inter-électrodes. Cependant, dans la littérature, certains auteurs considèrent ce champ variable radialement et le représente par différents types de profils [21] (en cloche ou en selle). Le choix du profil caractérisant le B<sub>AMF</sub> impacte fortement les résultats de la simulation. Si la variation du courant était prise en compte, nous devrions alors considérer un champ B<sub>AMF</sub> dépendant de la valeur de celle-ci.

Durant la phase d'arc, la composition du milieu évolue avec l'augmentation du courant et donc l'augmentation de la température et de la pression. Cependant, la grande majorité des auteurs fixe la valeur de la charge moyenne du plasma quel que soit le courant et les températures qui en découlent. Dans la section suivante une étude paramétrique sur l'impact de la variation de la charge Z est présentée.

## 2.5. Variation de la charge moyenne Z du plasma

La littérature considère par hypothèse que la charge globale du plasma est fixe et égale à 1.9 [26]. En réalité elle dépend de l'ablation de la cathode qui varie suivant l'intensité du courant. Nous proposons dans cette section, d'étudier l'impact de la variation de la charge moyenne du plasma notamment sur la température ionique. Afin d'estimer cette variation de charge nous utilisons des données issues de calculs de propriétés d'un plasma de cuivre à l'équilibre thermodynamique local (ETL) et hors ETL. Les calculs sont réalisés au sein de l'équipe AEPPT du Laplace à Toulouse. Nous calculons la variation de la composition d'un plasma de cuivre en fonction de la température du plasma pour une pression de 0.1 bar, pour plusieurs valeurs de θ :

$$ \theta = \frac{T_e}{T_i} \qquad (30) $$

Dans notre cas de référence nous retrouvons dans les champs de température figure 8 a) et b), des valeurs de θ qui varient entre 1.04 et 0.8. Il convient donc de considérer l'influence de ce déséquilibre sur la composition du plasma. Les ions de cuivre pris en compte dans ce calcul sont le Cu<sup>+</sup>, Cu<sup>2+</sup>, Cu<sup>3+</sup>, Cu<sup>4+</sup>, Cu<sup>-</sup>.

On observe sur la figure 11, la variation de la densité en m<sup>-3</sup> des particules en fonction de la température en eV pour P = 0.1 bar. Cette composition a été utilisée pour calculer la variation de la charge Z.
On remarque qu'à haute température, pour environ 3-4 eV, seuls les ions de cuivre Cu⁺, Cu²⁺, Cu³⁺, Cu⁴⁺, et les électrons restent dans des proportions élevées. A partir de ces données et celles hors ETL, nous pouvons calculer une estimation de la charge du plasma évoluant avec la température par la relation suivante :

$$Z = \frac{n_e}{n_i} = \frac{n_e}{n_{Cu+} + n_{Cu2+} + n_{Cu3+} + n_{Cu4+} + n_{Cu-}} \quad (31)$$

La figure 12 représente la variation de la charge moyenne en fonction de la température pour un plasma de cuivre pour un θ = 1 (plasma ETL), θ = 0.8 et θ = 1.04. Ces valeurs encadrent les valeurs de θ de la figure 8.

On remarque que la charge Z croit progressivement par paliers, relatifs aux différentes réactions d'ionisations du cuivre, jusqu'à environ Z = 4 pour une température de 4 eV (environ 50 000 K) quelles que soient les valeurs de θ et de la pression. Sur la figure 8, la plage de température ionique et électronique se trouve majoritairement entre 3 et 4 eV (zone grise sur la figure 12), ce qui ne correspond pas à une charge de Z = 1.9. Une valeur de Z = 4 serait donc plus appropriée. Cette estimation de la charge Z provient de calcul de composition qui est surtout légitime pour des pressions supérieures au dixième de bar. Quoi qu'il en soit la charge Z est loin d'être constante sur l'intervalle de température.

<figure><img src="image_13.png" /><figcaption>**Figure 11** : Densité d'espèces du cuivre en fonction de la température à l'ETL pour une pression de 0.1 bar.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>**Figure 12** : Variation de la charge moyenne du plasma Z en fonction de la température en eV pour P = 0.01 bar (θ=1, θ=0.8 et θ=1.04) et pour P = 0.1 bar (θ=1).</figcaption></figure>

Sur la figure 13, nous observons l'influence de la valeur de la charge du plasma Z sur les champs de température ionique pour un Z = 1.9 en a) et Z = 4 en b). Pour prendre en compte différentes valeurs de Z nous adaptons nos coefficients du tableau 1 à la valeur de Z considérée et les coefficients des équations 11 et 23, à l'aide des données proposées par Braginskii [22].

<figure><img src="image_15.png" /><figcaption>**Figure 13** Champs de la température ionique Tᵢ en eV a) pour Z = 1.9 et b) pour Z = 4 pour un courant de 15 kA, un champ magnétique appliqué B<sub>AMF</sub> de 75 mT.</figcaption></figure>
Nous remarquons que la température ionique est nettement plus élevée pour un Z = 4. En effet, l'augmentation de Z induit une augmentation de la densité d'ions positifs (voir figure 12) et donc une augmentation de la pression des ions (équation 16), et donc de leur température. Pour souligner son importance, la figure 14 représente les profils de température ionique pour un Z = 1.9 (en bleu), un Z = 3 (en vert) et un Z = 4 (en rouge) sur une ligne de coupe au centre du plasma.

<figure><img src="image_16.png" /><figcaption>**Figure 14** : Température ionique le long de ligne de coupure au centre du plasma pour Z=1, Z=0.8 et Z=1.4.</figcaption></figure>

On observe que la température ionique augmente avec l'augmentation de la charge Z. Une charge plus élevée induit donc un plasma plus ionisé et énergétique. Les collisions sont plus nombreuses et les échanges de chaleur sont plus élevés donc la température ionique augmente. On observe également que le terme d'échange d'énergie par collision élastiques $Q_i$ de l'équation 20 augmente pour Z = 4. Le nombre de collisions entre les ions de cuivre et les électrons augmentent, car la proportion d'ions $Cu^{4+}$ devient plus importante comme on peut le voir sur la figure 11. L'énergie transférée aux ions du fait de leurs collisions avec les électrons est donc plus importante.

## Conclusion

Un modèle hors équilibre thermique, basé sur les travaux de la littérature a été mis en place sous Fluent pour décrire un arc électrique présent dans les disjoncteurs à vide moyenne tension. Dans un mode d'arc subsonique et stationnaire, les températures, les vitesses et la pression absolue nous permettent de mieux appréhender les phénomènes physiques présents au sein du plasma dans le disjoncteur à vide. Cependant il reste encore des hypothèses et des approximations sur lesquelles revenir pour améliorer notre modèle numérique. La variation de la composition du milieu et donc de la charge Z influence les caractéristiques de l'arc. Nous avons calculé la composition d'un plasma de cuivre à l'ETL et hors ETL à P = 0.01 bar en fonction de la température et pour plusieurs valeurs de θ. Cela nous a permis d'estimer une valeur de Z qui coïncide avec nos températures trouvées dans la simulation. En prenant en compte un Z = 4 correspondant aux conditions physiques pour un courant de 15 kA, et un B$_{AMF}$ de 75 mT, on observe une augmentation significative de la température ionique et des grandeurs associées.
Références

[1] Y. Langlois, (2010), « Modélisation de l'arc électrique dans un disjoncteur à vide », thèse, Institut Jean Lamour, Nancy.
https://hal.univ-lorraine.fr/tel-01748841v1

[2] A. Palle, (2016), « Evolution dynamics of European electricity transmission networks: Technical networks between growth and degrowth », *Flux*, vol. 105, n° 3, p. 18-32.
https://shs.cairn.info/journal-flux-2016-3-page-18?lang=en

[3] P. Picot, (2000), "La coupure du courant électrique dans le vide" *Cahier technique n° 198*.
http://www.schneider-electric.com

[4] Yuqing Pei, (2014), « Computer Simulation of Fundamental Processes in High Voltage Circuit Breakers Based on an Automated Modelling Platform », The University of Liverpool (United Kingdom), Theses.

[5] GEIC, (2021), « Technical Summary, IPCC-ARG-WGI FullReport 2021 », https://doi.org/10.1017/9781009157896.002.

[6] « Site web Siemens Energy ». [En ligne]. Disponible sur : https://www.siemens-energy.com/global/en/home/products-services/product-offerings/circuit-breakers.html

[7] Siemens Energy Global GmbH & Co. KG, (2022), « High-Voltage Circuit Breakers - From 72.5 kV up to 1100 kV ».

[8] J. V. R. Heberlein et J. G. Gorman, (1980), « The High Current Metal Vapor Arc Column Between Separating Electrodes », https://doi.org/10.1109/TPS.1980.4317328

[9] B. Tezenas Du Moncel, (2021), « Étude expérimentale et modélisation d'un plasma de vapeurs métalliques dans une ampoule à vide pour le développement de disjoncteurs Haute-Tension », thèse, Institut Jean Lamour, Nancy.
https://hal.univ-lorraine.fr/tel-01925198v1

[10] E. Hantzspe, (2003), « Mysteries of the arc cathode spot: A retrospective glance », *IEEE Transactions on Plasma Science*, vol. 31, no. 5, pp. 799-808.
https://doi.org/10.1109/TPS.2003.818412

[11] Burkhard Jüttner, (2001), « Cathode spots of electric arcs », J. Phys. D: Appl. Phys. 34 R103
https://doi.org/10.1088/0022-3727/34/17/202

[12] C. W. Kimblin, (1973), « Erosion and Ionization in the Cathode Spot Region of Vacuum Arcs », *J. Appl. Phys*, vol. 44, n° 7, p. 3074-3080.
https://doi.org/10.1063/1.1662710

[13] R. L. Boxman, (1977), « Magnetic constriction effects in high-current vacuum arcs prior to the release of anode vapor », *Journal of Applied Physics*, vol. 48, n° 6, p. 2338-2345.
https://doi.org/10.1063/1.323994

[14] E. Schade et D. L. Shmelev, (2003), « Numerical simulation of high-current vacuum arcs with an external axial magnetic field », *IEEE Transactions on Plasma Science*, vol. 31, n° 5 I, p. 890-901.
https://doi.org/10.1109/TPS.2003.818436.

[15] E. Schade, (2005), « Physics of high-current interruption of vacuum circuit breakers », *IEEE Transactions on Plasma Science*, vol. 33, n° 5 I, p. 1564-1575
https://doi.org/10.1109/TPS.2005.856530

[16] M. Homma, M. Sakaki, E. Kaneko, et S. Yanabu, (2006), « History of vacuum circuit breakers and recent developments in Japan », *IEEE Transactions on Dielectrics and Electrical Insulation*, p. 85-92.
https://doi.org/10.1109/TDEI.2006.1593405
[17] W. Hartmann, A. Lawall, R. Renz, M. Römheld, N. Wenzel, et W. Wietzorek, (2011), « Cathode spot dynamics and arc structure in a dense axial magnetic-field-stabilized vacuum arc », *IEEE Transactions on Plasma Science*, p. 1324-1329
https://doi.org/10.1109/TPS.2011.2140402

[18] D. L. Shmelev et I. V. Uimanov, (2015), « Hybrid computational model of high-current vacuum arcs with external axial magnetic field », *IEEE Transactions on Plasma Science*, vol. 43, n° 8, p. 2261-2266
https://doi.org/10.1109/TPS.2015.2430372.

[19] B. Tezenas Du Montcel, P. Chapelle, C. Creusot, et A. Jardy, (2019), « Numerical Study of the Current Constriction in a Vacuum Arc at Large Contact Gap », *IEEE Transactions on Plasma Science*, vol. 47, n° 5, p. 2765-2774
https://doi.org/10.1109/TPS.2019.2909964.

[20] Y. Langlois, P. Chapelle, A. Jardy, et F. Gentils, (2011), « On the numerical simulation of the diffuse arc in a vacuum interrupter », *Journal of Applied Physics*, vol. 109, n° 11, https://doi.org/10.1063/1.3587180

[21] L. Wang, S. Jia, L. Zhang, D. Yang, Z. Shi, F. Gentils, B. Jusselin, (2008), « Current constriction of high-current vacuum arc in vacuum interrupters », *Journal of Applied Physics*, vol. 103, n° 6, 063301
https://doi.org/10.1063/1.2875813

[22] S. I. Braginskii, (1965), "Transport processes in a plasma", *Reviews of Plasma Physics*, p 205-311.

[23] J. Kutzner et H. C. Miller, (1989), « Ion Flux From the Cathode Region of a Vacuum Arc », *IEEE Transactions on Plasma Science*, vol. 17, n° 5, p. 688-694
https://doi.org/10.1109/27.41183

[24] J. Kutzner et H. C. Miller, (1992) « Integrated ion flux emitted from the cathode spot region of a diffuse vacuum arc », J. Phys. D: Appl. Phys., Vol 25, p 686
https://doi.org/10.1088/0022-3727/25/4/015

[25] J. E. Daalder, (1976), « Components of cathode erosion in vacuum arcs », J. Phys. D: Appl. Phys., vol. 9, n° 16, p. 2379-2395
https://doi.org/10.1088/0022-3727/9/16/009

[26] A. Anders, (1997), « Ion charge state distributions of vacuum arc plasmas: The origin of species », Phys. Rev. E 55, p 969-981
https://doi.org/10.1103/PhysRevE.55.969

[27] ANSYS, (2023), Ansys Fluent Theory Guide Release 2023 R2.