# Application de la méthode de Boltzmann sur réseau aux simulations instationnaires en mécanique des fluides

Simon Marié

► To cite this version:

Simon Marié. Application de la méthode de Boltzmann sur réseau aux simulations instationnaires en mécanique des fluides. Analyse numérique [math.NA]. Conservation National des Arts et Métiers, 2021. tel-04077809

HAL Id: tel-04077809

https://hal.science/tel-04077809

Submitted on 21 Apr 2023

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
CONSERVATOIRE NATIONAL DES ARTS ET MÉTIERS

ÉCOLE DOCTORALE SCIENCES DES MÉTIERS DE L'INGÉNIEUR

MÉMOIRE SCIENTIFIQUE PROVISOIRE

présenté par

SIMON MARIÉ

en vue d'obtenir :
L'HABILITATION À DIRIGER DES RECHERCHES

Application de la méthode de Boltzmann sur réseau aux
simulations instationnaires en mécanique des fluides.

LE 17 NOVEMBRE 2021

COMPOSITION DU JURY :

GARANT

Pr. Iraj MORTAZAVI Professeur des Universités au CNAM

RAPPORTEURS

Pr. Eric LAMBALLAIS Professeur à l'Université de Poitiers

Pr. Emmanuel LÉVÈQUE Directeur de Recherche CNRS au LMFA

Pr. Florian DEVUYST Professeur à l'Université Technologique de Compiègne

EXAMINATEURS

Pr. Maria Vittoria SALVETTI Professeur à l'Université de Pise (Italie)

Pr. Pierre SAGAUT Professeur à l'Université d'Aix Marseille

Pr. François DUBOIS Professeur des Universités au CNAM

Dr. Sébastien DECK Directeur de Recherche ONERA

Pr. Xavier GLOERFELT Professeur des Universités à l'ENSAM Paris (Membre invité)
<figure><img src="image_2.png" /></figure>
## II Etude des écoulements complexes en LBM 41

### 3 Etude LBM d'une interface à fort gradient de densité 43

3.1 Le modèle d'interface de He-Chen-Zhang 43

3.2 Filtrage adaptatif pour les forts gradients de densité 45

3.3 Application à l'instabilité de Rayleigh-Taylor 45

### 4 Modélisation de la diffusion multi-espèces en LBM 49

4.1 Modélisation LBM multi-espèces 49

4.1.1 Position du problème 49

4.1.2 Utilisation d'un terme de forçage 50

4.1.3 Validation du modèle sur des cas tests diffusifs et convectifs 51

4.2 Application à l'instabilité de Saffman-Taylor 52

4.2.1 Modélisation d'un milieu poreux 52

4.2.2 Effet du nombre de Péclet 53

4.2.3 Instabilité à trois espèces 53

4.2.4 Simulations à l'échelle d'un pore 55

## Références 61

## III Etudes aérodynamiques instationnaires 63

### 5 Etude des écoulements d'arrière corps 65

5.1 Etude des efforts instationnaires sur géométrie simplifiée 65

5.1.1 Reconstruction des efforts instationnaires 66

5.1.2 Reconstruction optimisée 67

5.2 Etude des configurations réelles 68

### 6 Etude des écoulements à thermodynamique complexe 71

6.1 Transition d'une couche de mélange supersonique 71

6.2 Modélisation des fluctuations de concentration 72

6.2.1 Validation sur un jet coaxial supersonique 73

6.2.2 Simulation d'une couche de mélange réactive 74

## Références 79

## Conclusion et Perspectives 83

## Suppléments 91

CV 91
# Présentation Générale
<figure><img src="image_3.png" /></figure>
# Présentation

## Parcours scientifique

Si l'on considère que mon parcours scientifique a débuté lors de mes années de master (DEA) à Poitiers, alors c'est incontestablement Yves GERVAIS qui a tout déclenché. En effet, c'est à l'issue de mon master en aérodynamique et aéroacoustique que l'évidence de poursuivre mon projet professionnel dans l'enseignement et la recherche s'est imposé. Après un stage au Centre National d'Etudes Spatiales (CNES) sur l'étude des écoulements d'arrière-corps, j'ai eu l'opportunité de poursuivre en thèse avec Pierre SAGAUT et Denis RICOT sur les thématiques des méthodes de Boltzmann sur Réseau qui commençaient juste à prendre de l'ampleur sur le plan industriel. Cette oscillation initiale entre deux thèmes de recherche relativement différents ne sera pas immédiatement amortie et occupera une bonne place dans ce mémoire comme nous allons le voir au cours de ces pages.

Au cours de ma thèse chez RENAULT, j'ai été plongé au coeur des problématiques industrielles à travers le projet MIMOSA qui visait à proposer des méthodes innovantes pour la modélisation des sources aéroacoustiques en collaboration avec ALSTOM, la SNCF et l'ECL. Le fil conducteur de ma thèse était de pouvoir caractériser de façon rigoureuse les capacités acoustiques du schéma LBM qui ne jouissait pas encore de la même littérature que les schémas Navier-Stokes sur ce sujet. Avec Denis RICOT et Pierre SAGAUT nous avons pu alors développer une méthode d'analyse basée sur la stabilité linaire pour mieux comprendre les apétences remarquables du schéma pour propager une information acoustique. Ces résultats se sont traduits également sur le plan industriel par la prise de conscience qu'il fallait une alternative au code commercial PowerFLOW qui était à l'époque le seul code disponible sur le marché utilisant la LBM. Quelques années après ma thèse, un consortium d'industriels mené par RENAULT a donc créé le code ProLB qui s'appelait à l'époque LABS et qui fait partie aujourd'hui des outils importants dans le secteur industriel de la LBM.

Après ma thèse j'ai eu l'opportunité d'approfondir les éléments que j'avais commencés lors de mon stage de fin d'études au CNES sur les instationnarités d'arrière-corps. J'ai donc pu retourner dans le milieu spatial à travers un Post-Doctorat à l'ONERA avec Sébastien DECK qui m'a beaucoup apporté sur le plan scientifique. L'idée était d'identifier des phénomènes locaux sur une marche axisymétrique permettant d'expliquer les efforts globaux mesurés expérimentalement. Cette identification permettait une réduction importante du nombre de capteurs nécessaire à la captation du phénomène. Puis, lors d'un post-doctorat au CNES
avec Hadrien LAMBARÉ, j'ai pu appliquer ces méthodes sur la géométrie réelle du lanceur Ariane5. J'ai ainsi pu avoir la chance de participer aux campagnes d'essais expérimentales d'un projet européen mené par le CNES, l'ESA, l'ONERA, TU-DELFT et ASTRIUM (Aujourd'hui ArianeGroup). L'idée du projet était de mieux comprendre les effets de géométrie réelle sur le buffeting de culot qui engendrait des efforts importants sur la tuyère du moteur Vulcain-2. Nous verrons dans ce mémoire que l'interaction entre le calcul et les essais a été très importante dans ce projet.

Avant d'arriver au Conservatoire des Arts et Métiers en 2011, j'étais déjà guidé par une envie forte de transmettre les connaissances et j'ai pu m'exercer à ce travail chez MECAVENIR puis à Paris 6 lors de mon dernier Post-doc à l'institut Jean le Rond d'Alembert en 2011 pour poursuivre mes activités sur les arrières-corps notamment avec Philippe DRUAULT qui m'a sensibilisé aux approches POD et DMD. J'ai ensuite eu l'opportunité de rejoindre le Conservatoire National des Arts et métiers (le CNAM) sur un poste de Maître de Conférence au laboratoire DynFLuid en cotutelle (ENSAM/CNAM). J'ai alors pu travailler avec Francesco GRASSO sur les thématiques de gaz chauds afin de poursuivre mes activités sur les arrières-corps qui n'avaient encore pas fait l'objet d'études en présence de gaz chauds propulsifs. Après une collaboration avec l'université de Rome La Sapienza pour essayer de comprendre les phénomènes intervenant dans la zone de transition d'une couche de mélange issue de deux couches-limites chaudes, nous avons lancer une thèse sur la modélisation de ce type d'écoulement en collaborant avec l'équipe de Pierre BRENNER chez ArianeGroup. L'idée était d'implémenter une modélisation statistique des fluctuations de concentration dans les outils industriels développés chez ArianeGroup pour permettre l'étude des écoulements d'arrières-corps en présence de jets chauds. Cette ambition a été un peu réorientée en raison des difficultés de prise en main de l'outil industriel. Peu après 2015, l'idée de revenir à la LBM s'est petit à petit imposée. J'ai alors pu améliorer les techniques de filtrage proposées pendant ma thèse en les rendant adaptatives lors d'une collaboration avec Xavier GLOERFELT. Par la suite, j'ai alors naturellement essayé d'apprivoiser la modélisation des écoulements multi-espèces en utilisant la LBM. J'ai alors pu lancer une thèse sur le sujet et j'ai pu développer un code me permettant aujourd'hui de poursuivre mes activités de façon autonome. Aujourd'hui une collaboration avec Iraj MORTAZAVI et Chloé MIMEAU a permis une comparaison des schémas LBM avec les méthodes Vortex et une collaboration avec l'ONERA est en cours sur les thématiques de couplage LBM/Navier-Stokes.

## Organisation du mémoire

Ce mémoire est composé de trois grandes parties, reprenant les différentes méthématiques de recherche présentées plus haut. La première partie présente les grandes lignes de la méthode de Boltzmann sur réseau ainsi que les avancées récentes sur le sujet. Puis cette première partie présentera mes activités de recherche liées à la stabilisation de la méthode par filtrage, à la comparaison de la LBM avec les autres approches ainsi qu'au couplage des méthodes LBM avec les solveurs
Navier-Stokes. Le fil conducteur de cette première partie sera principalement fo-
calisé sur la reconstruction des grandeurs macroscopiques à partir des grandeurs
mésoscopiques.

La deuxième partie sera dédiée aux applications de la LBM à la simulation des écoulements complexes. Nous aborderons d'abord les aspects diphasiques à travers l'application des filtres adaptatifs à la simulation de l'instabilité de Rayleigh-Taylor. Nous verrons ensuite comment les approches LBM peuvent décrire les mélanges et les interactions entre espèces lorsque la diffusion et la convection s'affrontent. Là encore, cette partie sera guidée par les relations entre les effets microscopiques locaux et les caractéristiques globales d'un mélange.

La troisième partie du mémoire traite des sillages instationnaires et des activités liées à l'aérodynamique instationnaire en général et aux écoulements d'arrière-corps. En particulier, nous présenterons les méthodes de reconstruction d'efforts sur des géométries plus ou moins simplifiées à partir de grandeurs locales et nous verrons les effets de la géométrie sur la dynamique des écoulements d'arrière-corps. Nous verrons également comment une dynamique de couche limite peut localement transitionner vers une dynamique de couche de mélange. Nous verrons enfin une tentative de construction d'approche hybride RANS/LES dédiée aux écoulements à plusieurs espèces. En particulier, il sera question des effets de diffusion massique à travers un nombre de Schmidt variable. Le fil conducteur de cette partie sera encore une fois piloté par les relations entre les effets de petites-échelles locales et de structure macroscopiques globales.
# Première partie

## La méthode de Boltzmann sur
Réseau : Principes et Comparaisons
<figure><img src="image_4.png" /></figure>
# Chapitre 1

## La méthode Boltzmann sur Réseau : Méthode Numérique en évolution

L'objectif de ce chapitre est de faire une présentation rapide des éléments importants de la méthode de Boltzmann sur réseau qui occupe une place importante dans mes travaux de recherche. L'idée est de présenter les aspects qui la caractérisent le plus, ainsi que de passer en revue les principales améliorations qui ont été apportées au schéma initial dans les dernières années. Je me bornerai dans ce chapitre aux écoulements chimiquement inertes et constitués d'une seule espèce. Les aspects multi-espèces seront abordés dans le chapitre 4. Le début du chapitre est consacré au rappels fondamentaux de physique statistique qui serviront à mieux comparer les modèles entre eux dans la partie numérique. Le lecteur désireux d'approfondir ces points exposés ici de façon sommaire pourra se référer aux nombreux ouvrages sur le sujet [64, 34, 29, 58].

### 1.1 Construction d'un réseau de vitesse

#### 1.1.1 L'équation de Boltzmann continue

Le premier élément important des méthodes LBM est le lien qui les connecte à la physique statistique dont l'élément descriptif fondamental est la fonction de distribution $f(x, c, t)$ qui décrit la densité de particules ayant une vitesse $c$ et présentent au point $x$ et au temps $t$. D'un point de vue descriptif, cela revient à décrire le comportement d'un fluide à l'aide de ses propriétés locales. A l'échelle microscopique on peut utiliser une approche Eulérienne pour décrire la trajetoire de chaque particule et reconstituer une dynamique à l'aide de moyennes d'ensembles. Cette approche est parfois appliquée dans les méthodes SPH ou les méthodes Vortex qui feront l'objet d'une comparaison avec la LBM dans la section 2.4. A l'échelle mésoscopique, on peut réaliser une moyenne d'ensemble sur un volume élémentaire en séparant les différentes vitesses admises par les particules pour aboutir à la définition de $f(x, c, t)$. Dès lors, afin de décrire la dynamique de ces échelles mésoscopiques, on considère une équation de transport standard de type transport/collision en écrivant :
$$ \frac{\partial f}{\partial t} + c_i \frac{\partial f}{\partial x_i} + \frac{F_i}{m} \frac{\partial f}{\partial c_i} = \left( \frac{\partial f}{\partial t} \right)_{\text{coll}} \qquad (1.1) $$

Cette équation de transport est appelée équation de Boltzmann [8, 11] et fut proposée par Ludwig Boltzmann en 1872. Dans cette équation $F/m$ désigne les forces extérieures exercées sur les particules par unité de masse. Nous verrons que ce terme de force s'avèrera très utile pour modéliser la présence de conditions aux limites complexes, des effets d'interaction moléculaire ou encore des tensions d'interfaces (Cf chapitre 4).

Pour retrouver l'échelle macroscopique qui est la plus souvent utilisée pour décrire la dynamique des fluides, on peut donc effectuer une intégration sur l'ensemble des valeurs de vitesse. Cet ensemble est généralement pris dans $\mathbb{R}^3$ bien que les vitesses admissibles par la physique restent bornées :

$$ \rho = \int_{\mathbb{R}^3} f d\mathbf{c} \qquad (1.2) $$

$$ \rho \mathbf{u} = \int_{\mathbb{R}^3} \mathbf{c} f d\mathbf{c} \qquad (1.3) $$

$$ \rho e + \frac{1}{2} \rho |\mathbf{u}|^2 = \frac{1}{2} \int_{\mathbb{R}^3} |\mathbf{c}|^2 f d\mathbf{c} \qquad (1.4) $$

Ces formulations sont appelées les moments de la fonction de distribution. La méthode de Boltzmann sur réseau consiste à trouver un ensemble discret de vitesse sur lequel les équations (1.2) (1.3) et (1.4) sont toujours valables. Le calcul des moments est donc un aspect important de la méthode, nous y reviendrons par la suite.

L'opérateur de collision $\left(\frac{\partial f}{\partial t}\right)_{\text{coll}}$ décrit l'effet des collisions entre particules et détermine donc le type de physique que l'on souhaite étudier. En effet, la description des collisions entre particules est intimement liée à celle de l'équilibre qui correspond à l'annulation du terme de collision. Ainsi, le choix d'un opérateur de collision implique très souvent un état d'équilibre dont quelques exemples sont donnés dans le tableau 1.1.

Les formulations d'équilibre sont établies en recherchant les distributions qui annulent l'opérateur de collision et en explicitant les invariants de collision, c'est à dire les moments qui sont conservés pendant la collision. En mécanique classique, l'élasticité des collisions implique la conservation des moments (1.2) (1.3) et (1.4) et engendre alors la distribution d'équilibre de Maxwell-Boltzmann :

$$ f^{\text{eq}}(\mathbf{x}, \mathbf{c}, t) = \rho \left(\frac{m}{2\pi k_B T}\right)^{3/2} \exp \left[-\frac{m}{2k_B T} (\mathbf{c} - \mathbf{u})^2\right] \qquad (1.5) $$

où $k_B$ est la constante de Boltzmann et $T$ la température. Pour définir un opérateur de collision simple, on peut constater que la quantité :

$$ \mathcal{H} = \int_{\mathbb{R}^3} f \ln(f) d\mathbf{c} \qquad (1.6) $$
<table><caption>TABLE 1.1 – *Quelques distributions d'équilibre*</caption><thead><tr><th>Type/Application</th><th>Formulation</th></tr></thead><tbody><tr><td>Maxwell-Boltzmann (Hydrodynamique)</td><td>ρ <span style="font-size: 2.5em; vertical-align: middle;">(</span><div style="display: inline-block; vertical-align: middle; text-align: center;"><span style="display: block; border-bottom: 1px solid black; padding: 0 0.2em;">m</span><span style="display: block; padding: 0 0.2em;">2πk<sub>B</sub>T</span></div><span style="font-size: 2.5em; vertical-align: middle;">)</span><sup>3/2</span> exp <span style="font-size: 2.5em; vertical-align: middle;">[</span><div style="display: inline-block; vertical-align: middle; text-align: center;"><span style="display: block; border-bottom: 1px solid black; padding: 0 0.2em;">m</span><span style="display: block; padding: 0 0.2em;">2k<sub>B</sub>T</span></div>(<b>c</b> - <b>u</b>)<sup>2</sup> <span style="font-size: 2.5em; vertical-align: middle;">]</span></td></tr><tr><td>Maxwell-Jutner (Méca relativiste) [51]</td><td><div style="display: inline-block; vertical-align: middle; text-align: center;"><span style="display: block; border-bottom: 1px solid black; padding: 0 0.2em;">γ<sup>2</sup>mvc</span><span style="display: block; padding: 0 0.2em;">k<sub>b</sub>TK<sub>2</sub></span></div> exp <span style="font-size: 2.5em; vertical-align: middle;">(</span><div style="display: inline-block; vertical-align: middle; text-align: center;"><span style="display: block; border-bottom: 1px solid black; padding: 0 0.2em;">mc<sup>2</sup></span><span style="display: block; padding: 0 0.2em;">k<sub>b</sub>T</span></div><span style="font-size: 2.5em; vertical-align: middle;">)</span></td></tr><tr><td>Bose-Einstein (Méca quantique) [65]</td><td><div style="display: inline-block; vertical-align: middle; text-align: center;"><span style="display: block; border-bottom: 1px solid black; padding: 0 0.2em;">Ξ</span><span style="display: block; padding: 0 0.2em;">exp <span style="font-size: 2.5em; vertical-align: middle;">(</span><div style="display: inline-block; vertical-align: middle; text-align: center;"><span style="display: block; border-bottom: 1px solid black; padding: 0 0.2em;">E<sub>i</sub> - μ<sub>i</sub></span><span style="display: block; padding: 0 0.2em;">k<sub>b</sub>T</span></div><span style="font-size: 2.5em; vertical-align: middle;">)</span> - 1</span></div></td></tr></tbody></table>

est une quantité décroissante. Ce théorème est connu sous le nom de théorème *H* de Boltzmann. Une version détaillée du théorème est exposée dans [62]. Ce théorème est à la base de la construction de l'opérateur BGK (Bhatnagar, Gross et Krook [6]) qui est très largement utilisé dans la plupart des modèles LBM. Cet opérateur considère que les collisions entre particules ont pour effet de modifier l'état d'équilibre 1.5 puis de relaxer vers cet état d'équilibre en un temps caractéristique donné $\lambda$, l'opérateur de collision peut prendre la forme suivante :

$$ \left(\frac{\partial f}{\partial t}\right)_{\text{coll}} = -\frac{1}{\lambda}[f - f^{\text{eq}}] \qquad (1.7) $$

Dans la suite, nous considérerons l'équation de Boltzmann munie de l'opérateur de collision BGK :

$$ \frac{\partial f}{\partial t} + c_i \frac{\partial f}{\partial x_i} = -\frac{1}{\lambda} [f - f^{eq}] \qquad (1.8) $$

L'équation (1.8), est une équation exacte et continue utilisant des variables statistiques. Il est important de noter ici que c'est une EDP linéaire comportant une partie convective pouvant être intégrée le long des vitesses caractéristiques *c*. Pour faire apparaître explicitement les non-linéarités il faut calculer les moments 1.2, 1.3, 1.4 en remplaçant *f* par l'équation (1.8) :

$$ \frac{\partial \rho}{\partial t} + \frac{\partial \rho u_i}{\partial x_i} = 0 \qquad (1.9) $$

$$ \frac{\partial \rho u_j}{\partial t} + \frac{\partial(\rho u_i u_j + P_{ij})}{\partial x_i} = 0 \qquad (1.10) $$

$$ \frac{\partial E}{\partial t} + \frac{\partial(u_i E + u_j P_{ij} + q_i)}{\partial x_i} = 0 \qquad (1.11) $$

L'équation (1.9) exprime la conservation de la masse, l'équation (1.10) la conservation de la quantité de mouvement dans laquelle $P_{ij}$ représente le tenseur des efforts :
$$P_{ij} = \int_{\mathbb{R}^3} (c_i - u_i)(c_j - u_j) f d\mathbf{c} \qquad (1.12)$$

et l'équation (1.11) traduit la conservation de l'énergie totale $E = \rho(\frac{1}{2}u^2 + e)$ dans laquelle $e$ représente l'énergie interne du gaz, définie par $e = \frac{D}{2}rT$ où $D$ est le nombre de degrés de liberté des particules, $r$ la constante réduite du gaz, $T$ la température et $q_i$ représente le flux de chaleur :

$$q_i = \frac{1}{2} \int_{\mathbb{R}^3} (c_i - u_i) |\mathbf{c} - \mathbf{u}|^2 f d\mathbf{c} \qquad (1.13)$$

Il est à noter que les expressions de $P_{ij}$ et $q_i$ dépendent ici de $f$. Pour lever cette indétermination, on utilise un développement asymptotique de la fonction de distribution selon le principe de Chapmann et Enskog [12]. L'idée réside dans le développement des fonctions de distribution en fonction du nombre de Knudsen. Ce nombre sans dimension, fondamental en physique statistique, traduit le rapport entre le libre parcours moyen $\mathcal{L}$ des particules, c'est-à-dire la distance sur laquelle ces particules ne subissent aucune collision, et la longueur caractéristique $L$ de l'écoulement :

$$\epsilon = \frac{\mathcal{L}}{L} = \frac{\lambda_0 \sqrt{rT_0}}{L} \qquad (1.14)$$

où $\sqrt{rT_0}$ a la dimension d'une vitesse et correspond à la vitesse d'agitation moyenne des particules et $\lambda_0$ représente le temps de relaxation moyen entre deux collisions. Dans toute la suite, les grandeurs seront réduites à leur valeur adimensionnelle sous la forme :

$$\tilde{t} = t \frac{\sqrt{rT_0}}{L} \qquad \tilde{x} = \frac{x}{L} \qquad \tilde{\lambda} = \frac{\lambda}{\lambda_0} \qquad \tilde{rT} = \frac{T}{T_0} \qquad \tilde{c} = \frac{c}{\sqrt{rT_0}} \qquad \tilde{f} = \frac{f}{\rho_0}$$

Ainsi, l'équation (1.8) prendra la forme :

$$\frac{\partial \tilde{f}}{\partial \tilde{t}} + \tilde{c}_i \frac{\partial \tilde{f}}{\partial \tilde{x}_i} = - \frac{1}{\epsilon \tilde{\lambda}} [\tilde{f} - \tilde{f}^{eq}] \qquad (1.15)$$

Dans la suite, pour des raisons de clarté, les variables adimensionnelles seront notées sans le tilde. Il est toutefois intéressant de remarquer que la forme sans dimension fait intervenir le nombre de Knudsen devant le temps de relaxation. Cela est dû au fait que le rapport $\sqrt{rT_0}/L$ apparaît lors de la normalisation.

Le développement de Chapman-Enskog, permet de trouver les expressions explicites des relations (1.12) et (1.13) en développant les fonctions de distribution de la façon suivante :

$$f = f^{(0)} + \epsilon f^{(1)} + \epsilon^2 f^{(2)} + \dots \qquad (1.16)$$

En injectant ce développement dans l'équation (1.15), on peut alors identifier les différents termes en fonction de leur ordre en $\epsilon$ pour obtenir :
$$
\begin{align}
\epsilon^{-1}: \quad & f^{(0)} &&= f^{eq} \nonumber \\
\epsilon^0: \quad & \frac{\partial f^{(0)}}{\partial t} + c_i \frac{\partial f^{(0)}}{\partial x_i} &&= -\frac{1}{\lambda} f^{(1)} \nonumber \\
\epsilon^1: \quad & \frac{\partial f^{(1)}}{\partial t} + c_i \frac{\partial f^{(1)}}{\partial x_i} &&= -\frac{1}{\lambda} f^{(2)} \tag{1.17}
\end{align}
$$

Ces relations permettent alors de calculer les moments de l'équation (1.15) à des
ordres différents en $\epsilon$. Pour retrouver les équations de conservation décrivant la
dynamique d'un fluide visqueux, il faut ainsi utiliser l'ordre 1 en epsilon, pour
obtenir :

$$
P_{ij}^{(1)} = P_{ij}^{0} + \epsilon \int_{\mathbb{R}^3} (c_i - u_i)(c_j - u_j) f^{(1)} d\mathbf{c} \quad (1.18)
$$

$$
q_i^{(1)} = q_i^0 + \epsilon \frac{1}{2} \int_{\mathbb{R}^3} (c_i - u_i) |\mathbf{c} - \mathbf{u}|^2 f^{(1)} d\mathbf{c} \quad (1.19)
$$

En utilisant la relation (1.17) à l'ordre 0 on obtient alors les expressions de $P_{ij}^{(1)}$ et $q_i^{(1)}$ en fonction des dérivées temporelles et spatiales des moments de la fonction d'équilibre. Les dérivées temporelles pouvant être exprimées en fonction des dérivées spatiales par l'intermédiaire des équations d'Euler, on aboutit finalement en posant $\tau = \epsilon\lambda$ à :

$$
P_{ij}^{(1)} = p\delta_{ij} - \tau\rho r T \left( \frac{\partial u_i}{\partial x_j} + \frac{\partial u_j}{\partial x_i} \right) - \frac{2}{3} \tau \rho r T \frac{\partial u_k}{\partial x_j} \quad (1.20)
$$

$$
q_i^{(1)} = -\frac{D+2}{2} \tau \rho r T \frac{\partial r T}{\partial x_i} \qquad (1.21)
$$

On retrouve alors une viscosité dynamique $\mu = \tau\rho r T$, une viscosité volumique
$\eta = 0$ et une conductivité thermique $\kappa = \frac{D+2}{2}\tau\rho r T$. Ce lien important entre
une équation statistique et des lois de conservation est un des points importants
qui distingue les méthodes LBM des méthodes de résolution numériques directes
des équations de Navier-Stokes. En effet, ce lien montre que lorsque le régime
des milieux continus est assuré, les solutions de l'équation de Boltzmann sont
asymptotiquement liées aux solutions de Navier-Stokes. Ainsi il devient intéressant
d'établir une méthode numérique permettant la résolution numérique de l'équation
(1.15) qui possède des caratéristiques convectives plus simples que les équations
Navier-Stokes.

1.1.2 Equation de Boltzmann à vitesses discrètes

Pour résoudre numériquement (1.15)), il faut d'abord discrétiser l'espace des
vitesses puis discrétiser l'espace et le temps. En effet, l'idée principale des méthodes
LBM est de considérer que les particules du fluide ne peuvent prendre qu'un en-
semble fini de vitesses ayant des directions particulières. Cet ensemble que l'on
appelle le « Réseau »de vitesses impose nécessairement des restrictions physiques
à la validité du modèle ([10, 23]). Afin de construire mathématiquement un réseau consistant, on impose le fait que les moments de la fonction de distribution doivent être conservés, c'est à dire que les équations (1.2), (1.3), (1.4) restent valables sur le réseau de vitesses. Plusieurs approches existent afin de construire systématiquement des réseaux de vitesses. De nombreux auteurs utilisent les quadratures de Gauss-Hermite pour assurer l'égalité des moments jusqu'à l'ordre 2 [30, 55, 56, 40] d'autres préfèrent établir des développements asymptotiques à partir de réseaux supposés connus et établir les équations macroscopiques équivalentes à un ordre donné [20, 19, 21, 49]. Dans les deux cas la fonction de distribution est exprimée à l'aide de développements polynomiaux. Le formalisme de ce manuscrit utilise les polynômes d'Hermite qui sont présentés succinctement ci-dessous.

## 1.1.3 Les polynôme d'Hermite

Les polynômes d'Hermite sont définis par la relation :

$$ \mathcal{H}^{(n)}(x) = (-1)^n e^{x^2/2} \frac{d^n}{dx^n} e^{-x^2/2} \qquad (1.22) $$

Ces polynômes sont utilisés en LBM pour écrire le développement de la fonction de distribution $f$ jusqu'à l'ordre $N$ [28, 27, 55, 56] :

$$ f^N = \omega(\mathbf{c}) \sum_{n=0}^{N} \frac{1}{n!} \mathcal{H}^{(n)}(\mathbf{c}) : \mathbf{a}^{(n)} \qquad (1.23) $$

ou $\mathbf{a}$ est un coefficient dépendant de la distribution $f$ :

$$ \mathbf{a}^{(n)} = \int f \mathcal{H}^{(n)} d\mathbf{c} \qquad (1.24) $$

Le calcul des premiers polynômes montre donc que les coefficients (1.24) sont reliés aux moments (1.2) et (1.4). Ainsi, le développement en polynômes d'Hermite des fonctions de distributions appliqué au développement (1.16) permet d'établir une contrainte d'égalité à l'ordre $N$ des moments continus et des moments discrets des fonctions de distributions.

## 1.1.4 Le modèle à 19 vitesses : D3Q19

Cette contrainte est alors satisfaite sur un réseau discret de vitesse en utilisant les quadratures de Gauss-Hermite [56]. Les réseaux les plus couramment utilisés sont donc les solutions de ces quadratures et donnent des vitesses basées sur un stencil simple utilisant les plus proches voisins soit 9 vitesses en 2D et 27 vitesses en 3D avec la vitesse nulle. Ces réseaux standards, par le faible nombre de vitesses discrètes utilisées, n'imposent l'égalité des moments que jusqu'à l'ordre 2. Ainsi, l'équation sur l'énergie n'est pas correctement conservée restreignant l'utilisation de ces réseaux à des applications isothermes. La fonction d'équilibre est développée à l'ordre 2 :
$$f_{\alpha}^{eq}(\mathbf{x}, t) = \rho\omega_{\alpha}\left(1 + \frac{\mathbf{u}\cdot\mathbf{c}_{\alpha}}{(rT)} + \frac{(\mathbf{u}\cdot\mathbf{c}_{\alpha})^2}{2(rT)^2} - \frac{|\mathbf{u}|^2}{2(rT)}\right) \quad (1.25)$$

et les contraintes du réseau imposent une température réduite $rT = 1/3$. De plus, le développement de Chapmann-Enskog de ces réseaux montre que le tenseur des contraintes visqueuses obtenu par (1.18) devient :

$$\tau_{ij} = \mu \left( \frac{\partial u_i}{\partial x_j} + \frac{\partial u_j}{\partial x_i} \right) - \tau \frac{\partial \rho u_i u_j u_k}{\partial x_k} \quad (1.26)$$

Cette erreur en $O(M^3)$ limite donc l'application de ce type de modèle aux écoulements à faible nombre de Mach. Des approches proposent de réduire notablement cette erreur et seront présentées dans les prochaines sections.

En pratique, pour limiter le temps de calcul et l'occupation mémoire, le modèle à 27 vitesses est souvent dégradé à 19 vitesses sans restriction supplémentaire. Les pondérations et vitesses du modèle sont données dans le tableau 1.2.

<figure><img src="image_5.png" /><figcaption>FIGURE 1.1 – Vitesses discrètes du modèle D3Q19</figcaption></figure>

<table><caption>TABLE 1.2 – Pondération de la fonction d'équilibre et vitesses discrètes du modèle D3Q19.</caption><thead><tr><th>α</th><th>1</th><th>2-7</th><th>8-19</th></tr></thead><tbody><tr><td rowspan="3">c̃<sub>α</sub></td><td rowspan="3">(0, 0, 0)</td><td>(0,0,±1)</td><td>(0,±1,±1)</td></tr><tr><td>(0,±1,0)</td><td>(±1,±1,0)</td></tr><tr><td>(±1,0,0)</td><td>(±1,0,±1)</td></tr><tr><td>ω<sub>α</sub></td><td>1/3</td><td>1/18</td><td>1/36</td></tr></tbody></table>
Le réseau D3Q19, sauf indication contraire, sera le modèle utilisé dans tout le mémoire. Il permet donc une résolution numérique de l'équation suivante :

$$ \frac{\partial f_{\alpha}}{\partial t} + c_{\alpha,i} \frac{\partial f_{\alpha}}{\partial x_i} = -\frac{1}{\tau} [f_{\alpha} - f_{\alpha}^{eq}] \quad (1.27) $$

sous les conditions d'écoulement isotherme à faible nombre de Mach.

## 1.2 Les différents schémas LBM

Une fois le réseau de vitesses et ces contraintes assumées, il est alors possible de construire un grand nombre de schémas différents possédant chacun des caractéristiques différentes en terme de stabilité et dissipation. Les différentes approches possibles sont présentées dans les sections suivantes et sont résumées dans le diagramme de la figure 1.2.

### 1.2.1 Le schéma LBM standard

S'il est possible de discrétiser (1.27) à l'aide de méthodes classiques comme les différences finies ou les volumes finis [57, 63], la désignation LBM fait cependant référence à une discrétisation particulière basée sur l'intégration le long de la caractéristique $c_\alpha$. En effet, en intégrant l'équation (1.27) le long de cette caractéristique, on peut alors écrire :

$$ f_{\alpha}(\mathbf{x} + c_{\alpha}\Delta t, t + \Delta t) - f_{\alpha}(\mathbf{x}, t) = -\frac{1}{\tau} \int_{0}^{\Delta t} [f_{\alpha}(\mathbf{x} + c_{\alpha}s, t + s) - f_{\alpha}^{eq}(\mathbf{x} + c_{\alpha}s, t + s)] ds \quad (1.28) $$

En intégrant le membre de droite avec la méthode des trapèzes, on obtient l'équation implicite :

$$ f_{\alpha}(\mathbf{x} + c_{\alpha}\Delta t, t + \Delta t) - f_{\alpha}(\mathbf{x}, t) = -\frac{\Delta t}{2\tau} [f_{\alpha}(\mathbf{x} + c_{\alpha}\Delta t, t + \Delta t) - f_{\alpha}^{eq}(\mathbf{x} + c_{\alpha}\Delta t, t + \Delta t) + f_{\alpha}(\mathbf{x}, t) - f_{\alpha}^{eq}(\mathbf{x}, t)] + \mathcal{O}(\Delta t)^3 \quad (1.29) $$

L'obtention de l'équation de Boltzmann sur Réseau explicite se fait alors grâce au changement de variable [18] :

$$ g_{\alpha}(\mathbf{x}, t) = f_{\alpha}(\mathbf{x}, t) + \frac{\Delta t}{2\tau} (f_{\alpha}(\mathbf{x}, t) - f_{\alpha}^{eq}(\mathbf{x}, t)) \quad (1.30) $$

qui permet alors d'écrire :

$$ g_{\alpha}(\mathbf{x} + c_{\alpha}\Delta t, t + \Delta t) = g_{\alpha}(\mathbf{x}, t) - \frac{\Delta t}{\tau_g} [g_{\alpha}(\mathbf{x}, t) - g_{\alpha}^{eq}(\mathbf{x}, t)] + \mathcal{O}(\Delta x^2, \Delta t^2) \quad (1.31) $$
avec $\tau_g = \tau + \frac{1}{2}$ et $g_{\alpha}^{eq} = f_{\alpha}^{eq}$.

Ainsi, la LBM possède une précision d'ordre 2 en espace et en temps. D'après l'équation (1.31), le réseau obtenu est donc constitué de mailles cubiques de côté $\Delta x$. Ainsi, les fonctions de distribution se déplaceront d'une maille $\Delta x$ en un temps donné $\Delta t$ soit : $\Delta x = |c_{\alpha}|\Delta t$. Cette relation particulière entre le pas de temps et le pas d'espace est souvent appelé *Acoustic Scaling* en référence au coefficient de proportionnalité entre les deux paramètres. D'autres auteurs [25] proposent d'utiliser une relation quadratique ou *Diffusive scaling* entre la pas de temps et le pas d'espace du type $\Delta t = \alpha \Delta x^2$ pour obtenir des modèles purement incompressibles. Ces approches ne seront pas discutées ici.

Pour des raisons pratiques, on considère souvent $\Delta x = \Delta t = 1$ dans l'algorithme (1.31) et les grandeurs calculées sont alors sans dimension (ou souvent appelées *Lattice Unit*). Les valeurs physiques peuvent alors être retrouvées en utilisant celles de $\Delta t$ et $\Delta x$. Ici, on notera les grandeurs sans dimension avec un chapeau. Par exemple, la vitesse du son adimensionnelle $\hat{c}_0$ sera définie par :

$$ \hat{c}_0 = \sqrt{r\hat{T}} = \frac{1}{\sqrt{3}} \qquad (1.32) $$

Immédiatement, la vitesse du son réelle, sera alors définie par :

$$ c_0 = \hat{c}_0 \frac{\Delta x}{\Delta t} \qquad (1.33) $$

Enfin, la viscosité du fluide simulé est reliée au temps de relaxation $\tau_g$ par :

$$ \hat{\tau}_g = \frac{\nu\sqrt{3}}{\Delta x c_0} + \frac{1}{2} \qquad (1.34) $$

Ainsi l'équation (1.31) est la version d'origine du schéma souvent désigné sous la forme LBM-BGK. Cette version souffre d'un grand nombre de limitations et n'est quasiment plus utilisée aujourd'hui. En particulier, le schéma devient instable lorsque la viscosité devient faible (Haut Reynolds). Pour palier à ces différents problèmes, plusieurs améliorations sont aujourd'hui disponibles et sont principalement classées en deux catégories : les approches proposant une modification de l'équilibre et les approches proposant une modification de la relaxation.

## 1.2.2 Les approches basées sur une modification de l'équilibre

### 1.2.2.1 Les approches régularisées

L'algorithme classique 1.31 est basé sur une écriture de l'équilibre développé à l'ordre deux à l'aide des polynômes d'Hermite 1.22. En écrivant les fonctions de distribution sous la forme $g_{\alpha} = g_{\alpha}^{eq} + g_{\alpha}^{neq}$ l'équation (1.31) devient :

$$ g_{\alpha}(\mathbf{x} + c_{\alpha}\Delta t, t + \Delta t) = g_{\alpha}^{eq}(\mathbf{x}, t) + \left(-\frac{\Delta t}{\tau_g}\right) g_{\alpha}^{neq}(\mathbf{x}, t) \quad (1.35) $$
Cette équation est en tout point équivalente à (1.31) à ceci près qu'elle utilise une formulation basée sur la partie hors équilibre $g^{neq}$. Il vient alors à l'esprit que cette quantité peut être recalculée en utilisant les polynômes d'Hermite (1.22) :

— Latt [37] propose de recalculer la partie hors équilibre après chaque collision en l'assimilant à $f^{(1)}$ dans le développement de Chapmann-Enskog. Elle est alors évaluée en utilisant l'estimation issue du développement :

$$f_{\alpha}^{(1)} = -\frac{\omega_{\alpha} \rho \tau_g \Delta t}{c_s^2} (c_{\alpha,i} c_{\alpha,j} - c_s^2 \delta_{ij}) S_{ij} \quad (1.36)$$

Ou $S_{ij}$, le tenseur des déformations est évalué par des différences finies centrées d'ordre 2.

— Malaspinas [40] propose une estimation de la partie hors équilibre basée sur une formulation récurente du développement en polynômes d'Hermite.

— Jacob [32] propose de coupler les deux approches précédentes avec un coefficient $\sigma$ permettant une évaluation hybride de $g^{neq}$.

Ces approches s'avèrent être bien plus stables que les modèles classiques mais peuvent cependant devenir très dissipatives en fonction de la valeur de $\sigma$.

### 1.2.2.2 Les approches entropiques

Les approches entropiques [54, 60] sont basées sur une relaxation locale en introduisant un paramètre de « sur-relaxation » $\beta$ calculé en résolvant un problème de minimisation basé sur le théorème H.

$$g_{\alpha}(\mathbf{x} + c_{\alpha} \Delta t, t + \Delta t) = g_{\alpha}(\mathbf{x}, t) - \frac{\beta \Delta t}{\tau_g} (g_{\alpha}(\mathbf{x}, t) - g_{\alpha}^{\text{neq}}(\mathbf{x}, t))(1.37)$$

Cette approche revient à utiliser un temps de relaxation local et à augmenter la dissipation lorsque les fonctions de distributions s'approchent de la valeur nulle. Les auteurs du modèle ont également montré le lien du modèle entropique avec les modèles de sous-maille de type Smagorinsky [3]. Une version utilisant plusieurs temps de relaxation a également été développée afin de palier aux lourdeurs engendrées par la résolution du problème de minimisation [9].

Les approches développées dans la section 1.4 sont également très similaires au modèle entropique en utilisant cependant des critères adaptatifs plus faciles et plus rapides à estimer.

## 1.2.3 Les approches basées sur une amélioration de la relaxation

Les modèles présentés précédemment sont basés sur une relaxation simple en supposant que chaque moment de la fonction de distribution est relaxé à l'équilibre avec le même temps caractéristique. On peut cependant considérer que les moments peuvent être relaxés avec des temps propres. Cette approche permet d'associer des temps de relaxation différents pour chaque moment physique. Cependant, pour que l'étape de collision soit consistante avec le réseau de vitesses, il faut donc choisir un nombre de moments identiques. Ainsi, en trois dimensions, on peut définir
facilement les cinq premiers moments (1.2) (1.3) et (1.4) mais il faut construire quatorze nouveaux moments dans le cadre d'un réseau à 19 vitesses. Le choix de ces nouveaux moments est particulièrement disputé et donne naissance à différentes formulations.

### 1.2.3.1 Le modèle MRT

La construction originale de d'Humière fut obtenue en construisant les moments à partir de combinaisons des moments physiques puis en orthogonalisant par une procédure de Gram-Schmidt. Une autre approche consiste à construire les moments à partir des polynômes d'Hermite via la relation 1.24. Dans tous les cas il faut également définir un paramètre de relaxation pour chaque moment afin de construire la matrice $S = diag[s_\alpha, \alpha = 1..19]$ avec $s_\alpha = 1/\tau_\alpha$.

Les coefficients du modèle de d'Humière ont été calculés pour une stabilité optimale [15, 36]. Ainsi, la matrice S proposée est la suivante :

$$ S = diag(0, s_1, s_2, 0, s_3, 0, s_3, 0, s_3, s_\nu, \quad (1.38) \\ s_2, s_\nu, s_2, s_\nu, s_\nu, s_\nu, s_4, s_4, s_4) $$

Une fois les moments **m** déterminés on peut procéder à l'étape de collision dans l'espace des moments en utilisant une matrice de passage **M** pour relier les fonctions de distribution aux moments **m** = **Mg**.

La collision du modèle s'écrit dans l'espace des moments :

$$ g(x + c, t + 1) = g(x, t) - M^{-1}S[m(x, t) - m^{eq}(x, t)] \quad (1.39) $$

Les coefficients de transport s'écrivent :

$$ \left\{ \begin{aligned} \hat{\nu} &= \frac{1}{3} \left( \frac{1}{s_{\nu}} - \frac{1}{2} \right) \\ \hat{\eta} &= \hat{\rho} \frac{2}{9} \left( \hat{\tau}_{1} - \frac{1}{2} \right) \end{aligned} \right. \quad (1.40) $$

Ainsi, la viscosité en volume $\eta$ peut être réglée indépendamment de $\nu$ et est arbitrairement élevée dans le modèle d'origine afin de garantir une meilleure stabilité. En revanche cela impose une mauvaise propagation des ondes acoustiques qui sont donc artificiellement dissipées par le modèle MRT.

### 1.2.3.2 Le modèle TRT

Le modèle TRT [33] (Two Relaxation Times) est une variante du modèle MRT utilisant seulement deux temps de relaxation. Le premier temps de relaxation est relié à la viscosité et concerne les moments d'ordre pair, le second temps de relaxation est un paramètre qui peut être choisi librement et concerne les moments d'ordre impair. En pratique il existe un choix judicieux du second temps de relaxation basé sur l'expression :
$$ \Lambda = \left( \frac{\tau_1}{\Delta t} - \frac{1}{2} \right) \left( \frac{\tau_2}{\Delta t} - \frac{1}{2} \right) \qquad (1.41) $$

qui montre que des valeurs particulières de $\Lambda$ (appelé *Magic parameter*) peuvent augmenter l'ordre de précision et la stabilité du schéma. Ainsi, en choisissant $\Lambda = \frac{1}{12}$, on peut annuler le terme d'ordre 3 intervenant dans les développements [16]. D'autres études montrent que $\Lambda = \frac{1}{4}$ donne des simulations plus stables [26, 35].

### 1.2.3.3 Les moments centraux et ses dérivées

Le principe des modèles à moments centraux est d'utiliser une définition des moments tenant compte de la composante principale de la vitesse. Cela revient à écrire les moments dans le référentiel mobile de l'écoulement. Ce décalage permet de régler certains problèmes d'invariance Galiléenne et permet une meilleure stabilité. Il fut introduit initialement par Geier [24] dans un formalisme assez différent des approches LBM puis fut ensuite étendu et repris par De Rosis [53, 52] puis par Février [22]. La relaxation des moments fut par la suite généralisée à celle des cumulants [25] dont les différences surviennent principalement à partir des moments d'ordre 4.

#### LBM Scheme Architecture

<figure><img src="image_6.png" /><figcaption>FIGURE 1.2 – Architecture d'un schéma LBM dans le diagramme (Equilibre, Collision, Post-Collision). Les flèches vertes correspondent au schéma standard (1.31), les flèches rouges aux schémas MRT et TRT (1.39) et les flèches bleues aux approches régularisées (1.35).</figcaption></figure>
## 1.3 Stabilisation par filtrage

Les travaux effectués lors de ma thèse [50] ont montrés que les instabilités inhérentes aux schémas LBM, en particulier dans sa version BGK, pouvaient être supprimées en filtrant l'opérateur de collision à l'aide des filtres sélectifs.

### 1.3.1 Filtrage sélectif simple

Un filtre sélectif consiste à introduire une viscosité artificielle sur une gamme de fréquence donnée. Le principe de ce genre de technique fut introduit initialement par von Neumann et Richtmyer [47]. Aujourd'hui les filtres sélectifs sont largement utilisés en simulations numériques pour la stabilisation de modes parasites. On peut écrire le filtrage $\langle \cdot \rangle$ d'une quantité $v$ sous la forme :

$$ \langle v(x) \rangle = v(x) - \sigma \sum_{j=1}^{D} \sum_{n=-N}^{N} d_n v(x + n \Delta x_j) \quad (1.42) $$

En pratique l'amplitude du filtre est souvent fixée à une valeur comprise entre $\sigma_0 = 0.0$ et $\sigma_0 = 0.1$. On peut alors augmenter la sélectivité en fréquence en augmentant l'ordre du filtre ce qui se traduit par l'utilisation d'un plus grand nombre de points. En pratique les simulations réalisées sur une double couche de mélange périodique ont montrés que le filtrage de la partie hors-équilibre $g_{\alpha}^{heq} = g_{\alpha} - g_{\alpha}^{eq}$ avec des filtres à 5 points donnait de bons résultats. On peut alors remarquer que cette méthode de filtrage peut être assimilée à l'utilisation d'un temps de relaxation modifié $\tau_g^f = \frac{\tau_g}{1 - \sigma F}$, où $F$ représente la fonction de transfert du filtre.

Ce type d'approche est donc relativement efficace mais augmente de façon arbitraire le stencil du schéma puisque le recours à des filtres à 5 points requiert l'utilisation des voisins d'ordre 2.

### 1.3.2 Filtrages adaptatifs

Une approche introduite par la suite pour améliorer l'efficacité des filtres fut de modifier localement la valeur du coefficient $\sigma$ [42]. En effet, les instabilités sont la plupart du temps liées à la présence de forts cisaillements, ainsi l'utilisation du filtrage dans des zones de vitesse uniforme devient superflu et peut dégrader localement la dissipation du schéma. Ainsi, on peut introduire une sélectivité supplémentaire basée sur le calcul du cisaillement local :

$$ \xi(\mathbf{x}) = \left(1 - e^{-\left(\frac{|S(\mathbf{x})|}{S_{max}}\right)^2}\right)^2 \quad (1.43) $$

Avec :

$$ S_{max} = \eta \frac{\sqrt{2U_0^2}}{2\nu} \sim \frac{Re_{\delta}U_0}{\delta} \quad (1.44) $$

Le critère proposé (1.44) est basé sur une maximisation du cisaillement estimé à partir des critères développés dans [54, 60] imposant la positivité des fonctions de distribution. On peut alors relier ce filtrage adaptatif aux approches entropiques
1.37 en utilisant $\beta = 1 - F(\xi)$ où $F(\xi)$ représente la fonction de transfert du filtre utilisé avec un coefficient $\sigma = \xi\sigma_0$.

Ces approches ont été testées et validées sur le tourbillon de Taylor-Green [42] et ont montré que la sélectivité spatiale permettait de diminuer l'ordre du filtre tout en gardant l'efficacité du filtrage. Ainsi, le fait de filtrer localement les zones cisaillées permet de gagner du temps de calcul en utilisant des filtres d'ordre plus faible.

<figure><img src="image_7.png" /><figcaption>FIGURE 1.3 – Comparaison des filtres simples (lignes pleines) et des filtres adaptatifs (pointillés) pour la simulation du Taylor-Green vortex sur deux grilles différentes : 96³ (à gauche) et 128³ (à droite). Les couleurs Rouge, bleu et verte représentent les filtres à 3-points, 5-points et 9-points respectivement.</figcaption></figure>

## 1.4 Temps de relaxation adaptatifs

On peut alors appliquer le senseur (1.43) à l'utilisation des modèles à temps de relaxation multiples. Il suffit de modifier la matrice des temps de relaxation $S$ de l'équation (1.38) en imposant :

$$s'_{i} = s_{i} + (1 - \xi)(s_{\nu} - s_{i}) \qquad (1.45)$$

L'interêt de cette approche est de pouvoir profiter de la stabilité des modèles MRT tout en gardant la faible dissipation des ondes acoustiques apportée par la relaxation du BGK. On a montré [41] que cette approche permettait de propager efficacement les ondes acoustiques issues de l'écoulement autour d'un cylindre même dans le cas de faibles résolutions.
<figure><img src="image_8.png" /><figcaption>FIGURE 1.4 – Pression acoustique obtenue avec les temps de relaxation adaptatifs (à gauche) et les temps de relaxation classiques (à droite). Niveaux symétriques de $-1.10^{-3}$ à $1.10^{-3}$</figcaption></figure>
# Chapitre 2

## Comparaison et couplage avec d'autres méthodes

Dans ce chapitre, je vais présenter mes activités liées à l'analyse et à la comparaison des schémas LBM avec les autres approches utilisées en CFD. En effet, la LBM présente des particularités qui nécessitent d'être rappelées dans le contexte des autres méthodes. Ainsi après avoir présenté le principe de l'analyse, nous verrons comment il devient possible d'envisager la construction de solveurs qui tireraient partie des avantages de chaque méthode pour la résolution efficace des problèmes de CFD.

### 2.1 Analyse de stabilité des schémas Boltzmann discrets

L'analyse de stabilité linéaire (LSA) consiste à linéariser un opérateur autour d'un état de référence souvent uniforme et homogène et d'analyser les valeurs propres de cet opérateur. Dans le cas de la LBM, on linéarise l'équation (1.31) en décomposant les fonctions de distribution en une partie uniforme $f_\alpha^{(0)}$ et une partie fluctuante $f'_\alpha$. Les termes non-linéaires de cette équation sont contenus dans la fonction d'équilibre qu'il faut alors linéariser :

$$f^{\text{eq}}(f_{\alpha}^{(0)} + f'_{\alpha}) = f_{\alpha}^{\text{eq},(0)} + \mathcal{D}_{\alpha\beta} f'_{\alpha} + o(f'_{\alpha}^2) \quad (2.1)$$

avec

$$\mathcal{D}_{\alpha\beta} = \left. \frac{\partial f_{\alpha}^{\text{eq}}}{\partial f_{\beta}} \right|_{f_{\beta}=f_{\beta}^{(0)}} \quad (2.2)$$

La forme exacte de $\mathcal{D}_{\alpha\beta}$ peut être évaluée analytiquement. En appliquant ces développements dans le cadre des schémas LBM on a une équation linéarisée du type :

$$e^{-i\omega f'} = M^{\text{coll}} f' \quad (2.3)$$
où $M^{\text{coll}}$ est une matrice dépendante du modèle de collision. On peut alors appliquer ce type d'analyse à l'ensemble des modèles de collision et résoudre le problème aux valeurs propres pour en tirer les propriétés dispersives et dissipatives des modèles de collision. Des études récentes ont montrée une comparaison intéressante des modèles de collisions dans le contexte de leur analyse de stabilité [61, 14]. La figure 2.1 présente les propriétés des principaux modèles de collision.

<figure><img src="image_9.png" /><figcaption>FIGURE 2.1 – Comparaison des analyses de stabilité linéaire des trois modèles de collision principaux pour le mode Acoustique principal. $M_a = 0.2$ et $\tau = 0.0025$.</figcaption></figure>

Au regard des résultats (figure 2.1), on remarque d'abord que pour les paramètres convectifs choisis, les 3 modèles présentent des caractéristiques dispersives similaires. Le modèle MRT est nettement plus dissipatif pour les ondes acoustiques car il est nativement construit avec une viscosité en volume plus importante. Le modèle HRR est quant à lui similaire au modèle BGK en terme de dissipation des ondes acoustiques. Une autre caractéristique des schémas LBM observée sur la figure 2.1, est la propagation des oscillations maille à maille ($\bar{k} = \pi$). En effet, pour les schémas classiques, les informations à 2 points par longueurs d'ondes ne sont pas propagées car l'erreur de dispersion pour ces ondes est de 100%. Pour la LBM, ces ondes ont une erreur de dispersion mais leur vitesse de propagation reste non nulle.

Ainsi, après avoir caractérisé la dispersion et la dissipation des modèles LBM, il serait intéressant de pouvoir comparer ces résultats à ceux obtenus avec les schémas numériques utilisés en aéroacoustique.

## 2.2 Comparaison avec les méthodes Navier-Stokes d'ordre élevé

Les analyses réalisées précédemment peuvent être appliquées aux schémas aux différences finies centrées d'ordre $N$ utilisant $2N + 1$ points. Pour cela on évalue les dérivées spatiales avec la discrétisation suivante :
$$ \frac{\partial \mathbf{U}}{\partial x_i}(x_i^0) = D_i(x_i^0) = \frac{1}{\Delta x_i} \sum_{j=-N}^{N} a_j \mathbf{U}(x_i^0 + j \Delta x_i) \quad (2.4) $$

et l'intégration temporelle à l'aide d'un schéma de Runge-Kutta à p-pas :

$$ \mathbf{U}^{n+1} = \mathbf{U}^n + \sum_{j=1}^{p} \gamma_j \Delta t^j F^j(\mathbf{U}^n) \quad (2.5) $$

Les coefficients $a_j$ dépendent de l'ordre du schéma. Pour l'aéroacoustique, ces schémas sont généralement d'ordre élevé (ordre 6 ou 7) pour satisfaire aux exigences d'un calcul acoustique. Pour les schémas standards, ils sont calculés par identification aux développements de Taylor. Il existe également des schémas dit optimisés [59, 7] qui calculent les coefficients $a_j$ de façon à minimiser l'erreur de dispersion commise. Les coefficients $\gamma_j$ ont ici le même rôle que les coefficients $a_j$ des différences finies. Ils sont calculés par identification aux développements de Taylor pour les schémas standards et par minimisation des erreurs de dispersion et de dissipation pour les schémas optimisés [5, 7].

<figure><img src="image_10.png" /><figcaption>FIGURE 2.2 – Comparaison des erreurs de dissipation des schémas LBM (ligne pleine) et Navier-Stokes (symboles) pour le mode acoustique principal convecté à $M_a = 0.2$. [43].</figcaption></figure>

La figure 2.2 montre les erreurs de dissipation sur le mode acoustique principal pour la LBM et pour quelques schémas différences finies optimisés. Ces comparaisons montrent que la méthode Boltzmann sur Réseau présente des capacités acoustiques bien supérieures aux schémas différences finies classiques, en particulier pour la propagation acoustique (Les détails de cette étude sont disponibles dans [43]). Il devient dès lors interessant de pouvoir coupler ces caractéristiques de la LBM aux performances proches-paroi des solveurs Navier-Stokes.
### 2.3 Couplage des méthodes LBM avec les méthodes Navier/Stokes

Comme nous l'avons vu précédemment, les maillages uniformes imposés par la construction des schémas LBM rendent les simulations de couches limites particulièrement inadaptées en LBM. En revanche ces propriétés acoustiques peu dissipatives associées à sa faible compressibilité en font un outil de premier choix pour les phénomènes aéroacoustiques. En revanche, les approches Navier-Stokes classiques sont particulièrement performantes dans le voisinage de la paroi de part l'utilisation de maillage à très grand rapport d'aspect. En effet, les maillages différences finies ou volumes finis permettent de garder une bonne précision à la paroi en utilisant des mailles longitudinales très étendues.

Dès lors, il devient intéressant d'étudier le couplage entre ces deux types d'approches qui présentent chacune des intérêts propres dans des zones bien distinctes. Ce travail est l'objet de la thèse d'Alexandre SUSS lancé en 2020 en collaboration avec l'ONERA. La première étape a été de mettre en place une bonne stratégie numérique pour permettre une transition consistante entre deux approches qui utilisent chacune des grandeurs différentes en particulier pour la résolution temporelle. De plus la transition entre deux solveurs engendre un grand nombre de problèmes inhérents à la simulation numérique :

- Tout d'abord, le problème du couplage peut être vu comme un problème de conditions aux limites et pose dès lors la question des quantités conservées à la frontière. Certains auteurs [38] ont étudié ce point dans une approche purement incompressible en focalisant leurs efforts sur les quantités à l'équilibre et obtiennent quelques résultats dans le cas d'un écoulement simple de type Poiseuille.

- D'autre part, si le problème du passage LB-NS semble bien posé en utilisant la définition des moments, le passage inverse posent des problèmes d'interpolation. En effet, le solveur NS ne donne que des informations macroscopiques (densité, vitesse, température) alors que le solveur LBM a besoin de quantités mésoscopiques (19 variables en 3D). Certains auteurs proposent de passer traditionnellement par la fonction d'équilibre [31] dans le cadre d'une équation d'advection diffusion. Cependant cette approche est mise en défaut en présence de termes convectifs et en particulier lorsqu'ils présentent de forts gradients de vitesse. De plus, la nature faiblement compressible de la LBM peut donner naissance à des ondes de réflexion qui peuvent se propager et perturber les solutions en aval de la zone LBM.

- Une autre difficulté liée au couplage des méthodes LBM et Navier-Stokes se trouve dans l'intégration temporelle. En effet les méthodes LBM, purement instationnaires, sont basées sur une formulation couplée en espace-temps à CFL constant. Dès lors il faut s'interroger sur le choix des méthodes d'intégration Navier-Stokes pertinentes pour le couplage [48].

- Ensuite, un point essentiel réside dans le choix du modèle de collision a utiliser pour le couplage. Les études existantes sont majoritairement basées
sur l'opérateur BGK connu également pour ses problèmes de stabilité nu-
mérique à grand nombre de Reynolds. Dès lors, l'utilisation des modèles de
collision alternatifs (MRT, Cumulant, Entropic, Régularisation) peut être
envisagée.

— Enfin, le raccordement de maillages de type couche limite à des maillages de type LBM pose également des problèmes de compatibilité de données. Une solution consiste à utiliser des recouvrements permettant une transition plus douce entre les deux méthodes de résolutions [48]. Une autre possibilité réside dans l'utilisation de maillages consistant avec les deux approches. Des maillages de type Octree ont par exemple été testés avec succès dans des cas simples [44]. La compatibilité des ces approches avec les outils existants reste a définir.

Des premiers résultats encourageant ont été obtenus en utilisant un solveur Navier-Stokes basé sur des schémas spatiaux de type senseurs [1] d'ordre 2 et des schémas temporelles de type Runge-Kutta d'ordre 3. La principale difficulté réside dans la reconstruction précise des fonctions de distributions à partir des informations macroscopiques issues de la partie Navier-Stokes. Pour cela les pistes privilégiées sont basées sur une reconstruction régularisée [32] de la partie hors-équilibre aux interfaces.

<figure><img src="image_11.png" /><figcaption>FIGURE 2.3 – Propagation d'un pulse de pression d'un domaine LBM au centre vers 7 domaines Navier-Stokes, profil de masse volumique obtenu après passage des interfaces.</figcaption></figure>

## 2.4 Comparaison avec les méthodes Vortex

Les méthodes vortex [2, 4] et les méthodes LBM incompressibles sont des méthodes qui s'écartent toutes les deux des méthodes standards basées sur Navier-Stokes. D'abord, elles opèrent toutes les deux dans un espace qui n'est pas complètement macroscopique, et peuvent partager dans certains cas des stencils identiques. En revanche, la construction hétérogène des méthodes vortex rend leur étude globale particulièrement difficile en particulier dans le cadre des approches
linéarisées. Dès lors, il devient intéressant de pouvoir comparer les caractéristiques relatives de ces deux méthodes afin de mieux comprendre leur comportement et de mieux les situer dans le paysage numérique. Dans le cadre de cette collaboration avec le laboratoire M2N [17, 46], les méthodes Vortex avec noyau de remaillage [39] ont été étudiées et comparées à une méthode LBM incompressible standard de type MRT sur des cas tests 3D bien documentés.

<figure><img src="image_12.png" /><figcaption>FIGURE 2.4 – Evolution temporelle de l'enstrophie du Taylor-Green Vortex. A gauche LBM, a droite VM.</figcaption></figure>

Les résultats montrent que la partie lagrangienne des méthodes vortex leur confère des propriétés dispersives plus intéressantes en particulier pour les faibles résolutions. De plus pour les faibles résolutions il apparaît que les méthodes VM présentent une sous-dissipation qui rend leur convergence atypique par rapport aux comportements observés en LBM et en Navier-Stokes (Cf Fig.2.4). Enfin, les simulations sur un cube 3D montrent des résultats similaires malgré l'importante sensibilité du cas test aux perturbations apportées par les différences de traitement des conditions de paroi dans les deux méthodes [45]. Ces comparaisons ont montré que les deux méthodes pouvaient être complémentaires et de futures investigations comme la construction d'une hybridation des deux schémas pourrait être envisagée [13].
Références

[1] Nicolas ALFEREZ. "Simulation des grandes échelles du processus de décrochage par éclatement de bulbe de décollement laminaire". Theses. ISAE-ENSMA Ecole Nationale Supérieure de Mécanique et d'Aérotechnique - Poids-tiers, mar. 2014 (cf. p. 33).

[2] C. ANDERSON et C. GREENGARD. "On vortex methods". In: *SIAM J. Num. Anal.* 22 (1985), p. 413-440 (cf. p. 33).

[3] S. ANSUMALI, I.V. KARLIN et S. SUCCI. "Kinetic theory of turbulence modeling: smallness parameter, scaling and microscopic derivation of Smagorinsky model." In: *Physica A* 338 (2004), p. 379-394 (cf. p. 22).

[4] J.T. BEALE et A. MAJDA. "Vortex methods II : high order accuracy in 2 and 3 dimensions". In : *Math. Comput.* 32 (1982), p. 29-52 (cf. p. 33).

[5] J. BERLAND, C. BOGEY et C. BAILLY. "Low-dissipation and low-dispersion fourth order Runge-Kutta algorithm". In : *Computers and Fluids* 10 (2006), p. 35 (cf. p. 31).

[6] P.L. BHATNAGAR, E.P. GROSS et M. KROOK. "A model for collision process in gases. I. Small amplitude process in charged and neutral one-component systems". In : *Phys. Rev.* 94(3) (1954), p. 511-525 (cf. p. 15).

[7] BOGEY, C. AND BAILLY, C. "A family of low dissipative explicit schemes for flow and noise computations". In : *J. Comput. Phys.* 194 (2004), p. 194-214 (cf. p. 31).

[8] BOLTZMANN, LUDWIG. "Weiter Studien über das Wärmegleichgewicht unter Gas-molekülen." In : *Wien. Ber.* 66 (1878), p. 275-370 (cf. p. 14).

[9] Fabian BÖSCH, Shyam S. CHIKATAMARLA et Ilya V. KARLIN. "Entropic multirelaxation lattice Boltzmann models for turbulent flows". In : *Phys. Rev. E* 92 (4 oct. 2015), p. 043309. DOI : 10.1103/PhysRevE.92.043309 (cf. p. 22).

[10] J.E. BROADWELL. "Study of rarified shear flow by the discrete velocity method." In : *Journal of Fluid Mechanics* 19 (1964), p. 401-414 (cf. p. 18).

[11] CERCIGNANI, C. *Theory and application of the Boltzmann Equation*. Elsevier., 1975 (cf. p. 14).

[12] CHAPMAN, S. AND COWLING, T.G. *The Mathematical theory of Non-Uniform gases*. Cambridge, 3rd edition., 1991 (cf. p. 16).
[13] Leitao CHEN et al. "Semi-Lagrangian implicit Bhatnagar-Gross-Krook collision model for the finite-volume discrete Boltzmann method". In: *Physical Review E* 101.6 (juin 2020). DOI: 10.1103/physreve.101.063301 (cf. p. 34).

[14] C. COREIXAS et al. "Impact of collision models on the physical properties and the stability of lattice Boltzmann methods". In: *Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences* 378.2175 (juin 2020), p. 20190397. ISSN: 1471-2962. DOI: 10.1098/rsta.2019.0397 (cf. p. 30).

[15] D. D'HUMIÈRE et al. "Multiple relaxation time lattice Boltzmann models in three dimensions". In: *Phil. Trans. R. Soc. Lon. A* 360 (2002), p. 437-451 (cf. p. 23).

[16] Dominique D'HUMIÈRES et Irina GINZBURG. "Viscosity independent numerical errors for Lattice Boltzmann models : From recurrence equations to "magic" collision numbers". In : *Computers & Mathematics with Applications* 58.5 (sept. 2009), p. 823-840. DOI: 10.1016/j.camwa.2009.02.008 (cf. p. 24).

[17] François DELASSAUX et al. "Sensitivity analysis of hybrid methods for the flow around the ahmed body with application to passive control with rounded edges". In : *Computers Fluids* 214 (jan. 2021), p. 104757. DOI: 10.1016/j.compfluid.2020.104757 (cf. p. 34).

[18] DELLAR, P.J. "Bulk and shear viscosities in lattice Boltzmann equations". In : *Phys. Rev. E* 64 (2003) (cf. p. 20).

[19] François DUBOIS. "Equivalent partial differential equations of a lattice Boltzmann scheme". In : *Computers & Mathematics with Applications* 55.7 (avr. 2008), p. 1441-1449. DOI: 10.1016/j.camwa.2007.08.003 (cf. p. 18).

[20] François DUBOIS. "Third order equivalent equation of lattice Boltzmann scheme". In : *Discrete and Continuous Dynamical Systems* 23.1/2 (sept. 2008), p. 221-248. DOI: 10.3934/dcds.2009.23.221 (cf. p. 18).

[21] François DUBOIS et Pierre LALLEMAND. "Towards higher order lattice Boltzmann schemes". In : *J. Stat. Mech.* 2009 (nov. 2008). DOI: 10.1088/1742-5468/2009/06/P06006 (cf. p. 18).

[22] Tony FÉVRIER. "Extension et analyse des schémas de Boltzmann sur réseau : les schémas à vitesse relative." Thèse de doct. Univeristé de Paris-Saclay, 2015 (cf. p. 24).

[23] R. GATIGNOL. *Théorie cinétique des gaz á répartition dicrète des vitesses*. T. 36. Springer, Berlin, 1975 (cf. p. 18).

[24] Martin GEIER, Andreas GREINER et Jan G. KORVINK. "Cascaded digital lattice Boltzmann automata for high Reynolds number flow". In: *Physical Review E* 73.6 (juin 2006). DOI: 10.1103/physreve.73.066705 (cf. p. 24).
[25] Martin GEIER et al. "The cumulant lattice Boltzmann equation in three dimensions : Theory and validation". In : *Computers and Mathematics with Applications* 70.4 (2015), p. 507-547. ISSN : 0898-1221. DOI : https://doi.org/10.1016/j.camwa.2015.05.001 (cf. p. 21, 24).

[26] Irina GINZBURG, Dominique D'HUMIÈRES et Alexander KUZMIN. "Optimal Stability of Advection-Diffusion Lattice Boltzmann Models with Two Relaxation Times for Positive/Negative Equilibrium". In : *Journal of Statistical Physics* 139.6 (avr. 2010), p. 1090-1143. DOI : 10.1007/s10955-010-9969-9 (cf. p. 24).

[27] H. GRAD. "Note on N-dimensional Hermite polynomials." In : *Communica-tion in Pure and Applied Mathematics* 2 (1949), p. 331 (cf. p. 18).

[28] H. GRAD. "On the Kinetic theory of rarified gases." In : *Communication in Pure and Applied Mathematics* 2 (1949), p. 325 (cf. p. 18).

[29] Zhaoli GUO et Chang SHU. *Lattice Boltzmann Method and Its Applications in Engineering*. WORLD SCIENTIFIC, 2013. DOI : 10.1142/8806. eprint : https://www.worldscientific.com/doi/pdf/10.1142/8806 (cf. p. 13).

[30] HE, X. AND LUO, L.S. "A priori derivation of the lattice Boltzmann equation". In : *Phys. Rev. E* 55 (1997), R6333 (cf. p. 18).

[31] H-B. HUANG, X-Y LU et M.C. SUKOP. "Numerical study of lattice Boltzmann methods for a convection–diffusion equation coupled with Navier–Stokes equations". In : *Journal of Physics A : Mathematical and Theoretical* 44 (2011) (cf. p. 32).

[32] Jerome JACOB, Orestis MALASPINAS et Pierre SAGAUT. "A new hybrid recursive regularised Bhatnagar-Gross-Krook collision model for Lattice Boltzmann method-based large eddy simulation". In : *Journal of Turbulence* 19.11-12 (2018), p. 1051-1076. DOI : 10.1080/14685248.2018.1540879 (cf. p. 22, 33).

[33] Siarhei KHIREVICH, Irina GINZBURG et Ulrich TALLAREK. "Coarse- and fine-grid numerical behavior of MRT/TRT lattice-Boltzmann schemes in regular and random sphere packings". In : *Journal of Computational Physics* 281 (jan. 2015), p. 708-742. DOI : 10.1016/j.jcp.2014.10.038 (cf. p. 23).

[34] Timm KRÜGER et al. *The Lattice Boltzmann Method : Principles and Practice*. Springer International Publishing, 2017. DOI : 10.1007/978-3-319-44649-3 (cf. p. 13).

[35] A. KUZMIN, I. GINZBURG et A.A. MOHAMAD. "The role of the kinetic parameter in the stability of two-relaxation-time advection–diffusion lattice Boltzmann schemes". In : *Computers & Mathematics with Applications* 61.12 (juin 2011), p. 3417-3442. DOI : 10.1016/j.camwa.2010.07.036 (cf. p. 24).

[36] LALLEMAND, P. AND LUO, L.S. "Theory of the lattice Boltzmann method : Dispersion, dissipation, isotropy, Galilean invariance and stability". In : *Ph. Review E* 61(06) (2000) (cf. p. 23).
[37] J. LATT et B. CHOPARD. "Lattice Boltzmann method with regularized pre-collision distribution functions". In: Mathematics and Computers in Simulation 72.2 (2006), p. 165-168 (cf. p. 22).

[38] Jonas LATT, Bastien CHOPARD et Paul ALBUQUERQUE. "Spatial Coupling of a Lattice Boltzmann fluid model with a Finite Difference Navier-Stokes solver". In: arXiv :physics 05112243 (2005) (cf. p. 32).

[39] A. MAGNI et G.-H. COTTET. "Accurate, non-oscillatory remeshing schemes for particle methods". In : *J. Comput. Phys.* 231(1) (2012), p. 152-172 (cf. p. 34).

[40] Orestis MALASPINAS. "Increasing stability and accuracy of the lattice Boltzmann scheme : recursivity and regularization". In : *arXiv e-prints*, arXiv :1505.06900 (mai 2015), arXiv :1505.06900. arXiv : 1505.06900 [physics.flu-dyn] (cf. p. 18, 22).

[41] S. MARIÉ. "On the use of adaptive relaxation times in Lattice Boltzmann Methods". In : *Accepted for publication in Progress in Computational Fluid Dynamics* (2021) (cf. p. 26).

[42] S. MARIÉ et X. GLOERFELT. "Adaptive filtering for the lattice Boltzmann method." In : *Journal of Computational Physics* 333C (2017), p. 212-226 (cf. p. 25, 26).

[43] S. MARIÉ, D. RICOT et P. SAGAUT. "Comparison between Lattice Boltzmann Method and Navier-Stokes high order schemes for Computational Aeroacoustics." In : *J. Comput. Phys.* 228 (4 2009), p. 1056-1070 (cf. p. 31).

[44] M. MEHL, T. NECKEL et Ph. NEUMANN. "Navier-Stokes and Lattice-Boltzmann on octree-like grids in the Peano framework". In: *International Journal for Numerical Methods in Fluids* 65.1 (2011), p. 67-86. DOI: 10.1002/fld.2469. eprint: https://onlinelibrary.wiley.com/doi/pdf/10.1002/fld.2469 (cf. p. 33).

[45] Chloé MIMEAU, Simon MARIÉ et Iraj MORTAZAVI. "A comparison of semi-Lagrangian Vortex method and Lattice Boltzmann method for incompressible flows". In: *Computers & Fluids* 224 (juin 2021), p. 104946. DOI: 10.1016/j.compfluid.2021.104946 (cf. p. 34).

[46] Chloé MIMEAU et Iraj MORTAZAVI. "A Review of Vortex Methods and Their Applications : From Creation to Recent Advances". In : *Fluids* 6 (fév. 2021), p. 68. DOI: 10.3390/fluids6020068 (cf. p. 34).

[47] J. von NEUMANN et R.D. RICHTMYER. "A method for the numerical calculation of hydrodynamics shocks." In : *Journ. of Applied Phys.* 21 (1950), p. 232-237 (cf. p. 25).

[48] Philipp NEUMANN. "On transient hybrid Lattice Boltzmann–Navier-Stokes flow simulations". In : *Journal of Computational Science* 17 (2016). Discrete Simulation of Fluid Dynamics 2015, p. 482-490. ISSN : 1877-7503. DOI : https://doi.org/10.1016/j.jocs.2016.02.003 (cf. p. 32, 33).
[49] Hiroshi OTOMO, Bruce M. BOGHOSIAN et François DUBOIS. "Two complementary lattice-Boltzmann-based analyses for nonlinear systems". In : *Physica A : Statistical Mechanics and its Applications* 486 (nov. 2017), p. 1000-1011. DOI : 10.1016/j.physa.2017.06.010 (cf. p. 18).

[50] D. RICOT, S. MARIÉ et P. SAGAUT. "Lattice Boltzmann Method with selective viscosity filter". In : *Journal of Computational Physics* 228-12 (2009), p. 4478-4490 (cf. p. 25).

[51] P. ROMATSCHKE, Miller MENDOZA et Sauro SUCCI. "Fully relativistic lattice Boltzmann algorithm". In : *Physical Review C - PHYS REV C* 84 (juin 2011). DOI : 10.1103/PhysRevC.84.034903 (cf. p. 15).

[52] Alessandro De ROSIS et Kai H. LUO. "Role of higher-order Hermite polynomials in the central-moments-based lattice Boltzmann framework". In : *Physical Review E* 99.1 (jan. 2019). DOI : 10.1103/physreve.99.013301 (cf. p. 24).

[53] Alessandro De ROSIS. "Nonorthogonal central-moments-based lattice Boltzmann scheme in three dimensions". In : *Physical Review E* (2017) (cf. p. 24).

[54] Ansumali S., Karlin I.V. et Ottinger H.C. "Minimal entropic kinetic models for hydrodynamics". In : *Europhys. Lett.* 63 (6 2003), p. 798-804 (cf. p. 22, 25).

[55] X. SHAN et X. HE. "Discretization of the velocity space in solution of the Boltzmann equation." In : *Phys Rev Lett* 80 (1998), p. 65 (cf. p. 18).

[56] X. SHAN, X. YUAN et H. CHEN. "Kinetic theory representation of hydrodynamics : a way beyond the Navier-Stokes equation." In : *Journal of Fluid Mechanics* 550 (2006), p. 413-441 (cf. p. 18).

[57] W. SHI et W.and Mei R. SHYY. "Finite-Difference-Based Lattice Boltzmann Method for inviscid compressible flows." In : *Numerical Heat Transfer Part B : Fundamentals* 40 (2001), p. 1-21 (cf. p. 20).

[58] Sauro SUCCI. *The Lattice Boltzmann Equation for Fluid Dynamics and Beyond*. Oxford University Press, 2001. ISBN : 978-0-19-850398-9 (cf. p. 13).

[59] TAM, C.K.W. AND WEBB, J.C. "Dispersion relation preserving finite difference schemes for computational acoustics." In : *J. Comput. Phys.* 107 (1993), p. 262-281 (cf. p. 31).

[60] F. TOSI et al. "Numerical stability of Entropic versus positivity-enforcing Lattice Boltzmann schemes". In : *Mathematics and Computers in Simulation* 72 (2006), p. 227-231 (cf. p. 22, 25).

[61] Gauthier WISSOCQ, Pierre SAGAUT et Jean-François BOUSSUGE. "An extended spectral analysis of the lattice Boltzmann method : modal interactions and stability issues". In : *J. Comput. Phys.* 380 (mar. 2019), p. 311-333. DOI : 10.1016/j.jcp.2018.12.015 (cf. p. 30).

[62] WOLF-GLADROW, D.A. *Lattice-Gas Cellular Automata and Lattice Boltzmann Models*. Springer, 2000 (cf. p. 15).
[63] H. XI, G. PEND et S.H. CHOU. "Finite Volume Lattice Boltzmann Method."
In: *Physical Review E* 59 (1999), p. 6202-6205 (cf. p. 20).

[64] Liming YANG et al. *Lattice Boltzmann and Gas Kinetic Flux Solvers*. WORLD SCIENTIFIC, 2020. DOI: 10.1142/11949. eprint: https://www.worldscientific.com/doi/pdf/10.1142/11949 (cf. p. 13).

[65] L. ZHONG et al. "Lattice Boltzmann schemes for the nonlinear Schrödinger equation." In: *Phys. Rev. E*. 74 (2006), p. 036704 (cf. p. 15).
# Deuxième partie

## Etude des écoulements complexes en LBM
<figure><img src="image_13.png" /></figure>
# Chapitre 3

## Etude LBM d'une interface à fort
gradient de densité

La simulation numérique des écoulements diphasiques se caractérise par la difficulté à calculer l'interface entre les deux phases avec précision. En effet, cette frontière entre deux états particuliers est caractérisée par une forte discontinuité (par exemple en masse volumique) et représente donc une difficulté numérique bien connue. Une description précise de l'interface est importante pour d'une part conserver la nature étroite de cette discontinuité mais aussi pour décrire de façon fiable l'évolution de cette interface voire à sa rupture lorsqu'elle est soumise à des efforts importants. Un type de méthodes utilisées pour simuler ce type d'écoulement consiste à ajouter aux équations de conservation une équation de transport relative à un paramètre de suivi d'interface qui peut être continu (level-set method) ou discontinu (VOF).

D'un point de vue mésoscopique les problématiques restent similaires et plussieurs approchent sont possibles comme les approches à pseudo-potentiel de Shan&Chen [32] et les approches à « énergie libre » de Swift [38]. La plupart des approches proposent l'adjonction d'une nouvelle fonction de distribution dédiée au suivi de l'interface dans l'esprit des méthodes level-set. En particulier le modèle de He, Chen et Zhang [17] permet la prise en compte de densités différentes via l'introduction d'une pression thermodynamique régit par une loi d'état. Cependant, lorsque les rapports de densité deviennent importants, des instabilités numériques apparaissent. Dès lors, il devient intéressant de vouloir tester l'efficacité des filtres adaptatifs présentés dans la première partie à ce type d'écoulement avec fort gradient de densité. Le but de ce chapitre sera donc de présenter l'utilisation des filtres adaptatifs dans le cadre de simulation LBM diphasique et de tester son application sur l'instabilité de Rayleigh-Taylor.

### 3.1 Le modèle d'interface de He-Chen-Zhang

Le modèle HCZ [17] considère donc deux fonctions de distributions sur un réseau de vitesses standard de type D2Q9. La première distribution $f_α$ permet de suivre l'évolution de l'interface. Cette distribution peut être considéré comme l'ana-
logue mésoscopique de la fonction d'index $\phi$ utilisée dans les méthodes level-set. La seconde fonction de distribution $g_\alpha$ permet de faire le lien entre l'interface et la pression hydrodynamique. L'évolution de ces distributions est donnée par le terme de collision classique BGK auquel on ajoute les effets des actions intermoléculaires comme les tensions de surfaces ou la gravité.

$$ 
\begin{aligned}
f_{\alpha}^{+} &= \mathbf{BGK}(f_{\alpha}) - \frac{2\tau - 1}{2\tau}(\mathbf{c}_{\alpha} - \mathbf{u}) \cdot \frac{\nabla\psi(\phi)}{RT} \Gamma_{\alpha}(\mathbf{u}) \\
g_{\alpha}^{+} &= \mathbf{BGK}(g_{\alpha}) - \frac{2\tau - 1}{2\tau}(\mathbf{c}_{\alpha} - \mathbf{u}) \cdot [\nabla\psi(\rho)(\Gamma_{\alpha}(\mathbf{u}) - \Gamma_{\alpha}(0)) - \Gamma_{\alpha}(\mathbf{u})(\mathbf{F}_{\mathbf{s}} + \mathbf{G})]
\end{aligned}
\quad (3.1) 
$$

Les états d'équilibres sont alors donnés par :

$$ 
\begin{aligned}
f_{\alpha}^{eq} &= \phi \Gamma_{\alpha}(\mathbf{u}) \\
g_{\alpha}^{eq} &= \omega_{\alpha} \left[ p + \rho RT \left( \frac{\mathbf{u} \cdot \mathbf{c}_{\alpha}}{RT} + \frac{(\mathbf{u} \cdot \mathbf{c}_{\alpha})^2}{2(RT)^2} - \frac{|\mathbf{u}|^2}{2RT} \right) \right]
\end{aligned}
\quad (3.2) 
$$

Avec :

$$ \Gamma_{\alpha}(\mathbf{u}) = \omega_{\alpha} \left( 1 + \frac{\mathbf{u} \cdot \mathbf{c}_{\alpha}}{RT} + \frac{(\mathbf{u} \cdot \mathbf{c}_{\alpha})^2}{2(RT)^2} - \frac{|\mathbf{u}|^2}{2RT} \right) \quad (3.3) $$

Les moments des fonctions $f$ et $g$ permettent alors de retrouver les quantités macroscopiques nécessaires à l'équilibre :

$$ 
\begin{aligned}
\phi &= \sum_{\alpha} f_{\alpha} \\
p_{h} &= \sum_{\alpha} g_{\alpha} - \frac{1}{2} \mathbf{u} \cdot \nabla \psi(\rho) \\
\rho R T \mathbf{u} &= \sum_{\alpha} \mathbf{c}_{\alpha} g_{\alpha} + \frac{RT}{2} (\mathbf{F}_{\mathbf{s}} + \mathbf{G})
\end{aligned}
\quad (3.4) 
$$

Il est alors important de noter ici que la fonction $\psi$ présente dans les relations précédentes est une fonction analytique d'une variable. Le modèle HCZ défini une fonction cubique pour $\psi$ et distingue l'utilisation de $\psi(\phi)$ avec la fonction de distribution $f_\alpha$ et de $\psi(\rho)$ avec la fonction de distribution $g_\alpha$ afin de pouvoir séparer la pression thermodynamique de la pression hydrodynamique :

$$ 
\begin{aligned}
\psi(\phi) &= p_{th} - \phi RT \\
\psi(\rho) &= p_{h} - \rho RT
\end{aligned}
\quad (3.5) 
$$

où la pression thermodynamique est calculée via l'équation d'état de Carnahan-Starling qui permet une pression identique pour deux densités différentes. Cette distinction entre pression hydrodynamique et pression thermodynamique est un aspect très important du modèle HCZ et permet une amélioration notable de la stabilité numérique pour les forts gradients de densité (~ 10).

Cependant des instabilités persistent lorsque le nombre de Reynolds devient important. Dès lors il devient intéressant de vouloir utiliser le filtrage adaptatif présenté en 1.3.2 pour tenter de réduire les instabilités du modèle HCZ à haut Reynolds.
## 3.2 Filtrage adaptatif pour les forts gradients de densité

L'idée d'utiliser le filtrage sélectif pour les écoulements diphasique a été présenté initialement par Kim [24] puis par Chao [9] en utilisant un modèle LBM dérivé du modèle HCZ. Les deux auteurs préconisent l'utilisation de filtres d'ordre élevé pour conserver les propriétés basses fréquences du schéma. Ceci peut-être vu comme une contrainte numérique pouvant augmenter artificiellement le stencil du schéma. Ainsi comme expliqué en 1.3.2 il peut être intéressant de rendre le filtre adaptatif pour localiser son action et permettre l'utilisation des filtres à faible stencil en conservant la sélectivité. Pour utiliser les filtres adaptatifs dans le cadre des écoulements diphasiques, il faut choisir une quantité pertinente pour la sensibilité du filtre. Une idée simple consiste à utiliser le gradient de densité présent dans le terme $\nabla\psi$ des équations 3.1. Le coefficient $\sigma$ devient alors :

$$ \sigma = \left(1 - e^{-(\nabla\psi/\psi_0)^2}\right)^2 \qquad (3.6) $$

$\psi_0$ est alors le seuil que l'on peut localement évaluer par $\psi_0 \sim \frac{\rho_h - \rho_l}{\eta}$ où $\eta$ est une grandeur caractéristique associée au gradient de densité. L'intérêt de choisir ce type de senseur c'est qu'il rajoute très peu de temps de calcul car $\nabla\psi$ est nécessaire à l'évaluation du terme source et peut être évalué par des schémas locaux utilisant les paramètres du réseau :

$$ \nabla\psi = \sum_{\alpha \neq 0} \frac{\omega_{\alpha} \mathbf{c}_{\alpha} \cdot \mathbf{i} [\psi(\mathbf{x} + \mathbf{c}_{\alpha}) - \psi(\mathbf{x} - \mathbf{c}_{\alpha})]}{2c_0^2} \qquad (3.7) $$

On peut alors ensuite filtrer classiquement la quantité $\phi$ qui sera donc affectée aux endroits ou le gradient de densité sera important. L'algorithme devient alors :

- Calcul de $\nabla\psi(\rho)$ avec (3.7)

- Mise à jour de $\sigma_d$

- Étape de collision

- Étape de propagation

- Calcul des moments $(\phi_h, \Gamma)$ Filtragede$\phi$

- Calcul de l'état d'équilibre $g_{\alpha}^{eq}$ et $f_{\alpha}^{eq}$

Ces filtres adaptatifs sur le gradient de densité ont été testé sur une instabilité de Rayleigh-Taylor pour différents rapports de densité et différents nombres de Reynolds.

## 3.3 Application à l'instabilité de Rayleigh-Taylor

L'instabilité de Rayleigh-Taylor consiste à placer initialement un fluide lourd de densité $\rho_h$ au dessus d'un fluide plus léger de densité $\rho_l$ et de déformer initialement l'interface entre les deux fluides pour la déstabiliser. Ce type d'instabilité devient particulièrement difficile à simuler lorsque le rapport de densité et le nombre de Reynolds deviennent important. La dimension caractéristique de cette instabilité
est basée sur la largeur du domaine W qui correspond à la largeur initiale de
la perturbation. L'échelle de vitesse est alors définie par la quantité $\sqrt{Wg}$ où g représente l'accélération gravitationnelle imposée ici a une valeur telle que $\sqrt{Wg} = 0.04$. L'échelle de temps est définie par $T = \sqrt{W/g}$. On peut alors définir le nombre de Reynolds :

$$Re = \frac{\sqrt{WgW}}{\nu} \qquad (3.8)$$

et le rapport de densité est caractérisé par le nombre d'Atwood :

$$A_t = \frac{\rho_h - \rho_l}{\rho_h + \rho_l} \qquad (3.9)$$

Les simulations sont réalisées sur une grille rectangulaire $L_{xy} = W \times 4W$. Les résultats montrent un effet du filtrage sur la stabilité numérique permettant une augmentation notable du Reynolds pour un grand nombre d'Atwood (figure 3.1). Typiquement, on peut obtenir $Re_{max} = 4096$ pour $At = 0.9$ contre $Re_{max} = 632$ pour le modèle HCY original. En revanche, les instabilités persistent lorsque l'on souhaite augmenter le rapport de densité.

D'autre part, nous avons pu dans cette étude, effectuer une comparaison avec des calculs Navier-Stokes afin de vérifier que les simulations filtrées à haut Reynolds conservaient de bonnes propriétés. Les calculs Navier-Stokes ont été effectué à l'aide d'un code basé sur des schémas de type TVD avec un limiteur Super-Bee [30] sur le même maillage et le même CFL. Les résultats montrent une bonne adéquation jusqu'à $Re = 2048$ mais ont permis de mettre en évidence un déficit de vitesse dans les simulations LBM filtrées lorsque le nombre d'Atwood devient important [27].

La figure 3.2 représente l'état de l'interface au même instant pour différentes

<figure><img src="image_14.png" /><figcaption>FIGURE 3.1 – Position de l'interface à t/T ~ 2.2. De gauche à droite : ($A_t$, $R_e$) = (0.8, 1024) (0.8, 2048) (0.9, 1024) et (0.9, 4096). – Navier-Stokes, – LBM.</figcaption></figure>
valeurs de Reynolds et d'Atwood. Les valeurs extrêmes de la figure représentent des cas instables dont aucun calcul de référence n'a pu être établi dans la littérature jusqu'à présent.

Des améliorations du modèle HCZ basées sur un état d'équilibre ou sur l'utilisation de lois d'état différentes ont été réalisées et donnent des résultats prometteurs [7, 33, 37, 25, 44].
<figure><img src="image_15.png" /><figcaption>FIGURE 3.2 – Evolution de l'interface à l'instant $tA_t/T = 2.2$ en fonction des nombres de Reynolds et d'Atwood.</figcaption></figure>
# Chapitre 4

## Modélisation de la diffusion multi-espèces en LBM

Ce chapitre s'inscrit dans la continuité des travaux de recherche entrepris précédemment. Après avoir étudié la LBM dans sa forme canonique à une seule espèce avec les problématiques de stabilité et l'amélioration des modèles de collision, après avoir fait un détour par les écoulements d'arrière corps et tenté d'étudier les écoulements réactifs multi-espèces dans le sillage de jets propulsifs, il m'a semblé intéressant de revenir aux approches mésoscopiques en essayant de comprendre comment les interactions entre espèces pouvaient être modélisées dans les approches mésoscopiques. Ainsi, après la modélisation diphasique du chapitre précédent, nous allons aborder dans ce chapitre la modélisation des écoulements multi-espèces par les approches LBM. Nous allons surtout nous concentrer sur les écoulements miscibles en focalisant notre attention sur la diffusion complexe qui peut survenir en présence de plusieurs espèces. Ce travail est essentiellement issu du travail de thèse de Lucien VIENNE que j'ai pu encadrer entre 2016 et 2019. Après un tour d'horizon des approches possibles en LBM, nous présenterons l'approche que nous avons privilégiée pour ensuite s'intéresser à l'instabilité de Saffman-Taylor.

### 4.1 Modélisation LBM multi-espèces

#### 4.1.1 Position du problème

La principale difficulté dans le cadre de la simulation des écoulements engendrant plusieurs espèces réside dans la description des interactions entre les espèces. En effet, comme nous l'avons vu dans le chapitre 6.2 si le mélange d'espèces peut engendrer des réactions chimiques, il faut prendre en compte les taux de réaction dans les équations de conservations. Si toutefois les espèces sont considérées non réactives, elles peuvent cependant se mélanger à des vitesses différentes. Ce processus de diffusion à plusieurs espèces est très souvent simplifié dans les approches macroscopiques en considérant un mélange stable et établi pour lequel on peut définir une vitesse globale et un coefficient de viscosité global. Cette simplification est tout à fait acceptable lorsque les espèces sont parfaitement mélangées et que
leur mouvement relatif peut se confondre avec un mouvement moyen mais devient discutable dans les zones ou le mélange reste naissant et transitionnel lentement vers un état d'équilibre global. De plus la modélisation de la diffusion en présence de plusieurs espèces doit être prise en compte par une loi de type Maxwell-Stephan qui permet d'étendre la validité de la loi de Fick à plus de deux composants.

La LBM offre une manière alternative et pratique de modéliser les écoulements multi-espèces. L'approche la plus courante consiste à utiliser une fonction de distribution pour chaque espèce. Cependant la principale difficulté est de choisir un opérateur de collision approprié et en particulier un état d'équilibre cohérent. Une approche possible consiste à séparer les collisions entre molécules d'une même espèce (self-collision) et les collisions entre molécules d'espèces différentes (cross-collision) [26, 34, 20, 40]. D'autres auteurs [5, 6, 35, 2, 21] choisissent d'utiliser le formalisme de BGK en reconsidérant l'écriture de l'équilibre et distinguent les phases d'équilibre local des espèces et l'équilibre global du mélange. D'autres approches sont basées sur un seul terme de collision global tel que le modèle d'Asinari [4], dérivé de l'opérateur de collision AAP-BGK [1]. Les références [3, 8, 46] ont légèrement modifié la formulation de l'état d'équilibre et la définition de la vitesse d'équilibre.

L'approche choisie dans le cadre de la thèse de Lucien VIENNE a été de considérer un équilibre basé sur la vitesse propre de chaque espèce et une interaction inter espèces modélisée par un terme de forçage sur un réseau de vitesse isotherme du type D2Q9 ou D3Q19.

### 4.1.2 Utilisation d'un terme de forçage

Un mélange est composé de plusieurs espèces et chaque espèce est définie par sa propre fonction de distribution qui est gouvernée par sa propre équation cinétique. Pour des raisons de simplicité, on considère uniquement un opérateur de collision de type BGK mais les opérateurs présentés dans les sections 1.2.2 et 1.2.3 peuvent aussi être implémentés. Soit les indices $m$ et $n$ correspondant aux différentes espèces ($m, n = 1, 2, ..., N$; $N$ étant le nombre total d'espèces). La fonction de distribution de l'espèce $m$, $f_{\alpha}^{m}$, obéit à l'équation cinétique discrète suivante :

$$f_{\alpha}^{m}(\mathbf{x} + \mathbf{e}_{\alpha}\delta_{t}, t + \delta_{t}) = f_{\alpha}^{m}(\mathbf{x}, t) - \frac{\delta_{t}}{\tau_{m}} [f_{\alpha}^{m}(\mathbf{x}, t) - f_{\alpha}^{m(eq)}(\mathbf{x}, t)] \\ + (1 - \frac{\delta_{t}}{2\tau_{m}})\delta_{t}S_{\alpha}^{m}(\mathbf{x}, t) \quad (4.1)$$

La fonction de distribution à l'équilibre, $f_{\alpha}^{m(eq)}$, est donnée par la formulation polynomiale usuelle

$$f_{\alpha}^{m(eq)} = \rho_{m}\omega_{\alpha}\left[1 + \frac{\mathbf{u}_{m} \cdot \mathbf{e}_{\alpha}}{c_{s}^{2}} + \frac{(\mathbf{u}_{m} \cdot \mathbf{e}_{\alpha})^{2}}{2c_{s}^{4}} - \frac{\mathbf{u}_{m} \cdot \mathbf{u}_{m}}{2c_{s}^{2}}\right]. \quad (4.2)$$

dans laquelle $\mathbf{u}_m$ désigne la vitesse de l'espèce $m$. $S_{\alpha}^m$ est le terme source pour chaque espèce similaire au terme de forçage de Guo [16], couramment utilisé pour
inclure des forces dans l'algorithm LBM,

$$S_{\alpha}^{m} = \omega_{\alpha} \left[ \frac{\mathbf{e}_{\alpha} - \mathbf{u}_{m}}{c_{s}^{2}} + \frac{(\mathbf{e}_{\alpha} \cdot \mathbf{u}_{m})\mathbf{e}_{\alpha}}{c_{s}^{4}} \right] \cdot F_{m}, \quad (4.3)$$

$F_m$ étant la force agissant sur la m-ième espèce qui est dérivée dans la suite pour prendre en compte les interactions entre les espèces.

Les quantités macroscopiques, notamment la densité et la quantité de chaque espèce, sont obtenues en calculant les différents moments des fonctions de distribution.

$$\rho_m = \sum_{\alpha} f_{\alpha}^{m}, \quad \rho_m \mathbf{u}_m = \sum_{\alpha} f_{\alpha}^{m} \mathbf{e}_{\alpha} + \frac{\delta_t}{2} F_m. \quad (4.4)$$

Afin de prendre en compte l'interaction entre les différentes espèces miscibles, nous introduisons les forces de diffusions ou forces de friction selon l'approche de Maxwell-Stefan au transfert de masse :

$$F_{D,m} = -p \sum_{n=1}^{N} \frac{x_m x_n}{D_{mn}} (\mathbf{u}_m - \mathbf{u}_n), \quad (4.5)$$

et $F_m$ devient

$$F_m = F_{D,m} + F_{B,m}, \quad (4.6)$$

où $F_{B,m}$ est une force volumique. Ainsi les équations cinétiques discrètes (4.1) pour les différentes espèces sont couplées par l'intermédiaire de $F_m$. Puisque les forces de diffusion $F_{D,m}$ dépendent de la vitesse, de la pression totale $p$, des fraction molaires $x_m$, et des coefficients de diffusion (Maxwell-Stefan) $D_{mn}$, un système linéaire doit être résolu à chaque pas de temps afin de calculer la quantité de mouvement de chaque espèce grâce à l'équation (4.4). Cette force également appelée force de friction intermoléculaire dépend de la vitesse relative de chaque espèce. L'idée d'inclure les effets de diffusion comme une force agissant sur les molécules date des premiers travaux de la théorie cinétique des gaz par Maxwell [29]. La même expression a été rigoureusement dérivée plus tard par Chapman et Cowling [10], une discussion détaillée est donnée par Hirschfelder, Curtiss, et Bird [18] et Kerkhof et Geboers présentent une dérivation plus récente dans la référence [22].

Les coefficients de viscosité du modèle sont recalculés en fonction de leur concentration et permettent de fixer les temps de relaxation associés à chaque espèce. Les coefficients de diffusion de (4.5) sont également recalculés en fonction des masses molaires [41].

### 4.1.3 Validation du modèle sur des cas tests diffusifs et convectifs

Le modèle proposé a été testé et validé sur plusieurs cas tests documentés. En particulier le tube de Loschmidt dans lequel deux mélanges de même composition mais avec différente proportion sont mis en contact, a permis de valider la diffusion complexe simulée dans le cas de composants ayant une masse molaire différente.
<figure><img src="image_16.png" /><figcaption>FIGURE 4.1 – (Gauche) Evolution des fractions molaires et de la vitesse sur une ligne médiane de deux jets opposés. (Droite) Evolution temporelle des fractions molaires de l'Argon, l'Hydrogène et du Méthane dans le tube de Loschmidt. Les symboles représentent les points expérimentaux [2].</figcaption></figure>

D'autre part, la simulation de la collision de deux jets opposés composés de 4 espèces a permis de valider la partie convective du modèle ainsi que le calcul des coefficients de transports en comparant les résultats à ceux obtenus dans le cadre d'une tabulation des coefficients (Cf figure 4.1).

## 4.2 Application à l'instabilité de Saffman-Taylor

L'instabilité de Saffman-Taylor se produit dans un milieu poreux lorsqu'un fluide moins visqueux s'écoule vers un fluide plus visqueux. Lorsque les deux fluides ne se mélangent pas (non miscibles) l'instabilité est essentiellement pilotée par ses caractéristiques convectives. En revanche si les fluides en présence sont miscibles, une nouvelle échelle caractéristique intervient liée à la diffusion d'un fluide dans l'autre. Cette nouvelle échelle peut s'avérer plus complexe en présence de mélange pré-établi à plusieurs espèces. La simulation de ce type d'instabilité en LBM requiert une modélisation fiable d'un milieu poreux ainsi que la prise en compte de la diffusion complexe de type Maxwell-Stefan.

### 4.2.1 Modélisation d'un milieu poreux

La modélisation d'un milieu poreux en LBM a fait l'objet de nombreuses études et passe par la prise en compte de la porosité soit sous forme de conditions aux limites directes de type Bounce-Back soit sous la forme d'une force externe faisant intervenir les grandeurs associées aux milieux poreux.

La première approche utilise une formulation modifiée du Bounce-Back consistant à ne réfléchir qu'une partie des distributions. Cette approche appelée Gray-LBM (GLBM) a été développée sous plusieurs formes en fonction de la quantité réfléchie. Walsh propose [42] de réfléchir les distributions pré-collision, Zhu [45] les distributions post-collision, Yoshida et Hayashi [43] propose les distributions post-collision des points voisins. Toutes ces formulations visent à reproduire le
comportement décrit par les équations de Darcy-Brinkman [11] construisant le lien entre gradient de pression et gradient de vitesse dans le milieu poreux. La deuxième approche, appelée modèle de force de Brinkman (ou Brinkman Force), consiste à prendre en compte la résistance du milieu poreux par un terme de forçage supplémentaire[31, 13, 15, 28, 36, 14] :

$$ \mathcal{F}_{porous,m} = -\frac{\mu_m}{K} \mathbf{u}_m. \qquad (4.7) $$

Cette stratégie est préférée ici car elle est locale et a un coût de calcul négligeable puisque nous implémentons déjà un terme source en raison des forces de frottement intermoléculaires. Comme mentionné dans les références [39, 23, 12], la perméabilité numérique dépend légèrement de la viscosité pour les modèles de Brinkman et Gray-LBM. Cette variation non physique peut être atténuée en garantant constant le paramètre Λ (Cf, 1.2.3.2). Nous utilisons ici le modèle TRT avec $\Lambda = 3/16$ correspondant à $s_q = 8^{2-1/\tau_m/8-1/\tau_m}$.

<figure><img src="image_17.png" /><figcaption>FIGURE 4.2 – (Gauche) Évolution des longueurs interfaciales ($l_{inter}$) et de mélange ($l_{mix}$) de l'instabilité. (Droite) Évolution de l'interface à $t^* = 0.15, 0.46, 1.07, 3.05, 6.10$.</figcaption></figure>

### 4.2.2 Effet du nombre de Péclet

Le nombre de Péclet ($P_e = UL/D$) permet de caractériser l'importance relative de la convection et de la diffusion. C'est donc un nombre important pour l'étude de la dynamique de ce type d'instabilité pour laquelle les deux phénomènes sont en concurrence. Par exemple, les taux de croissances observés montrent une dépendance au nombre de Péclet permettant de confirmer une valeur limite en dessous de laquelle l'instabilité n'est pas observée. La figure 4.3 présente l'évolution de l'instabilité avec un premier régime diffusif d'autant plus long que le Péclet est grand et un régime convectif qui permet la croissance rapide des franges instables. La courbe neutre permet également d'observer l'enrichissement fréquentiel attendu lorsque le nombre de Péclet augmente.

### 4.2.3 Instabilité à trois espèces

Des simulations ont également été réalisées en présence d'une espèce « catalytique » permettant la déstabilisation d'un cas défavorable par le mécanisme
<figure><img src="image_18.png" /><figcaption>FIGURE 4.3 – [Left] : Évolution temporelle de la longueur de mélange normalisée pour différents nombres de Péclet. [Right] : Courbe neutre pour différents Péclet</figcaption></figure>

de retro-diffusion. L'objectif était de montrer que l'instabilité peut se déclencher même dans un cas défavorable (la viscosité du premier mélange étant plus grande). Dans ce cas, nous avons montré que le modèle était capable de décrire correctement la rétro-diffusion et ainsi de prédire l'apparition d'une instabilité lorsque le transfert de viscosité était suffisant pour porter la viscosité du premier mélange en dessous de celle du second mélange. Pour ce dernier cas, il a été possible d'identifier 3 mécanismes distincts :

- La rétro-diffusion permettant le transfert de fluide moins visqueux dans le premier mélange,

- L'apparition du régime diffusif lié à la nouvelle viscosité du premier mélange,

- Le régime convectif contenant la forte croissance des franges instables.

<figure><img src="image_19.png" /><figcaption>FIGURE 4.4 – Courbe de dispersion des 3 espèces de $t^* = 0.005$ à $t^* = 0.1$ avec un $\Delta t^* = 0.05$ entre chaque courbe (le temps est représenté par le gradient de couleur). ($n_x = 4000$, $n_y = 4000$).</figcaption></figure>

Les courbes neutres réalisées pour les 3 espèces (figure 4.4), montrent une forte croissance de la troisième espèce qui joue le rôle de « catalyseur » et se diffuse rapidement dans le premier mélange. Les autres espèces ont d'abord un taux de croissance négatif qui converge ensuite vers une valeur unique pour les 3 espèces.
Ce type de comportement a été observé par Hota et Mishra [19] pour un profil de viscosité non monotone en ayant recours à une analyse non modale.

4.2.4 Simulations à l’échelle d’un pore

Si l'on considère qu'un milieu poreux peut se caractériser par son échelle ca-
ractéristique *p*, les simulations réalisées précédemment supposaient $\Delta x \gg p$ et
prenait en compte le milieu poreux par l'intermédiaire d'un modèle (Brinkmann
Force). On peut cependant se passer de modèle si l'on représente le milieu poreux
explicitement en imposant $\Delta x < p$.

<figure><img src="image_20.png" /><figcaption>FIGURE 4.5 – *Simulation explicite de la digitation visqueuse* (*p* = 16Δ*x*)</figcaption></figure>

Les premiers résultats (figure 4.5) montrent que l’instabilité est bien observée lorsque le rapport de viscosité est favorable (ici R=3). Dans le cas étudié, le milieu poreux est constitué de pores cylindriques régulièrement et symétriquement espacés. Pourtant l’instabilité se développe et les premières franges observées sont très nettement non-symétriques. Ces premiers résultats devront être approfondis et seront l’objet de futures investigations numériques.
Références

[1] Pierre ANDRIES, Kazuo AOKI et Benoit PERTHAME. In : *Journal of Statistical Physics* 106.5/6 (2002), p. 993-1018. DOI : 10.1023/a:1014033703134 (cf. p. 50).

[2] S. ARCIDIACONO et al. "Lattice Boltzmann model for the simulation of multicomponent mixtures". In : *Phys. Rev. E* 76 (4 oct. 2007), p. 046703. DOI : 10.1103/PhysRevE.76.046703 (cf. p. 50, 52).

[3] P. ASINARI. "Lattice Boltzmann scheme for mixture modeling : Analysis of the continuum diffusion regimes recovering Maxwell-Stefan model and incompressible Navier-Stokes equations". In : *Physical Review E* 80 (5 2009), p. 056701. DOI : 10.1103/PhysRevE.80.056701 (cf. p. 50).

[4] Pietro ASINARI. "Multiple-relaxation-time lattice Boltzmann scheme for homogeneous mixture flows with external force". In : *Phys. Rev. E* 77 (5 mai 2008), p. 056706. DOI : 10.1103/PhysRevE.77.056706 (cf. p. 50).

[5] Pietro ASINARI. "Semi-implicit-linearized multiple-relaxation-time formulation of lattice Boltzmann schemes for mixture modeling". In : *Phys. Rev. E* 73 (5 mai 2006), p. 056705. DOI : 10.1103/PhysRevE.73.056705 (cf. p. 50).

[6] Pietro ASINARI et Li-Shi LUO. "A consistent lattice Boltzmann equation with baroclinic coupling for mixtures". In : *Journal of Computational Physics* 227.8 (avr. 2008), p. 3878-3895. DOI : 10.1016/j.jcp.2007.12.001 (cf. p. 50).

[7] Amir BANARI, Christian F. JANSSEN et Stephan T. GRILLI. "An efficient lattice Boltzmann multiphase model for 3D flows with large density ratios at high Reynolds numbers". In : *Computers Mathematics with Applications* 68.12, Part A (2014), p. 1819-1843. ISSN : 0898-1221. DOI : https://doi.org/10.1016/j.camwa.2014.10.009 (cf. p. 47).

[8] Sam BENNETT, Pietro ASINARI et Paul J. DELLAR. "A lattice Boltzmann model for diffusion of binary gas mixtures that includes diffusion slip". In : *International Journal for Numerical Methods in Fluids* 69.1 (mar. 2011), p. 171-189. DOI : 10.1002/fld.2549 (cf. p. 50).

[9] Jianghui CHAO et al. "A filter-based, mass-conserving lattice Boltzmann method for immiscible multiphase flows". In : *International Journal for Numerical Methods in Fluids* 66.5 (2011), p. 622-647. ISSN : 1097-0363. DOI : 10.1002/fld.2276 (cf. p. 45).
[10] S. CHAPMAN et T. G. COWLING. *The Mathematical Theory of Non-uniform Gases*. 3e éd. Cambridge University Press, 1970. ISBN : 052140844X (cf. p. 51).

[11] Chen CHEN et al. "Chapman-Enskog Analyses on the Gray Lattice Boltzmann Equation Method for Fluid Flow in Porous Media". In: *Journal of Statistical Physics* 171.3 (mar. 2018), p. 493-520. DOI: 10.1007/s10955-018-2005-1 (cf. p. 53).

[12] Irina GINZBURG. "Comment on An improved gray Lattice Boltzmann model for simulating fluid flow in multi-scale porous media: Intrinsic links between LBE Brinkman schemes". In : *Advances in Water Resources* 88 (fév. 2016), p. 241-249. DOI: 10.1016/j.advwatres.2014.05.007 (cf. p. 53).

[13] Irina GINZBURG. "Consistent lattice Boltzmann schemes for the Brinkman model of porous flow and infinite Chapman-Enskog expansion". In : *Physical Review E* 77.6 (juin 2008). DOI: 10.1103/physreve.77.066704 (cf. p. 53).

[14] Patrick GROSFILS et al. "Structural and dynamical characterization of Hele-Shaw viscous fingering". In : *Philosophical Transactions of the Royal Society of London. Series A : Mathematical, Physical and Engineering Sciences* 362.1821 (août 2004). Sous la dir. de M. AL-GHOUL, J. P. BOON et P. V. COVENEY, p. 1723-1734. DOI: 10.1098/rsta.2004.1398 (cf. p. 53).

[15] Zhaoli GUO et T. S. ZHAO. "Lattice Boltzmann model for incompressible flows through porous media". In : *Physical Review E* 66.3 (sept. 2002). DOI: 10.1103/physreve.66.036304 (cf. p. 53).

[16] Zhaoli GUO, Chuguang ZHENG et Baochang SHI. "Discrete lattice effects on the forcing term in the lattice Boltzmann method". In : *Phys. Rev. E* 65 (4 avr. 2002), p. 046308. DOI: 10.1103/PhysRevE.65.046308 (cf. p. 50).

[17] X. HE, S. CHEN et R ZHANG. "A Lattice Boltzmann Scheme for Incompressible Multiphase Flow and Its Application in Simulation of Rayleigh-Taylor Instability". In : *Journal of Computational Physics* 152 (1999), p. 642-663 (cf. p. 43).

[18] J. O. HIRSCHFELDER, C. F. CURTISS et R. B. BIRD. *Molecular Theory of Gases and Liquids*. Second printing. John Wiley & Sons, 1964. ISBN : 9780471400653 (cf. p. 51).

[19] Tapan Kumar HOTA et Manoranjan MISHRA. "Non-modal stability analysis of miscible viscous fingering with non-monotonic viscosity profiles". In : *Journal of Fluid Mechanics* 856 (2018), 552-579. DOI: 10.1017/jfm.2018.674 (cf. p. 55).

[20] Abhijit S JOSHI et al. "Lattice Boltzmann method for continuum, multi-component mass diffusion in complex 2D geometries". In : *Journal of Physics D : Applied Physics* 40.9 (avr. 2007), p. 2961-2971. DOI: 10.1088/0022-3727/40/9/044 (cf. p. 50).
[21] Jinfen KANG, Nikolaos I. PRASIANAKIS et John MANTZARAS. "Thermal multicomponent lattice Boltzmann model for catalytic reactive flows". In : *Phys. Rev. E* 89 (6 juin 2014), p. 063310. DOI : 10.1103/PhysRevE.89.063310 (cf. p. 50).

[22] Piet J. A. M. KERKHOF et Marcel A. M. GEBOERS. "Toward a unified theory of isotropic molecular transport phenomena". In : *AIChE Journal* 51.1 (2004), p. 79-121. DOI : 10.1002/aic.10309 (cf. p. 51).

[23] Siarhei KHIREVICH, Irina GINZBURG et Ulrich TALLAREK. "Coarse- and fine-grid numerical behavior of MRT/TRT lattice-Boltzmann schemes in regular and random sphere packings". In : *Journal of Computational Physics* 281 (jan. 2015), p. 708-742. DOI : 10.1016/j.jcp.2014.10.038 (cf. p. 53).

[24] S.H. KIM et H. PITSCH. "On the lattice Boltzmann method for multiphase flows". In : *Center for Turbulence Research, Annual Research Briefs* (2009), p. 377-388 (cf. p. 45).

[25] You LI et al. "An interfacial lattice Boltzmann flux solver for simulation of multiphase flows at large density ratio". In : *International Journal of Multiphase Flow* 116 (2019), p. 100-112. ISSN : 0301-9322. DOI : https://doi.org/10.1016/j.ijmultiphaseflow.2019.04.006 (cf. p. 47).

[26] Li-Shi LUO et Sharath S. GIRIMAJI. "Theory of the lattice Boltzmann method : Two-fluid model for binary mixtures". In : *Phys. Rev. E* 67 (3 mar. 2003), p. 036302. DOI : 10.1103/PhysRevE.67.036302 (cf. p. 50).

[27] S. MARIÉ et V. DARU. "On the use of adaptive filtering with multiphase lattice Boltzmann models." In : 13th ICMMES conference. Hamburg, Germany, 2016 (cf. p. 46).

[28] Nicos S. MARTYS. "Improved approximation of the Brinkman equation using a lattice Boltzmann method". In : *Physics of Fluids* 13.6 (juin 2001), p. 1807-1810. DOI : 10.1063/1.1368846 (cf. p. 53).

[29] James Clerk MAXWELL. "IV. On the dynamical theory of gases". In : *Philosophical Transactions of the Royal Society of London* 157 (1867), p. 49-88. DOI : 10.1098/rstl.1867.0004 (cf. p. 51).

[30] A. MONGRUEL et al. "Early post-impact time dynamics of viscous drops onto a solid dry surface". In : *Physics of Fluids* 21.3 (2009) (cf. p. 46).

[31] N. RAKOTOMALALA, D. SALIN et P. WATZKY. "Miscible displacement between two parallel plates : BGK lattice gas simulations". In : *Journal of Fluid Mechanics* 338 (mai 1997), p. 277-297. DOI : 10.1017/s0022112097004928 (cf. p. 53).

[32] X. SHAN et H. CHEN. "Lattice Boltzmann model for simulating flows with multiple phases and components". In : *Physical Review E* 47 (1993), p. 1815 (cf. p. 43).

[33] J.Y. SHAO et al. "Free-energy-based lattice Boltzmann model for the simulation of multiphase flows with density contrast". In : *Physical Review E* 89 (2014), p. 033309 (cf. p. 47).
[34] Lawrence SIROVICH. "Kinetic Modeling of Gas Mixtures". In: *Physics of Fluids* 5.8 (1962), p. 908. DOI: 10.1063/1.1706706 (cf. p. 50).

[35] Victor SOFONEA et Robert F. SEKERKA. "BGK models for diffusion in isothermal binary fluid systems". In: *Physica A: Statistical Mechanics and its Applications* 299.3-4 (oct. 2001), p. 494-520. DOI: 10.1016/S0378-4371(01)00246-1 (cf. p. 50).

[36] Michael A. A. SPAID et Frederick R. PHELAN. "Lattice Boltzmann methods for modeling microscale flow in fibrous porous media". In: *Physics of Fluids* 9.9 (sept. 1997), p. 2468-2474. DOI: 10.1063/1.869392 (cf. p. 53).

[37] Christopher D. STILES et Yongqiang XUE. "High density ratio lattice Boltzmann method simulations of multicomponent multiphase transport of H2O in air". In: *Computers & Fluids* 131 (2016), p. 81-90 (cf. p. 47).

[38] M.R. SWIFT et al. "Lattice Boltzmann simulations of liquid-gaz and binary fluid systems". In: *Physical Review E* 54 (1996), p. 5041 (cf. p. 43).

[39] L. TALON et al. "Assessment of the two relaxation time Lattice-Boltzmann scheme to simulate Stokes flow in porous media". In: *Water Resources Research* 48.4 (avr. 2012). DOI: 10.1029/2011WR011385 (cf. p. 53).

[40] Zi-Xiang TONG et al. "A multi-component lattice Boltzmann method in consistent with Stefan-Maxwell equations". In: *Computers & Fluids* 105 (déc. 2014), p. 155-165. DOI: 10.1016/j.compfluid.2014.09.022 (cf. p. 50).

[41] Lucien VIENNE. "Simulation of multi-component flows by the lattice Boltzmann method and application to the viscous fingering instability". Thèse de doct. Conservatoire National des Arts et Métiers, 2019 (cf. p. 51).

[42] Stuart D.C. WALSH, Holly BURWINKLE et Martin O. SAAR. "A new partial-bounceback lattice Boltzmann method for fluid flow through heterogeneous media". In: *Computers & Geosciences* 35.6 (juin 2009), p. 1186-1193. DOI: 10.1016/j.cageo.2008.05.004 (cf. p. 52).

[43] Hiroaki YOSHIDA et Hidemitsu HAYASHI. "Transmission-Reflection Coefficient in the Lattice Boltzmann Method". In: *Journal of Statistical Physics* 155.2 (fév. 2014), p. 277-299. DOI: 10.1007/s10955-014-0953-7 (cf. p. 52).

[44] Ye ZHOU et al. "Rayleigh-Taylor and Richtmyer-Meshkov instabilities". In: *Physica D: Nonlinear Phenomena* (2021), p. 132838. ISSN: 0167-2789. DOI: https://doi.org/10.1016/j.physd.2020.132838 (cf. p. 47).

[45] Jiujiang ZHU et Jingsheng MA. "An improved gray lattice Boltzmann model for simulating fluid flow in multi-scale porous media". In: *Advances in Water Resources* 56 (juin 2013), p. 61-76. DOI: 10.1016/j.advwatres.2013.03.001 (cf. p. 52).
[46] Jens ZUDROP et al. "A robust lattice Boltzmann method for parallel simulations of multicomponent flows in complex geometries". In: *Computers & Fluids* 153 (août 2017), p. 20-33. DOI: 10.1016/j.compfluid.2017.04.021 (cf. p. 50).
Troisième partie

Etudes aérodynamiques
instationnaires
<figure><img src="image_21.png" /></figure>
# Chapitre 5

## Etude des écoulements d'arrière corps

En phase transsonique de vol (entre mach 0.5 et mach 1.5), le lanceur Ariane V est soumis à des efforts aérodynamiques instationnaires particulièrement importants sur la partie arrière (buffeting de culot)[14, 5]. Ces efforts peuvent être critiqués pour la tenue mécanique du divergent du moteur Vulcain. La compréhension des écoulements d'arrière corps est donc fondamentale pour d'une part améliorer la maîtrise des spécifications de vol (prédiction) et d'autre part, pour définir des moyens de réduction du buffeting de culot (anticipation//contrôle). Les problématiques liées au buffeting sont d'autant plus critiques que le moteur Vulcain 2 du lanceur Ariane 5+ a une tuyère allongée pour améliorer sa poussée. [30, 11, 33]

Les approches numériques [12, 4] et expérimentales [13] pour tenter de comprendre les phénomènes mis en jeux, ont d'abord fait l'objet d'études sur une géométrie représentative mais grandement simplifiée comme la marche axisymétrique. Nous allons voir dans la section 5.1 que cette simplification peut permettre de dégager des principes généraux de ce type d'écoulement mais nous verrons ensuite dans la section 5.2 que les effets géométriques sont importants et qu'il est aujourd'hui capital d'étudier ce type d'écoulement en configuration réelle.

### 5.1 Etude des efforts instationnaires sur géométrie simplifiée

Les efforts aérodynamiques représentent des quantités globales directement mesurables et dont les effets sont visibles à grande échelle. Dans notre cas, c'est la tuyère du moteur qui est directement impactée, et si on décide de simplifier la géométrie, on parlera plutôt de l'émergence. Pourtant la compréhension des phénomènes sous-jacents passe nécessairement par l'étude des effets locaux responsables des effets à grande échelle. En ingénierie, on se borne souvent à mesurer ou à calculer un coefficient de proportionnalité entre les effets locaux (par exemple la pression) et leurs conséquences globales (par exemple la trainée). Cependant l'étude fine de la dynamique des grandeurs locales permet une compréhension bien
plus fondamentale des mécanismes en présence.

5.1.1 Reconstruction des efforts instationnaires

Pour étudier les efforts aérodynamiques sur la géométrie simple axisymétrique présentée sur la figure 5.1.1, on peut partir de la formulation générale :

$$
\mathbf{F}(t) = \int_{0}^{2\pi} \int_{0}^{L} p(x, \phi, t) r(x) \, \mathbf{n} dx \, d\phi \equiv \sum_{i=1}^{N_c} p_i(t) \mathbf{dS}_i \quad (5.1)
$$

Dans laquelle $p$ représente la pression et $r$ le rayon de l'émergence à une position donnée (ici $r(x) = R = Cste$). La pression est alors mesurée expérimentalement à l'aide de capteurs instationnaires ou calculée numériquement.

<figure><img src="image_22.png" /><figcaption>FIGURE 5.1 – (a) Géométrie axysymétrique. $D = 0.1$ m, $L = 1.2$ D, $d = 0.4$ D $U_0 = 237$ m/s. (b) Phénomène important de l'écoulement : 1. Couche de mélange, 2. Zone de recirculation, 3. Point de recollement, 4. Recirculation fluide, 5. Écoulement de coin, 6. Sillage turbulent.</figcaption></figure>

On peut alors s'intéresser au contenu spectral de ces efforts pour caractériser
les fréquences dominantes. Pour cela, il suffit de choisir une répartition de points
de mesures de la pression (position des capteurs) et d'estimer correctement les
éléments de surface de la relation 5.1. Dans le cadre de simulations numériques
cette répartition correspond aux points de maillages qui est souvent très fine. En
revanche, la difficulté réside dans la précision des schémas utilisés pour résoudre les
équations de conservation. Les approches utilisées ici sont basées sur la ZDES et
développées par l'équipe de Sébastien DECK à l'ONERA [10, 34]. Dans le cadre
expérimental, en revanche, le nombre de capteurs reste limité pour des raisons
d'encombrement et le choix du positionnement devient capital. Dès lors, il devient
intéressant de vouloir reconstruire les efforts en minimisant le nombre de capteurs
utilisés.
5.1.2 Reconstruction optimisée

Pour cela on peut décomposer la pression en mode de Fourier à condition
de satisfaire la périodicité spatiale nécessaire à cette décomposition. C'est le cas
ici pour une géométrie axisymétrique. Cette approche proposé dans le cadre des
charges latérales internes [23, 10], a été adapté ici pour les efforts externes. On
écrit donc la pression :

$$
\bar{p}(\phi, t) = a_0 + \sum_{k=1}^{\infty} [a_k(t)\cos(k\phi) + b_k(t)\sin(k\phi)] \quad (5.2)
$$

et on regarde l'évolution des coefficients $a_m$ et $b_m$. Pour une géométrie axisymé-
trique, la dynamique est pilotée par les premiers coefficients en raison du caractère
symétrique de l'écoulement. Pour les écoulements internes, le premier mode $m = 1$
peut être reconstruit en utilisant uniquement les pressions diamétralement oppo-
sées comme proposé par Nguyen [23]. En revanche pour les écoulements externes,
supposer que ce mode est uniformément réparti sur toute la longueur de l'émer-
gence conduit à une forte surestimation des efforts. Ainsi en regardant la cohérence
des signaux de pressions opposées (figure 5.2), on s'aperçoit que ce mode est par-
ticulièrement efficace sur une partie très localisée de l'émergence correspondant en
moyenne à la zone de recollement des couches de mélange issues du culot.

<figure><img src="image_23.png" /><figcaption>FIGURE 5.2 – Maximum du module de cohérence entre des points diamétralement opposés situés à deux positions azimutales : (gauche) $\phi_0 = 0$ et (droite) $\phi_0 = \pi/2$.</figcaption></figure>

Ainsi, les efforts peuvent être reconstruits en utilisant uniquement des cap-
teurs situés dans la zone d'intérêt mise en évidence sur la figure 5.2. Cette ap-
proche permet de réduire considérablement le nombre d'informations nécessaires
à la reconstruction d'un phénomène complexe.

Des auteurs ont récemment regardé plus précisément la dynamique de l'écou-
lement [24, 25] en utilisant une décomposition dynamique (DMD[27]) basée sur
une décomposition temporelle des coefficients de (5.3), ou en utilisant des mesures
PIV plus résolues [31]. Ces études ont permis d'identifier une structure en hélice
responsable du mode m = 1 au culot d'une géométrie axisymétrique.
<figure><img src="image_24.png" /><figcaption>FIGURE 5.3 – Comparaison des DSP d'efforts reconstruits avec le modèle Nguyen de base et le modèle optimisé.</figcaption></figure>

## 5.2 Etude des configurations réelles

Les problématiques liées au buffeting de culot sur le lanceur Ariane5 ont poussé les industriels européens à étudier en profondeur les paramètres importants dans ce type de phénomène. En particulier les effets de géométrie ont été très vite pointés du doigt. Dans le cadre du TRP entre l'ESA, le CNES, l'ONERA et ArianGroup [19, 17, 9], des calculs ZDES ont montré que l'écoulement sur une configuration dépourvue d'attaches arrières (DAAR), engendrait des efforts latéraux bien moins importants [32]. Ainsi une campagne expérimentale visant à caractériser les efforts sur deux types de configuration (figure 5.4) a été lancée.

<figure><img src="image_25.png" /><figcaption>FIGURE 5.4 – Maquette Ariane 5 à l'échelle 1/60 dans la veine du NLR.</figcaption></figure>

Une façon plus systématique d'étudier les pressions en parois dans le cas de géométries quelconques, consiste à calculer la matrice des corrélations point à point et de regarder s'il est possible d'en dégager une base particulière qui permettrait une représentation plus efficace des phénomènes. En pratique, cela revient à écrire la relation (5.2) sans supposer de base particulière [21, 29] :
<figure><img src="image_26.png" /><figcaption>FIGURE 5.5 – Corrélation entre les 4 premiers modes POD de la vitesse et les efforts suivant z reconstruits avec les 2 premiers modes de pression.</figcaption></figure>

$$p(\mathbf{X}, t) = \sum_{k=1}^{\infty} a_k(t) \phi_k(\mathbf{X}), \qquad (5.3)$$

L'explicitation des $a_k$ et $\phi_k$ se fait alors par la recherche des valeurs et vecteurs propres de la matrice des corrélations. Cette procédure est souvent appelée POD (Proper Orthogonal Decomposition) en mécanique des fluides et correspond thématiquement à une analyse en composantes principales.

Un des principaux résultats dans le cadre de la reconstruction des efforts latéraux sur les géométries choisies, est que le développement (5.3) n'est nécessaire que jusqu'à $k = 2$. En effet, en analysant la forme des $\phi_k$ il est apparu que sur l'extrémité de la tuyère (zone la plus en aval) on pouvait identifier $\phi_1(\mathbf{X}) \sim \cos(\phi + \delta)e_x + \sin(\phi + \delta).e_y$ et $\phi_2(\mathbf{X}) \sim \cos(2\theta + \delta)e_x + \sin(2\theta + \delta).e_y$, le déphasage $\delta$ étant dû à la dissymétrie des attaches. Cette identification devient globalement fausse en amont de l'arrière corps mais a permis de comprendre que les deux premiers modes POD de la pression pariétale permettait d'expliquer le comportement des deux modes de Fourier ($m = 1$ et $m = 2$) comme le montre les DSP de la figure 5.6.

D'autre part, les essais réalisés sur les deux configurations ont permis de confirmer les prédicitions du calcul ZDES qui avait mis en évidence l'impact des attaches arrières sur les efforts latéraux. Pour mieux comprendre les relations entre l'écoulement et sa signature de pression sur la tuyère, des mesures PIV résolues en temps ont été réalisées dans la soufflerie du NLR dans le plan normal du lanceur (XZ) [28].

Les analyses POD du champ de vitesse ont montré que la configuration A
<figure><img src="image_27.png" /><figcaption>FIGURE 5.6 – Contribution des modes POD aux efforts latéraux.</figcaption></figure>

présentait d'avantage de structures cohérentes rendant l'écoulement plus efficace et se traduisant par des efforts plus importants lors du recollement sur la tuyère. De plus, la corrélation entre les modes de vitesse et les modes de pression (figure 5.5) semble indiquer que les structures cohérentes générées par la présence des attaches sont responsables des efforts observés sur la tuyère. Des techniques de reconstruction de la pression à partir des données PIV [31] pourront permettre d'approfondir les relations entre le champ de pression de sillage et la signature en paroi.
# Chapitre 6

## Etude des écoulements à thermodynamique complexe

Ce chapitre présente les études réalisées dans le cadre d'une meilleure compréhension des mécanismes intervenant dans les écoulements d'arrière-corps en présence de jets propulsifs. La première étude se place dans le cadre d'une couche de mélange faisant intervenir un rapport de densité et la seconde présente les résultats obtenus au cours de la thèse de Lorris CHARRIER réalisée en collaboration avec ArianeGroup sur la modélisation hybride RANS/LES des écoulements réactifs.

### 6.1 Transition d'une couche de mélange superso-nique

Cette étude est issue d'une collaboration avec l'université de la Sapienza à Rome avec Matteo Bernardini et Sergio Pirozzoli. Le but est d'étudier la zone de transition d'une couche de mélange issue de deux couches limites compressibles en reproduisant numériquement la configuration expérimentale de Goebel & Dutton [18] avec un nombre de Mach convectif de 0.45. Pour les simulations, le code utilisé est basé sur des différences finies d'ordre 6 avec un WENO d'ordre 5 et une intégration Runge-Kutta d'ordre 4 [26, 20]. La couche de mélange est initialisée à partir des deux couches limites qui se développent spatialement sur une plaque plane explicitement prises en compte en utilisant les techniques de *recycling/rescaling* pour obtenir des couches limites turbulentes en sortie de plaque plane. Le calcul est dimensionné pour résoudre les structures associées à un nombre de Reynolds de 14000 et 17000 dans les deux couches limites et comporte ainsi plus de 1.4 milliard de points de calcul.

Les principaux résultats montrent un très bon accord avec les statistiques observées expérimentalement et une croissance rapide de la couche de mélange dans les 3 premières longueurs (basées sur l'épaisseur moyenne des couches limites en sortie de plaque). Cette première zone montre des tensions de Reynolds importantes et est principalement causée par la turbulence des couches limites amonts. D'autre part, une région intermédiaire a été mise en évidence jusqu'à une trentaine
<figure><img src="image_28.png" /><figcaption>FIGURE 6.1 – Mode DMD principal de la pression dans le plan xz</figcaption></figure>

de longueurs dans laquelle les structures cohérentes de pression et de température semblent inclinées d'un angle de 70 degrés par rapport à la direction de l'écoulement. Enfin, la zone de mélange observée après 35 longueurs se développe progressivement vers un état de couche de mélange développée sans pour autant atteindre d'état auto-similaire. Une décomposition dynamique de l'écoulement dans les 3 régions caractéristiques (figure 6.1) a confirmé ces observations et a également permis d'identifier les structures cohérentes qui pourraient expliquer la formation des rouleaux de vorticité classiquement observés dans le sillage lointain des couches de mélange [3, 7].

## 6.2 Modélisation des fluctuations de concentration

Ce travail s'inscrit directement dans la suite des activités liées aux écoulements d'arrières-corps. En effet, même si les mécanismes instationnaires se produisant dans la zone décollée ont été clairement identifiés en l'absence de jets propulsifs, l'impact des jets exothermiques sur la dynamique d'arrières-corps n'est pas encore identifié. D'autre part, dans le contexte du développement d'Ariane6 qui a vu sont arrière-corps largement caréné mais reprenant certains éléments de l'architecture d'Ariane5, les effets thermiques générés par l'environnement chaud des jets propulsifs deviennent importants. Ainsi l'étude de ce type d'écoulement dans un cadre industriel est un point important et les outils existants nécessitaient d'être adaptés aux besoins. Ainsi, nous avons lancé une thèse en partenariat avec ArianeGroup pour développer les outils numériques et les rendre compatibles avec l'étude des écoulements d'arrières-corps en présence de jets propulsifs réactifs. Ces travaux sont structurés en deux grandes parties :

- La modélisation acceptable des flux thermiques et massiques turbulents.

- L'implémentation d'un modèle d'interaction chimie/turbulence polyvalent et pas trop coûteux en ressources informatiques.

Le premier point est lié à la modélisation turbulente des fluctuations de concentration et de température des espèces composant le mélange. Les approches standards imposent des valeurs de Prandtl et de Schmidt turbulent pour calculer les coefficients de diffusion turbulente $\alpha_t$ et $D_t$ :

$$ \alpha_t = \nu_t / Pr_t \quad D_t = \nu_t / Sc_t \quad (6.1) $$
Dans le cas de la modélisation RANS, de nombreux modèles existent et pro-
posent de calculer directement les diffusivités turbulentes $\alpha_t$ et $D_t$ en rajoutant
des équations de transport [36] :

$$
\alpha_t = C_{\alpha} f_{\alpha} k \tau_{m,\alpha} \quad D_t = C_D f_D k \tau_{m,D} \qquad (6.2)
$$

Dans l'optique d'une modélisation RANS/LES, il a été proposé dans le cadre
de la thèse de Lorris CHARRIER d'utiliser des équations de transport basées sur
les variances et les taux de dissipation turbulente des espèces. L'hybridation est
alors classiquement basée sur la dissipation des fluctuations de température ou de
concentration :

$$
\epsilon_{\alpha} = k_{\alpha} \frac{\epsilon^{1/3}}{(l_{\alpha})^{2/3}} \qquad (6.3)
$$

Avec, $l_\alpha = \min(l_\alpha^{\text{LES}}, l_\alpha^{\text{RANS}})$ et la variance $k_\alpha = \overline{T'^2}$. La modélisation des variances proposée est basée sur les approches de Brinckmann [2] pour la température et celle de Xiao [35] pour la concentration. Dans ce cas, les nombres de Schmidt et de Prandtl turbulent peuvent être différents en fonction de la zone de l'écoulement.

Le second point concerne la modélisation turbulente des taux de réaction $\overline{\dot{\omega}_k}$ dans le cadre des mélanges réactifs. Les approches existantes en RANS utilisent des distributions de probabilités (PDF) mais restent très couteuses pour atteindre une bonne précision [1, 15, 16]. Des approches simplifiées ont été implémentées dans le cadre de la thèse en supposant la décorrélation des fluctuations de concentration et de température :

$$
\overline{\dot{\omega}_k} = W_{ki} \overline{c_i^z k(T)} \quad (6.4)
$$

où $\overline{k(T)}$ est développé en série de Taylor. Toutefois cette hypothèse s'est avérée
globalement fausse en raison des aspects fortement non-linéaires intervenant dans
les types d'écoulements considérés. Une amélioration proposée consiste à assumer
une échelle de mélange plus grossière regroupant des espèces initialement pré-
mélangées sous la forme de « blends »[6]. Là encore, les simulations n'ont pas
montré de résultats probants permettant de valider clairement ces approches.

Une implémentation RANS/LES des aspects présentés plus haut a été implé-
menté dans le code FLUSEPA d'ArianeGroup et deux simulations ont été réalisées
sur un jet coaxial supersonique et une couche de mélange réactive.

6.2.1 Validation sur un jet coaxial supersonique

Le jet de Cutler [8] consiste à injecter un jet d'hélium à Mach = 1.8 dans un jet d'air au même Mach. La vitesse du son étant plus faible dans l'air que dans l'Hélium, le Mach convectif est de 0.7. Des simulations ont été réalisées en reproduisant l'intégralité du dispositif expérimental et en utilisant les versions RANS/LES des modèles proposés (figure 6.2).
<figure><img src="image_29.png" /><figcaption>FIGURE 6.2 – Fraction massique de l'hélium à différentes positions du jet. Les symboles représentent les données expérimentales.</figcaption></figure>

Les résultats montrent un léger effet de la modélisation des fluctuations de concentration (Schmidt variable) en particulier dans le cadre des simulations RANS. D'autre part, l'hybridation de ces approches reste fortement sensible à la résolution et les calculs hybrides réalisés sur maillage fin ne montrent pas de dépendance au nombre de Schmidt.

### 6.2.2 Simulation d'une couche de mélange réactive

Afin de réaliser des simulations mettant en jeux des réactions chimiques, le calcul de la couche de mélange réactive de Miller [22] a été réalisé pour différentes résolutions. Cette configuration confronte de l'air chaud supersonique avec un mélange subsonique d'azote et d'hydrogène froid. Les deux gaz se rencontrent dans une veine d'essai de section réduite en aval d'une plaque séparatrice.

Les résultats montrent une bonne reproduction de l'écoulement moyen ainsi qu'une bonne description des réseaux de choc en aval de la plaque (figure 6.3). Les simulations ont également permis d'extraire certaines informations complémentaires comme l'apparition transitoire de tourbillons bidimensionnels à la pointe du séparateur et des oscillations périodiques entre les régimes supersoniques et subsoniques. Ces oscillations disparaissent après un long temps d'intégration rendant les simulations couteuses pour atteindre un état comparable à la configuration expérimentale.

La restitution des quantités fluctuantes dans le cadre des calculs hybrides souffre des hypothèses faites sur les termes de taux de réaction et ne donnent pas de résultats suffisamment fiables pour valider les approches proposées. De plus certains détails de la configuration expérimentale rendent les comparaisons diffi-
<figure><img src="image_30.png" /><figcaption>FIGURE 6.3 – Comparaison des Schlieren expérimentaux (Haut) et numériques de la couche de mélange de Miller obtenue pour différents maillages des simulations RANS/LES.(13, 7 et 4 million de points.)</figcaption></figure>

ciles à exploiter. Ainsi, de nouvelles investigations seront nécessaires pour améliorer
la modélisation des taux de réaction dans le cadre des approches RANS/LES en
incluant des corrélations concentration/température.
Références

[1] Roland BORGHI. "Chemical Reactions Calculations in Turbulent Flows : Application to a Co-Containing Turbojet Plume". In : *International Union of Theoretical and Applied Mechanics and International Union of Geodesy and Geophysics - Turbulent Diffusion in Environmental Pollution, Proceedings of a Symposium held at Charlottesville*. Elsevier, 1975, p. 349-365. DOI : 10.1016/s0065-2687(08)60592-0 (cf. p. 73).

[2] Kevin W. BRINCKMAN, William H. CALHOON et Sanford M. DASH. "Sca-lar Fluctuation Modeling for High-Speed Aeropropulsive Flows". In : *AIAA Journal* 45.5 (mai 2007), p. 1036-1046. DOI : 10.2514/1.21075 (cf. p. 73).

[3] G.L. BROWN et A. ROSHKO. "On density effects and large structure in tur-bulent mixing layers." In : *Journ.Fluid.Mech* 64 (1974), p. 775-816 (cf. p. 72).

[4] S. CAMARRI et al. "Large-eddy simulation of a bluff-body flow on unstructu-red grids". In : *International Journal for Numerical Methods in Fluids* 40.11 (2002), p. 1431-1460. DOI : 10.1002/fld.425 (cf. p. 65).

[5] R. CAMUSSI et al. "Wall pressure fluctuations induced by transonic boundary layers on a launcher model." In : *Aerospace Science and Tech.* 11 (2007), p. 349-359 (cf. p. 65).

[6] L. CHARRIER. "Développement de modèles Hybrides RANS/LES avancés pour la simulation numérique instationnaire de jets propulsifs supersoniques et réactifs". Thèse de doct. Ecole Nationale Supérieur d'Arts et Métiers, 2018 (cf. p. 73).

[7] N. T. CLEMENS et M. G. MUNGAL. "Large-scale structure and entrainment in the supersonic mixing layer". In : *Journal of Fluid Mechanics* 284 (fév. 1995), p. 171-216. DOI : 10.1017/s0022112095000310 (cf. p. 72).

[8] A. D. CUTLER et al. "Supersonic Coaxial Jet Experiment for Computational Fluid Dynamics Code Validation". In : *AIAA Journal* 44.3 (mar. 2006), p. 585-592. DOI : 10.2514/1.5781 (cf. p. 73).

[9] S. DAVID et S. RADULOVIC. "Prediction of Buffet Loads on the Ariane 5 Af-terbody". In : *6th Symposium on Launcher Technologies*. Munich, Germany, août 2005 (cf. p. 68).

[10] S. DECK et A.T. NGUYEN. "Unsteady Side-Loads in a Thrust-Optimized Contour Nozzle at Hysteresis Regime." In : *AIAA Journ.* 42.9 (sept. 2004), p. 1878-1888 (cf. p. 66, 67).
[11] S. DECK et P. THORIGNY. "Unsteadiness of an axisymmetric separating-reattaching flow : Numerical investigation." In : *Phys. Fluids* 19.065103 (2007), p. 1-20 (cf. p. 65).

[12] S. DECK et P. THORIGNY. "Unsteadiness of an axisymmetric separating-reattaching flow." In : *Phys. Fluids*. 19 (2007), p. 065103 (cf. p. 65).

[13] D. DEPRÉS, P. REIJASSE et J.P. DUSSAUGE. "Analysis of Unsteadiness in Afterbody Transonic Flows". In : *AIAA Journ.* 42.12 (déc. 2004), p. 2341-2550 (cf. p. 65).

[14] G.E. DUMNOV. "Unsteady Side-Loads Acting on the Nozzle with Developed Separation Zone." In : *AIAA-paper* 96-3220 (1996) (cf. p. 65).

[15] R.L. GAFFNEY et al. "Modeling temperature and species fluctuations in turbulent, reacting flow". In : *Computing Systems in Engineering* 5.2 (avr. 1994), p. 117-133. DOI : 10.1016/0956-0521(94)90044-2 (cf. p. 73).

[16] Peter GERLINGER et Manfred AIGNER. "Assumed PDF Modeling with Detailed Chemistry". In : *High Performance Computing in Science and Engineering '01*. Springer Berlin Heidelberg, 2002, p. 317-328. DOI : 10.1007/978-3-642-56034-7_31 (cf. p. 73).

[17] E. GEURTS. "Steady and unsteady pressure measurements on the rear section of various configurations of the Ariane 5 launch vehicle". In : *6th Symposium on Launcher Technologies*. Munich, Germany, août 2005 (cf. p. 68).

[18] Steven G. GOEBEL et J. Craig DUTTON. "Experimental study of compressible turbulent mixing layers". In : *AIAA Journal* 29.4 (avr. 1991), p. 538-546. DOI : 10.2514/3.10617 (cf. p. 71).

[19] K. HANNEMANN et al. *Launch vehicle base Buffeting : Recent experimental and numerical investigations*. Rapp. tech. Unsteady Subscale Force Measurements within a Launch Vehicle Base Buffeting Environment, 2011 (cf. p. 68).

[20] Christopher a KENNEDY et Andrea GRUBER. "Reduced aliasing formulations of the convective terms within the Navier-Stokes equations for a compressible fluid". In : *J. Comput. Phys.* 227.3 (2008), p. 1676-1700. ISSN : 00219991. DOI : 10.1016/j.jcp.2007.09.020 (cf. p. 71).

[21] J.L. LUMLEY. "The structure of inhomogeneous turbulent flows". In : *Atm. Turb. and Radio wave Prop.*. Sous la dir. d'YAGLOM et TATARSKY. 1967, p. 166-178 (cf. p. 68).

[22] M. F. MILLER, C. T. BOWMAN et M. G. MUNGAL. "An experimental investigation of the effects of compressibility on a turbulent reacting mixing layer". In : *Journal of Fluid Mechanics* 356 (fév. 1998), p. 25-64. DOI : 10.1017/s002211209700791x (cf. p. 74).

[23] A.T. NGUYEN. "Décollement instationnaire et charges latérales dans les tuyères propulsives". Thèse de doct. Université de Poitiers, 2003 (cf. p. 67).
[24] Romain PAIN, Pierre-Élie WEISS et Sebastien DECK. "Zonal Detached Eddy Simulation of the Flow Around a Simplified Launcher Afterbody". In: *AIAA Journal* 52 (sept. 2014), p. 1967-1979. DOI: 10.2514/1.J052743 (cf. p. 67).

[25] Romain PAIN et al. "Large scale dynamics of a high Reynolds number axisymmetric separating/reattaching flow". In: *Physics of Fluids* 31 (déc. 2019), p. 125119. DOI: 10.1063/1.5121587 (cf. p. 67).

[26] Sergio PIROZZOLI. "Numerical Methods for High-Speed Flows". In: *Annu. Rev. Fluid Mech.* 43.1 (jan. 2011), p. 163-194. ISSN: 0066-4189. DOI: 10.1146/annurev-fluid-122109-160718 (cf. p. 71).

[27] PETER J. SCHMID. "Dynamic mode decomposition of numerical and experimental data". In: *Journal of Fluid Mechanics* 656 (juil. 2010), p. 5-28. DOI: 10.1017/s0022112010001217 (cf. p. 67).

[28] F.F.J. SCHRIJER et al. "Experimental investigation of base flow buffeting on the Ariane 5 launcher using high speed PIV." In: *7th European Symposium on Aerothermodynamics for Space Vehicles*. Bruges, Belgium, 2011 (cf. p. 69).

[29] L. SIROVICH. "Turbulence and the dynamics of coherent structures. Part I : Coherent structures". In: *Q. Appl. Math* XLV (1987), p. 561-571 (cf. p. 68).

[30] F. GARÇON et J.P. DREVET. "Unsteady load measurements on the main engine nozzle of the Ariane 5 launch vehicle". In: *Third European Symposium on Aerothermodynamics for Space Vehicles*. ESA. ESTEC, Noordwijk, Yche Netherlands, déc. 1998, p. 623 (cf. p. 65).

[31] Paul VAN GENT et al. "Comparative assessment of pressure field reconstructions from particle image velocimetry measurements and Lagrangian particle tracking". In: *Experiments in fluids* 58:33 (mar. 2017). DOI: 10.1007/s00348-017-2324-z (cf. p. 67, 70).

[32] P.E WEISS et S. DECK. "Zonal Detached Eddy Simulation of the flow dynamics on an Ariane 5-type afterbody". In: *4th Eur. Conf. for Aerospace Sciences*. St-Petersbourg, Russia, avr. 2011 (cf. p. 68).

[33] P.E. WEISS et al. "On the dynamics of axisymmetric turbulent separating/reattaching flows." In: *Physics of Fluids* 21 (2009), p. 075103 (cf. p. 65).

[34] Pierre-Élie WEISS et Sébastien DECK. "Control of the antisymmetric mode ($m=1$) for high Reynolds axisymmetric turbulent separating/reattaching flows". In: *Physics of Fluids* 23.9 (sept. 2011), p. 095102. DOI: 10.1063/1.3614481 (cf. p. 66).

[35] X. XIAO, H. A. HASSAN et R. A. BAURLE. "Modeling Scramjet Flows with Variable Turbulent Prandtl and Schmidt Numbers". In: *AIAA Journal* 45.6 (juin 2007), p. 1415-1423. DOI: 10.2514/1.26382 (cf. p. 73).

[36] Dennis A. YODER. "Comparison of Turbulent Thermal Diffusivity and Scalar Variance Models". In: *54th AIAA Aerospace Sciences Meeting*. American Institute of Aeronautics et Astronautics, jan. 2016. DOI: 10.2514/6.2016-1561 (cf. p. 73).
Conclusion et Perspectives
<figure><img src="image_31.png" /></figure>
Conclusion

Nous avons donc fait dans ce mémoire une synthèse des travaux entrepris de-
puis une dizaine d'années. Nous avons d'abord rappelé le contexte théorique des
méthodes LBM ainsi que les avancées récentes sur le sujet permettant de mieux
comprendre les travaux réalisés par la suite. En effet, les approches développées au
cours de mes travaux basées sur des corrections locales de la dissipation s'appuient
principalement sur les méthodes locales ou entropiques visant à corriger le schéma
ou à modifier localement la relaxation. La particularité des approches développées
réside principalement dans la simplification des critères sélectifs entraînant une
modification de la dissipation visqueuse sans toucher à la dissipation acoustique.
Ces types de corrections couplées aux méthodes de régularisation paraissent pro-
metteuses pour envisager des couplages avec les approches Navier-Stokes ou dans
une moindre mesure avec les approches Vortex remaillées qui présentent des carac-
téristiques lagrangiennes intéressantes pour les faibles résolutions. Nous avons vu
que le point clé dans la comparaison et le couplage avec d'autres méthodes résidait
dans la traduction des quantités mésoscopiques en quantité macroscopiques via le
calcul des moments et la traduction inverse via le calcul de l'équilibre.

Dans la deuxième partie, nous nous sommes focalisés sur la modélisation des
écoulements complexes en LBM en traitant d'abord des interfaces. L'application
des filtrages adaptatifs aux modèles LBM des écoulements diphasiques, proches
des méthodes level-set, ont permis une légère augmentation des rapports de den-
sité et des nombres de Reynolds dans les simulations. Enfin nous avons présenté
un nouveau modèle d'interaction entre espèces dans les simulations LBM via l'ad-
jonction d'un terme de forçage. L'étude de l'instabilité de Saffman-Taylor avec ce
modèle a permis de mettre en évidence la phase de transition de l'instabilité en
présence de trois espèces dans un cas initialement défavorable.

Dans la dernière partie nous avons vu comment la projection des données de
pression sur des bases appropriées permettaient une meilleure compréhension des
phénomènes de recollement instationnaires sur les arrières-corps de lanceurs. Puis,
après avoir mis en évidence une zone de transition entre le développement de
couches limites et une couche de mélange, nous avons tenté d'explorer quelques
pistes pour l'hybridation des modèles de turbulence en présence de réactions chi-
miques à plusieurs espèces. Ce dernier point était destiné à mieux comprendre
l'impact des jets propulsifs sur la dynamique et la thermique des arrières-corps de
lanceurs.
# Perspectives

Pour finir, nous allons ici dresser un rapide aperçu des activités de recherches que je souhaite mener dans les prochaines années. Ces activités s'articulent autour de 3 axes principaux :

1. Utilisation de la LBM pour la simulation acoustique.

2. Les écoulements multi-espèces en LBM.

3. La dynamique de foule.

Ces trois axes sont présentés dans les sections suivantes.

## Utilisation de la LBM pour la simulation acoustique

### Simulation et Etudes des sources acoustiques non linéaires dans les instruments à vents

L'utilisation de l'acoustique non linéaire dans les guides cylindriques est un domaine largement étudié dans les applications industrielles. En particulier, le contrôle non-destructif (CND) pour la détection de fissures, l'échographie médicale pour l'identification des tumeurs et le filtrage des nano-particules ont fait l'objets de récentes études [13] [11] [16]. L'enjeu majeur dans la maîtrise des outils faisant intervenir la propagation d'ultrason se situe dans la description précise du champ acoustique rayonné ainsi que dans la description précise des matériaux [6]. En effet, les champs acoustiques non linéaires, de part leur forte amplitude et leur propagation déformée, reste difficile à étudier de façon précise aussi bien d'un point de vue expérimental que numérique.

La propagation acoustique non linéaire intervient également dans le design des instruments à vent, en particulier dans la maîtrise et l'optimisation de leur jouabilité [4]. L'étude numérique de ce type d'écoulement est encore très peu étudié par les approches directes (temporelles) qui est pourtant nécessaire en raison de la mise en défaut des approches modales (fréquentielles) pour les régimes non linéaires.

Comme nous l'avons vu dans ce mémoire, les caractéristiques numériques de la LBM en font un très bon candidat pour la simulation des fronts d'ondes acoustiques dans les guides cylindriques [1, 9] et en particulier pour la propagation des hautes fréquences et des régimes à forte amplitude. Ainsi un de mes prochains axes de recherche sera consacré à l'utilisation des outils développés jusqu'à présent pour
la modélisation des interactions non linéaires intervenant dans la génération des ondes acoustiques et leur propagation non linéaire dans les guides cylindriques.

Cet axe de recherche fera l'objet d'une collaboration avec Jean-Baptiste DOC [3] du laboratoire LMSSC qui effectue ses recherches principalement sur l'acoustique des instruments à vent et en particulier sur les embouchures de trompettes. L'objectif sera de lancer une thèse sur le sujet.

## Couplage de la LBM avec d'autres schémas pour la propagation acoustique

Les activités initiées et présentées dans la section 2.3 seront poursuivies dans le cadre de la thèse d'Alexandre SUSS. Les premiers résultats encourageants devront permettre de réaliser dans les deux prochaines années des simulations de taille industrielle dans lesquelles la résolution fine des couches limites assurée par la partie Navier-Stokes permettra d'alimenter les parties détachées du calcul responsable de la propagation d'ondes acoustiques assurées par le solveur LBM. Ainsi, les écoulements de sillages étudiés dans le chapitre 5 pourraient être simulés dans le cadre d'une formulation hybride faisant intervenir le solveur NS dans les couches limites attachées et un solveur LBM dans la zone de sillage turbulent.

## Etude des écoulements multi-espèces en LBM

Les écoulements multi-espèces étudiés au chapitre 4 ont montré leur intérêt pour l'étude des écoulements faisant intervenir de la diffusion entre plusieurs espèces et une convection globale du mélange. Cependant, afin de tirer profit des avantages d'une approche indépendante par espèce, une étude détaillée de la dynamique du mélange dans la zone de transition devra être menée. Ainsi, à la manière des études de la section 6.1 réalisées dans la zone de transition d'une couche de mélange, il sera intéressant de comparer les résultats obtenus par le modèle proposé à ceux de simulations classiques réalisées par des codes multi-espèces utilisant une approche parfaitement mélangée pour décrire la quantité de mouvement. En effet, le transport de la quantité de mouvement puis le mélange de deux espèces différentes ayant une viscosité différente, doit se traduire par une zone de transition dans laquelle le transfert des quantités de mouvement de chaque espèce ne peut pas être décrit par un seul coefficient de viscosité. Les résultats de ce type de simulation devront permettre de valider ou d'invalider un des deux types d'approche dans cette zone particulière de l'écoulement.

## Etude LBM de la dynamique de foule

Depuis les travaux de Helbin et Molnár [7], la simulation et l'étude de la dynamique de foule a fait l'objet de nombreuses études utilisant des approches différentes. Si le modèle de « Social Force » développé par ces auteurs est aujourd'hui
largement employé, d'autres approches existent. En particulier, les approches utilisant des automates cellulaires ont récemment été revisitées [15, 5] et présentent des similarités avec les modèles LBM,[12], que certains utilisent comme algorithme de contrôle pour la dynamique de foule [14]. Dès lors, il devient intéressant de pouvoir relier les deux approches en identifiant les éléments importants dans la dynamique de foule et développer un modèle de collision compatible avec cette dynamique [8].

<figure><img src="image_32.png" /><figcaption>FIGURE 6.4 – Comparaison du modèle de Yamamoto et du modèle réalisé dans le cadre du stage de C.Pruvot avec un automate cellulaire basé sur les niveaux d'énergies. Temps d'évacuation en fonction du nombre de personnes.</figcaption></figure>

Une première approche a été réalisée lors du stage de Colin Pruvot réalisé en 2020. L'objectif était de reconstruire un automate cellulaire simple capable de reproduire les résultats obtenus par Yamamoto en 2007 [15] en se basant sur les interactions à courtes distances et les contraintes de longues portées (figure 6.4). La prochaine étape de cet axe d'étude sera d'utiliser cet automate cellulaire pour générer des données de simulations qui pourront servir de données d'entrainement à des algorithmes d'identification de système dynamique [10] afin de pouvoir construire à terme un modèle de collision fiable pour la simulation des mouvements de foule. A plus long terme, les modèles à plusieurs composants développés dans le chapitre 4 pourront être intégrés dans les modèles de simulation des dynamiques de foules à multiples contraintes [2].
RÉFÉRENCES
Références

[1] James M. BUICK et al. "Investigation of non-linear acoustic losses at the open end of a tube". In : *The Journal of the Acoustical Society of America* 129.3 (mar. 2011), p. 1261-1272. DOI : 10.1121/1.3543987 (cf. p. 85).

[2] Xu CHEN et al. "A multiagent-based model for pedestrian simulation in subway stations". In : *Simulation Modelling Practice and Theory* 71 (fév. 2017), p. 134-148. DOI : 10.1016/j.simpat.2016.12.001 (cf. p. 87).

[3] Tom COLINOT et al. "Multiple two-step oscillation regimes produced by the alto saxophone". In : *The Journal of the Acoustical Society of America* 147.4 (avr. 2020), p. 2406-2413. DOI : 10.1121/10.0001109 (cf. p. 86).

[4] Jean-Baptiste DOC et Christophe VERGEZ. "Oscillation regimes produced by an alto saxophone : Influence of the control parameters and the bore inharmonicity". In : *The Journal of the Acoustical Society of America* 137.4 (avr. 2015), p. 1756-1765. DOI : 10.1121/1.4916197 (cf. p. 85).

[5] Claudio FELICIANI et Katsuhiro NISHINARI. "An improved Cellular Automata model to simulate the behavior of high density crowd and validation by experimental data". In : *Physica A : Statistical Mechanics and its Applications* 451 (2016), p. 135-148. ISSN : 0378-4371. DOI : https://doi.org/10.1016/j.physa.2016.01.057 (cf. p. 87).

[6] Stephan GOLLOB et Georg Karl KOCUR. "Analysis of the wave propagation paths in numerical reinforced concrete models". In : *Journal of Sound and Vibration* 494 (mar. 2021), p. 115861. DOI : 10.1016/j.jsv.2020.115861 (cf. p. 85).

[7] D. HELBING et P. MOLNÁR. "Social force Model for Pedestrian Dynamics." In : *Physical Review E* 51(5) (mai 1995), p. 4282-4286 (cf. p. 86).

[8] L. HENDERSON. "The Statistics of Crowd Fluids." In : *Nature* 229 (1971), 381-383 (cf. p. 87).

[9] Chenzhen JI et Dan ZHAO. "Lattice Boltzmann investigation of acoustic damping mechanism and performance of an in-duct circular orifice". In : *The Journal of the Acoustical Society of America* 135.6 (juin 2014), p. 3243-3251. DOI : 10.1121/1.4876376 (cf. p. 85).

[10] Jean-Christophe LOISEAU, Bernd R. NOACK et Steven L. BRUNTON. "Sparse reduced-order modelling : sensor-based dynamics to full-state estimation". In : *Journal of Fluid Mechanics* 844 (avr. 2018), p. 459-490. DOI : 10.1017/jfm.2018.147 (cf. p. 87).
[11] Claudio PECORARI et Milan POZNIĆ. "Nonlinear acoustic scattering by a partially closed surface-breaking crack". In : *The Journal of the Acoustical Society of America* 117.2 (fév. 2005), p. 592-600. DOI : 10.1121/1.1850052 (cf. p. 85).

[12] Y.-H. QIAN, D. D'HUMIÈRES et P. LALLEMAND. "Lattice BGK models for Navier-Stokes equation." In : *Europhys.Lett.* 17 (1992), p. 479-484 (cf. p. 87).

[13] Feng SHAN et al. "Acoustic field of an ultrasonic cavity resonator with two open ends : Experimental measurements and lattice Boltzmann method modeling". In : *Journal of Applied Physics* 121.12 (mar. 2017), p. 124502. DOI : 10.1063/1.4978013 (cf. p. 85).

[14] Yiran XUE et al. "Abnormal Prediction of Dense Crowd Videos by a Purpo-
seâDriven Lattice Boltzmann Model". In : *International Journal of Applied
Mathematics and Computer Science* 27.1 (2017), p. 181-194 (cf. p. 87).

[15] Kazuhiro YAMAMOTO, Satoshi KOKUBO et Katsuhiro NISHINARI. "Simulation for pedestrian dynamics by real-coded cellular automata (RCA)". In : *Physica A : Statistical Mechanics and its Applications* 379.2 (2007), p. 654-660. ISSN : 0378-4371. DOI : https://doi.org/10.1016/j.physa.2007.02.040 (cf. p. 87).

[16] W. T. YUEN et al. "The Use of Nonlinear Acoustics as an Energy-Efficient Technique for Aerosol Removal". In : *Aerosol Science and Technology* 48.9 (août 2014), p. 907-915. DOI : 10.1080/02786826.2014.938800 (cf. p. 85).
Adresse

Laboratoire DynFLuid,
151 boulevard de l'Hopital
75013 Paris, France

Téléphone

+33 6-25-50-74-15

Email

simon.marie@lecnam.net

Web

hpp.education

Date de Naissance

02/10/1981

Nationalité

Française

## Mots Clés

Boltzmann sur Réseau, mécanique des fluides numérique, schémas numérique, simulation, aérodynamique, aéroacoustique, dynamique de mélange, écoulement multi-espèces.

## Situation actuelle

### Maître de Conférence

CONSERVATOIRE NATIONAL DES ARTS ET MÉTIERS

Laboratoire DynFLuid

Cours, TD et TP en *Simulation numérique (M1)*, *LBM (M2)*, mécanique des fluides (L2), Aérody-namique (L2), Mécanique du Vol (L3), Mécanique spatiale(M2), Aéroacoustique (M2).

## Etudes

**2008** DOCTORAT "Etude de la méthode de Boltzmann sur Réseau pour les simulations en aéroacoustique."
INSTITUT JEAN LE ROND D'ALEMBERT - UNIVERSITÉ PIERRE ET MARIE CURIE, Paris

**2004** DEA - ACTA Aérodynamique et Aéroacoustique,
UNIVERSITÉ DE POITIERS

**2001-2004** ENSIP option Acoustique
ECOLE NATIONALE SUPÉRIEURE D'INGÉNIEURS DE POITIERS

**1999-2001** Classes préparatoires MPSI/MP au lycée LaFayette à Clermont-Ferrand

## Expérience Professionnelle

**2009-2011** Post-Doc - CENTRE NATIONAL D'ETUDES SPATIALES, "Etudes des écoulements d'arrière-corps."

**2008-2009** Post-Doc - ONERA - DAAP, "Etude numérique d'une configuration d'arrière corp à partir de simulation ZDES"

**2005-2008** Thèse de Doctorat chez RENAULT.
Groupe Acoustique de la Direction de la Recherche au Technocentre.

**2004** Stage de fin d'étude au CNES sur l'étude des écoulements d'arrière-corps.
# Projets de recherche

**CARNOT** 2021-2023 *Caractérisation numérique et expérimentale du rayonnement de sources acoustiques guidée à forte amplitudes.*, collaboration avec le LMSSC.

**GENCI** 2018-2020 *Simulation des écoulements multi-espèces par la méthode de Boltzmann sur réseau*, 2.1 Millions d'heures CPU sur Occigen (Grant 2018-A0052A10636, 2019-A0072A10636, 2020-A0092A10636)

**PRACE** 2014-2015 *Direct Numerical Simulation of a Spatially Developing Mixing Layer* en collaboration avec la Sapienza, Rome: 5 Millions d'heures CPU sur JUQUEEN (grant-2013081603)

**TRP-ESA** 2010-2013 "Experimental study of unsteady Ariane-5 base flow" 2010-2013 (TU-DELFT, ONERA, CNES, ESA, ArianeGroup)

**ATAC** 2008-2011 "Aérodynamique des Tuyères et des Arrières-Corps" 2009-2011 (CNES, ONERA)

**MIMOSA** 2006-2009 "Méthodes Innovantes pour la Modélisation de Sources Aéroacoustiques" 2006-2009 (Rnault, SNCF, ALSTOM, ECL, Paris-6)

# Encadrement

## Thèses

**2020-2023** A. SUSS "Construction et analyse d'un solver hybride Lattice Boltzmann / Navier-Stokes pour les simulations en aerodynamique" Thèse ONERA

**2016-2019** L. Vienne "Etude des écoulements multi-espèces par la méthode de Boltzmann sur Réseau" Thèse CNAM

**2014-2018** L. Charrier "Etude des écoulements multi-espèces réactifs par approches RANS/LES pour la simulation des écoulements d'arrière-corps.", CIFRE ARIANEGROUP

## Stage de Master

**2020** C. Pruvot "Etude de la dynamique de foule par un automate cellulaire".

**2020** A. Suss "Etude du couplage LBM / Navier stokes pour les simulations en aérodynamique".

**2016** B. Aladdin "Simulation Boltzmann sur Réseau de l'instabilité de Rayleigh-Taylor"

**2014** X. Dongmin "Etude d'une couche de mélange réactive"

# Activités administratives

**2021-2023** Responsable de la 2ème année de la formation IAS aéronautique et espace.

**2020-2022** Représentant CNAM du laboratoire DynFluid.

**2012-2021** Membre du conseil de laboratoire de DynFluid depuis 2012.

**2016-2020** Membre élu du conseil de l'EPN4.

# Animation et valorisation de la recherche

**2020-2021** Développement des simulateurs de vol et de lanceur spatial pour la formation IAS-Aéronautique et espace et la formation HTT.

**2016-2018** Création d'un groupe de travail pour l'utilisation d'un serveur JupyterHUB de notebook interactifs pour l'enseignement et la recherche.
**2015-2019** Co-responsable du groupe de Travail Numérique du laboratoire DynFLuid et développement du site intranet du laboratoire.

**2010-2021** Revues d'articles pour des journaux internationaux (Journal of computational Physics, Computers & Fluids, Journal of the American Society of Acoustics...)

# Publications

## Journaux à comité de lectures

[1] S. Marié. "On the use of adaptive relaxation times in Lattice Boltzmann Methods". In: *Progress in Computational Fluid Dynamics, Accepted for publication* (2021)

[2] L. Vienne and S. Marié. "A Lattice Boltzmann study of miscible viscous fingering for binary and ternary mixtures". In: *Physical Review F* 6 (5 2021), p. 053904

[3] Chloé Mimeau, Simon Marié, and Iraj Mortazavi. "A comparison of semi-Lagrangian Vortex method and Lattice Boltzmann method for incompressible flows". In: *Computers & Fluids* 224 (2021), p. 104946

[4] L. Vienne, S. Marié, and F. Grasso. "Lattice Boltzmann model for miscible gases: a forcing term approach". In: *Physical Review E* 100 (2) (2019), p. 023309

[5] Lorris Charrier, Mathieu Jubera, Grégoire Pont, Simon Marié, Pierre Brenner, and Francesco Grasso. "Simulations of reactive supersonic/subsonic flow interactions for space launcher applications". In: *International Journal of Numerical Methods for Heat & Fluid Flow* 31.10 (Aug. 2021)

[6] S. Marié and X. Gloerfelt. "Adaptive filtering for the lattice Boltzmann method". In: *Journal of Computational Physics* 333C (2017), pp. 212-226

[7] S. Pirozzoli, M. Bernardini, S. Marié, and F. Grasso. "Early evolution of the compressible mixing layer issued from two turbulent streams". In: *Journal of Fluid Mechanics* 777 (2015), pp. 196-218

[8] S. Marié, P. Druault, H. Lambaré, and F. Schrijer. "Experimental analysis of the pressure-velocity correlations of external unsteady flow over rocket launchers". In: *Aerospace Science and Technology* 30-1 (2013), pp. 83-93

[9] S. Marié, S. Deck, and P-E. Weiss. "From pressure fluctuations to dynamic loads on axisymmetric step flows with minimal number of kulites". In: *Computers and Fluids* 39-5 (2010), pp. 747-755

[10] D. Ricot, S. Marié, and P. Sagaut. "Lattice Boltzmann Method with selective viscosity filter". In: *Journal of Computational Physics* 228-12 (2009), pp. 4478-4490

[11] S. Marié, D. Ricot, and P. Sagaut. "Comparison between lattice Boltzmann method and Navier-Stokes high order schemes for computational aeroacoustic". In: *Journal of Computational Physics* 228-4 (2009), pp. 1056-1070

[12] Y. Dong, P. Sagaut, and S. Marié. "Inertial consistent subgrid model for Large Eddy Simulation based on Lattice Boltzmann method". In: *Physics of Fluids* 20-3 (2008), p. 035104

## Chapitre de Livre

[1] L. Charrier, G. Pont, S. Marié, P. Brenner, and F. Grasso. "Hybrid RANS/LES simulation of a supersonic coaxial He/air jet experiment at various turbulent Lewis numbers". In: *Notes on Numerical Fluid Mechanics and Multidisciplinary Design* 137 (2018), pp. 337-346

## Conférences avec Actes

[1] Lucien Vienne., Simon Marié., and Francesco Grasso. "Simulation of Viscous Fingering Instability by the Lattice Boltzmann Method". In: vol. AIAA-paper 2019-3432. AIAA AVIATION CFD. Dallas, USA, 2019
[2] L. Charrier, G. Pont, S. Marié, P. Brenner, and F. Grasso. "Simulation of a supersonic coaxial He/Air jet experiment with a hybrid RANS/LES variable turbulent Schmidt number model". In: vol. AIAA-paper 2017-4280. 23rd AIAA Computational Fluid Dynamics Conference. Denver, CO, USA, 2017

[3] S. Marié and X. Gloerfelt. "Dynamical selective filtering for the Lattice Boltzmann Method". In: vol. AIAA-paper 2015-2514. 22nd AIAA Computational Fluid Dynamics Conference. Dallas (USA), 2015

[4] S. Marié, D. Ricot, and P. Sagaut. "Accuracy of the Lattice Boltzmann method for aeroacoustic simulations". In: vol. AIAA-paper 2007-3515. 13th AIAA Aeroacoustic conference. Rome (Italy), 2007

## Conférences sans Actes

[1] Lucien Vienne, Simon Marié, and Francesco Grasso. "Lattice Boltzmann model for miscible mixtures: application to the viscous fingering instability". In: ICMMES. Edimbourgh, UK, 2019

[2] Simon Marié, and Xavier Gloerfelt. "Adaptive relaxation times for aeroacoustic simulations." In: ICMMES 2019. Edimbourgh, UK, 2019

[3] S. Marié, L. Vienne, and X. Gloerfelt. "Consistent filtering for lattice Boltzman simulation of aerodynamic noise." In: 27th International Conference on Discrete Simulation of Fluid Dynamics. Worcester, USA, 2018

[4] L. Vienne, S. Marié, and F. Grasso. "Lattice Boltzmann method for the simulation of miscible gases." In: 27th International Conference on Discrete Simulation of Fluid Dynamics. Worcester, USA, 2018

[5] L. Vienne, S. Marié, and F. Grasso. "Lattice Boltzmann simulation of a multicomponent viscous fingering instability at high Schmidt number". In: 14th ICMMES conference. Nantes, France, 2017

[6] L. Charrier, G. Pont, S. Marié, P. Brenner, and F. Grasso. "Hybrid RANS/LES simulation of a supersonic coaxial He/Air jet experiment with variable turbulent Prandtl and Schmidt numbers". In: 6th HRLM Symposium. Strasbourg, France, 2016

[7] S. Marié and V. Daru. "On the use of adaptive filtering with multiphase lattice Boltzmann models." In: 13th ICMMES conference. Hamburg, Germany, 2016

[8] L. Charrier, M. Jubera, G. Pont, S. Marié, F. Grasso, and P. Brenner. "Simulations of reactive supersonic/subsonic flow interactions for space launcher applications on FLUSEPA solver." In: vol. FP32-2016. 51st 3AF International Conference on Applied Aerodynamics. Strasbourg, France, 2016

[9] S. Marié, H. Lambaré, and P. Druault. "A modal overview of the unsteady loads induced by the Ariane-5 base-flow". In: 4th European Conference for Aerospace Sciences. St-Petersbourg (Russia), 2011

[10] S. Marié and H. Lambaré. "On the unsteady loads induced by the bluff body wake of the Ariane 5 rocket". In: 6th IUTAM Symposium on Bluff Body Wakes and Vortex-Induced Vibrations. Capri (Italy), 2010

[11] Y-H. Dong, S. Marié, and P. Sagaut. "A study of subgrid models in lattice Boltzmann-based large eddy simulation". In: vol. Mechanics Down Under - 11443. 22nd International Congress of Theoretical and Applied Mechanics. Adelaide (Australia), 2008

[12] D. Ricot, S. Marié, and P. Sagaut. "Grid refinement in LBM based on continuous distribution functions". In: 5th ICMMES conference. Amsterdam (Netherlands), 2008

[13] D. Ricot, S. Marié, and P. Sagaut. "Lattice Boltzmann Equation with selective viscosity filters". In: 4th ICMMES conference. Munich (Germany), 2007