# Modélisation des écoulements chargés en particules solides

Christophe Henry

► To cite this version:

Christophe Henry. Modélisation des écoulements chargés en particules solides. Physique Numérique [physics.comp-ph]. Université Cote d'Azur; ED SFA, 2023. tel-04178774

HAL Id: tel-04178774

https://inria.hal.science/tel-04178774

Submitted on 8 Aug 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
<figure><img src="image_2.png" /></figure>

# HABILITATION A DIRIGER LES RECHERCHES

Présentée par **Christophe Henry**

## Modélisation des écoulements chargés en particules solides

Le Jeudi 28 mars 2023 à 10h00

Devant le Jury composé de :

<table><tr><td>Eric Climent</td><td>Institut de Mécanique des Fluide de Toulouse</td><td>Rapporteur</td></tr><tr><td>Rodney Fox</td><td>Iowa State University</td><td>Rapporteur</td></tr><tr><td>Elisabeth Guazzelli</td><td>Université Paris Diderot</td><td>Rapporteur</td></tr><tr><td>Sergio Chibbaro</td><td>Université Paris Saclay</td><td>Examineur</td></tr><tr><td>Hervé Guillard</td><td>Centre Inria d'Université Côte d'Azur</td><td>Examineur</td></tr><tr><td>Nicolas Mordant</td><td>Université Grenoble Alpes</td><td>Examineur</td></tr><tr><td>Aurore Naso</td><td>École Centrale Lyon</td><td>Examineur</td></tr></table>
Remerciements

Ce mémoire est consacré à mes travaux de recherche. La plupart des travaux que je résume ici a été réalisée avec des collègues/partenaires mais également, plus récemment, avec des docto-rant(e)s et post-doctorant(e)s que j'ai eu le plaisir d'encadrer.

Je tiens à remercier toutes les personnes avec qui j'ai eu la chance d'échanger ces dernières années, à la fois sur des aspects scientifiques et aussi sur des notions philosophiques & éthiques. J'ai une pensée toute particulière pour Jean-Pierre Minier, Mireille Bossy, Jérémie Bec et Jacek Pozorski. Plus généralement, je souhaite remercier les laboratoires et instituts qui m'ont accueilli ces dernières années lors de mes pérégrinations scientifiques : le laboratoire MFEE à EDF R&D Chatou et le LECIME à Chimie ParisTech entre 2009 et 2012, le laboratoire Écoulements Multiphasiques à l'Institut des Machines à Écoulement Fluide (IMP PAN) à Gdansk entre 2012 et 2016, puis le laboratoire Lagrange à l'Observatoire de la Côte d'Azur (OCA) à Nice entre 2016 et 2018 et enfin l'équipe Tosca ainsi que l'équipe CaliSto à l'Inria Sophia-Antipolis Méditerranée depuis 2018.

Je tiens également à remercier des personnes que je n'ai pas encore eu la chance de côtoyer, mais dont les écrits m'ont inspiré de nouvelles manières d'être et de penser. J'ai une pensée particulière pour Edgar Morin, Jean-Baptiste Morizot et Alain Damasio. Je continuerai à faire évoluer ces réflexions à la fois personnelles et scientifiques au cours de ces prochaines années et j'espère, en toute humilité, réussir à les transmettre à d'autres personnes pour qu'elles puissent à leur tour s'en inspirer pour tracer leur propre chemin.

Enfin, je remercie ma famille et mes amis (avec ou sans patte) pour tous ces moments de bonheur et de complicité passés ensemble. J'ai une pensée particulière pour mes parents et pour Arnaud, mon frère jumeau.
## Préface

Avant de s'embarquer dans une aventure scientifique, je commence ici par quelques considérations méta-physiques. En bon mordu des jeux sur les mots, leur orthographe et leur sens, je mélange ces différentes notions pour proposer une interprétation à une question qui taraude nombre de scientifiques : quel sens donner à la science ?

Malgré son apparente simplicité, cette question n'a pas de réponse universelle et ce même au sein de la communauté scientifique. Et la raison est simple : cette question soulève un problème qui n'est pas purement scientifique mais aussi méta-physique ou philosophique à. Ainsi, il n'y a pas une seule et unique réponse à cette question, mais plutôt un ensemble de réponses qui varient selon les personnalités et sensibilités de chacun(e). Dans les paragraphes qui suivent, je partage ma vision actuelle sur ce sujet (libre à chacun de s'en inspirer).

Revenons d'abord sur les origines étymologiques du terme<sup>b</sup>. Le mot “science” vient du latin “scientia”, signifiant connaissance. Cette relation avec le terme connaissance se retrouve dans la définition même du terme “science”. Ainsi, selon le dictionnaire Larousse, la science correspond à un “ensemble cohérent de connaissances relatives à certaines catégories de faits, d'objets ou de phénomènes obéissant à des lois et/ou vérifiés par les méthodes expérimentales”. Sciences et connaissances semblent donc intimement liées, ce que je me permets ici d'imager en jouant sur l'orthographe pour former le mot connait-science. L'intérêt de cette ré-écriture est, selon moi, de mettre en valeur le fait que connaissance et science sont le résultat d'une action de la part d'un être vivant. En effet, toujours selon le Larousse, la connaissance est “l'action, le fait de comprendre, de connaître les propriétés, les caractéristiques, les traits spécifiques de quelque chose”. En pratique, ces connaissances s'acquieren grâce à un échange (que ce soit la perception de l'objet d'intérêt ou bien la communication orale/écrite à propos de cet objet) et font aussi appel à des fonctions cognitives (comme la mémoire, le raisonnement, la prise de décision). D'ailleurs, la cognition fait référence à l'ensemble des processus mentaux qui se rapportent à la fonction de connaissance. Étymologiquement parlant, le terme “cognition” est emprunté au latin classique “cognitio”, signifiant “action d'apprendre à connaître”.

<footnote>a. Référence à l'ouvrage d'Edgar Morin, "Science avec conscience", édition Seuil (2017).</footnote>

<footnote>b. Les étymologies citées dans la suite sont tirées du Centre National de Ressources Textuelles et Lexicales, CNTRL (lien vers le site internet).</footnote>
Si l'on aborde maintenant le sujet d'un point de vue historique, on se rend compte que la
“*science*” est intimement liée à la compréhension de la nature. Cela est d'autant plus vrai lorsque
l'on s'intéresse aux modèles scientifiques, c'est-à-dire aux descriptions de phénomènes naturels
via des outils théoriques ou pratiques (je reviendrai sur les modèles dans l'introduction puisqu'il
s'agit du cœur de mes recherches). Plusieurs exemples attestent de ce lien étroit. Un premier
exemple concerne l'histoire des conceptions de l'Univers (cf. les écrits d'Arthur Koestler à ce
sujet c). Cette discipline a commencé dès l'Antiquité avec les premières observations des objets
célestes. Cela a donné lieu à diverses descriptions tentant d'expliquer certains des phénomènes
observés. Dans la directe lignée de ces modèles décrivant les mouvements de corps célestes, on
peut citer un second exemple notable avec les travaux d'Isaac Newton. Ses mémoires s'intitulent
*Philosophiæ naturalis principia mathematica*, ce qui signifie “*Principes mathématiques de la
Philosophie Naturelle*”. Cet ouvrage a posé les fondements de la mécanique dite Newtonienne (aussi
appelée mécanique classique), qui décrit le mouvement d'objets macroscopiques. Un troisième
exemple, plus général, est lié à l'étymologie du terme “*physique*” : ce dernier est issu du grec
“φυσική” (prononcé *phusikós*) signifiant “*Science de la nature*”. Ainsi, les sciences physiques ont
un ancrage sémantique avec la nature.

Au-delà des liens avec la nature, la science pose également des questions plus philosophiques. C'est d'ailleurs le sens premier du terme “*métaphysique*”, qui signifie “*après les choses de la nature*”. On trouve encore des traces de ce lien entre la science et la philosophie actuellement dans la vie académique. Par exemple, le grade universitaire de docteur est souvent abrégé en Ph.D., pour *Philosophy Doctor* en anglais. Un autre exemple est le titre de l'ouvrage d'Isaac Newton, faisant clairement référence à la philosophie. La référence à la philosophie est en fait liée à l'étymologie du mot “*philosophe*”, qui signifie “*celui qui aime la sagesse*” en grec ancien. En ce sens, la science peut être vue comme le fruit d'une pratique (parfois longue et laborieuse) donnant lieu à un savoir.

En abordant ce sujet d'un point de vue de l'individu/du scientifique se posent alors des
questions mêlant philosophie et histoire des sciences comme : quelle(s) motivation(s) pousse(nt)
un individu à étudier un phénomène naturel ? Quel est son objectif ? Dans l'exemple précédent
sur l'histoire des conceptions de l'Univers, certains ont étudié ce sujet suite à la fascination et
l'attirance que certains phénomènes naturels exercent sur l'Homme depuis l'Antiquité. Certains
voulaient prédire ces phénomènes dans un but précis (par exemple pour connaître les saisons
adéquates pour planter certaines graines/semances). Néanmoins, ces exemples montrent aussi les
limites et écueils possibles de ces descriptions : autant que possible, elles se doivent de rester
fidèles aux phénomènes observés sans être influencées par des croyances ou des dogmes. Le
meilleur exemple est l'influence qu'a exercé l'école de Platon sur le mouvement des planètes :
“*la forme du monde doit être une sphère parfaite et tout mouvement doit s'effectuer en cercles parfaits à une vitesse uniforme*”<sup>c</sup>. Cela a été l'un des fondements du système de Ptolémé, qui
a été adopté pendant plusieurs centaines d'années comme le système de référence décrivant le
mouvement des planètes. Malheureusement, la justification de cette hypothèse n'est aucunement
physique mais purement métaphysique, puisque le cercle et la sphère étaient alors considérés
comme des formes parfaites. Le même genre de considérations peuvent amener à sélectionner
une théorie physique plutôt qu'une autre sur des critères plus subjectifs (comme l'esthétique)<sup>d</sup>.

<footnote>c. Référence à l'ouvrage d'Arthur Koestler, "Les Somnanbules : essai sur l'histoire des conceptions de l'Univers", édition Les Belles Lettres (2017).</footnote>

<footnote>d. Référence à l'ouvrage de Sabine Hossenfelder, "Lost in Math : How Beauty Leads Physics Astray.", édition Hachette (2018).</footnote>
En lien avec les paragraphes précédents, je propose une autre manière de ré-écrire le mot co-
naissance. Certes, d'un point de vue purement étymologique, la racine des deux mots n'est pas la
même : le mot “connaissance” est dérivé du latin “nōscō” (qui signifie “familiarisé avec”)
que le mot “naissance” est dérivé du latin “nascor” (qui signifie “naître”). Néanmoins, l'étymologie du
terme “physique” est profondément liée au vivant. En effet, la racine indo-européenne signifiant
“croître, naître, être” a donné le verbe grec ancien “φυω” (prononcé), dont est tiré le terme
“φυσις” (prononcé phúsis) signifiant “nature”. C'est en faisant ce rapprochement étymologique
que m'est venu l'idée de ré-écrire le mot co-naissance.

Ce jeu sur les sens ouvre la réflexion suivante : la connaissance pourrait être vue comme issue
de la naissance simultanée d'un individu et de sa relation au monde qui l'entoure. D'autant plus
qu'un individu (au sens d'une entité isolée) nait nescient, c'est-à-dire sans connaissance. L'intérêt
de cette vision est qu'elle apporte alors un éclairage différent. Dans ce cas, la science (notamment
les sciences physiques) permet de décrire non pas l'Univers comme une entité bien identifiée,
mais bien la perception qu'un individu/un groupe d'individus a de cet univers avec lequel il/elle
interagit. L'observation d'un objet et sa description en terme scientifique devient alors intimement
liée à l'interaction entre l'observateur et l'objet observé. Dans ce sens, le scientifique travaille
pour comprendre et décrire l'ensemble de ces relations et de ces interdépendances e. Pour le dire
autrement, je dirai que la science ne consiste pas essentiellement “à découvrir des faits, mais à
découvrir de nouvelles manières de les penser”f.

En revenant à la question d'origine sur le sens de la science, l'une des réponses que je suis
tenté d'apporter aujourd'hui est la suivante : la science consiste à faire avancer nos connaissances
et donc à repousser les limites de nos méconnaissances. On peut alors se demander : comment ?
Ce à quoi je répondrai : grâce à un dialogue, un échange avec le monde qui nous entoure. Et
surtout en restant ouvert d'esprit, pour accueillir toute la complexité de ces inter-actions g et sans
essayer de les restreindre à des principes génériques (qui, malgré leur aspect parfois séduisant
mathématiquement et/ou rationellement parlant, s'avèrent souvent simplificateurs d'une réalité
bon plus complexe). Cette façon d'observer et de comprendre les choses permet aussi de rester
conscient que nous faisons partie intégrante de l'Univers qui nous entoure, nous en tant qu'êtres
vivants et conscients comme tant d'autres. Ainsi, la formulation d'une théorie ou d'un modèle
relève d'une certaine forme d'art, dont je détaille certaines des facettes dans le chapitre introductif
sur l'art de la modélisation (cf. Chapitre 1).

Enfin, à travers cette notion de dialogue avec la nature (qui a lieu depuis des milliers d'années),
il est beaucoup plus naturel d'admettre que les théories scientifiques évoluent constamment et
ne sont pas figées. C'est l'un des grands thèmes abordés par Karl Popper, selon qui une théorie
scientifique est réfutable/faillible, c'est-à-dire qu'elle n'est considérée comme acceptable que tant
qu'elle n'a pas été réfutée (à l'inverse des dogmes qui sont eux incontestables). Par ailleurs, la
notion même de dialogue avec la nature permet de mieux reconnaître le caractère fini de notre
savoir (cf. théorème d'incomplétude de Kurt Gödel).

Pour résumer cette réflexion étymologique, lexicale, philosophique et métaphysique sur le sens
de la science, je dirais : connaissance rime avec co-nait-sciences.

<footnote>e. Référence à l'ouvrage de Baptiste Morizot, "Manières d'être vivant", éditions Actes Sud (2020).</footnote>

<footnote>f. Référence à l'ouvrage de Lawrence Bragg, "A Short History of Science : Origins and Results of the Scientific Revolution", (1959).</footnote>

<footnote>g. Je ré-utiliserai souvent cette orthographe pour refléter l'importance de l'échange dans le préfixe "inter".</footnote>
Jeu sur les sons et les sens :
Les sens de la science et l'essence de la science

Pour finir ce préambule, je propose cette fois-ci de jouer sur les sons et les sens. A la ma-
nière d'Alain Damasio, je me permets de ré-orthographier certains mots, puis de les mettre en
regard/en relation afin d'ouvrir des pistes de réflexion sur les connexions entre la physique et la
métaphysique (voire la Méta-physique h) :

Conscience, con-science<sup>i</sup>
Décence, des-sens, des-sciences
Déficience, défi & sciences
Croissance, crois-science, crois-science
Prescience, pré-science, presse-science<sup>j</sup>
Patience, pas-science
Puissance, puits-sens, puisse-science

L'essence de la science, c'est de donner du sens à des connaissances qui ont déjà été recensées. Mais sans les encenser! Et en évitant les contre-sens. Sans quoi c'est la déliquescence des sens. En réalité, en sciences, chacun(e) consent à chercher de manière sensée et sensible une forme de quintessence du sens. Ainsi naît science initialement nesciente.

<footnote>h. Référence à l'entreprise Meta (anciennement Facebook).</footnote>

<footnote>i. Référence à l'ouvrage de Rabelais, "Pantagruel" (1532) : "Science sans conscience n'est que ruine de l'âme"</footnote>

<footnote>j. Référence à l'ouvrage de Victor Hugo, "Actes et Paroles, III" (1876) : "Pour moi, l'homme n'est complet que s'il réunit ces trois conditions, science, prescience, conscience".</footnote>
## Organisation du mémoire

J'ai choisi d'organiser ce mémoire de la manière suivante. Le premier chapitre introduit toutes les notions clés pour la compréhension de mes travaux de recherche. Ce 1<sup>er</sup> chapitre est écrit de manière à rendre abordable les sujets traités dans ce mémoire à un plus large public. Il introduit notamment ce que j'appelle par la suite "L'art de la modélisation", c'est-à-dire le processus qui permet d'établir un modèle décrivant un phénomène physique puis de l'implémenter dans un logiciel informatique afin de réaliser des simulations numériques. Les chapitres suivants s'appuient sur les notions introduites dans ce premier chapitre. Chacun des chapitres suivants illustre la démarche de modélisation en l'appliquant à une thématique spécifique parmi celles que j'ai abordées. Par souci de concision, les chapitres 2 à 4 sont des résumés de mes travaux de recherche (les détails étant disponibles dans les annexes qui rassemblent les articles publiés à ce jour).

Le plan détaillé du mémoire est résumé ci-dessous :

**Le chapitre 1** est un chapitre introductif sur les écoulements chargés en particules, leur modélisation et les simulations associées. Il contient une description détaillée de : (a) les phénomènes physique en jeu, (b) l'art de la modélisation (où l'enjeu est de développer des modèles rendant compte de ces phénomènes physiques), (c) les problématiques liées à l'implémentation de ces modèles dans des codes (programmes informatiques) et (d) l'analyse des résultats.

**Le chapitre 2** résume les modèles développés à partir de travaux expérimentaux et théoriques pour rendre compte des phénomènes ayant lieu en proche paroi, notamment le transport près des parois et la remise en suspension de particules déposées sur des parois.

**Le chapitre 3** illustre une démarche de modélisation reposant sur des simulations multi-échelle pour la dynamique de particules non-sphériques, c'est-à-dire où les simulations microscopiques servent à la fois à l'élaboration de modèles macroscopiques (avec plusieurs niveaux) et à leur validation.

**Le chapitre 4** montre comment cette démarche de modélisation peut s'appliquer aux modèles pour l'agglomération/fragmentation de particules (en prenant en compte plusieurs niveaux de description).

**Les annexes** comportent les principaux articles donnant plus de détails sur les méthodologies développées dans les chapitres 2 à 4.
Table des matières

# 1 Introduction aux écoulements à phase dispersée : de la physique à la modélisation 13

1.1 Les phénomènes et mécanismes physiques en jeu 14

1.2 Modélisation de ces phénomènes 21

1.2.1 Introduction à l'art de la modélisation 21

1.2.2 Formulation d'un modèle physique 22

1.2.3 Implémentation numérique 39

1.2.4 Simulations et analyse des résultats 44

1.2.5 Le cycle modèle / implémentation / analyse 46

# 2 Suspensions en milieux confinés : resuspension et phénomènes en proche paroi 49

2.1 Introduction 49

2.2 Phénomènes de resuspension de particules 52

2.3 Dynamique en proche paroi : conditions de bord 57

2.4 Perspectives 60

# 3 Agglomération et fragmentation de particules 63

3.1 Introduction 63

3.2 Simulations directes de l'agglomération 65

3.2.1 Méthodes Monte-Carlo 65

3.2.2 Simulations DNS 67

3.2.3 Méthode à probabilité de collisions 68

3.3 Modèles requérant une fréquence de collision 71

3.3.1 Dynamique de population (PBE) 71

3.3.2 Méthode à probabilité de collisions 73

3.4 Perspectives 76

# 4 Dynamique lagrangienne de particules non-sphériques 79

4.1 Introduction 79

4.2 Simulations microscopiques 82
<table>
  <tbody>
    <tr>
      <td>4.3 Modèle réduit</td>
      <td>87</td>
    </tr>
    <tr>
      <td>4.3.1 Chaînes de billes connectées</td>
      <td>87</td>
    </tr>
    <tr>
      <td>4.3.2 Ellipsoïdes et sphéroïdes</td>
      <td>89</td>
    </tr>
    <tr>
      <td>4.4 Perspectives</td>
      <td>92</td>
    </tr>
    <tr>
      <td><b>5 Conclusion générale</b></td>
      <td><b>95</b></td>
    </tr>
    <tr>
      <td><b>Bibliographie</b></td>
      <td><b>95</b></td>
    </tr>
    <tr>
      <td><b>Appendices</b></td>
      <td><b>111</b></td>
    </tr>
    <tr>
      <td><b>A CV</b></td>
      <td><b>111</b></td>
    </tr>
    <tr>
      <td><b>B Non-sphériques</b></td>
      <td><b>113</b></td>
    </tr>
    <tr>
      <td>B.1 Article sur la dynamique de fibres en turbulence</td>
      <td>113</td>
    </tr>
    <tr>
      <td>B.2 Article sur la dynamique de chaînes inextensibles</td>
      <td>114</td>
    </tr>
    <tr>
      <td><b>C Paroi</b></td>
      <td><b>115</b></td>
    </tr>
    <tr>
      <td>C.1 Article sur les conditions limites et la consistance des modèles</td>
      <td>115</td>
    </tr>
    <tr>
      <td><b>D Resuspension</b></td>
      <td><b>117</b></td>
    </tr>
    <tr>
      <td>D.1 Article décrivant l'implémentation numérique d'un modèle de resuspension</td>
      <td>117</td>
    </tr>
    <tr>
      <td>D.2 Article indiquant les limites des modèles de resuspension</td>
      <td>118</td>
    </tr>
    <tr>
      <td>D.3 Article montrant le rôle des collisions entre particules dans la resuspension</td>
      <td>119</td>
    </tr>
    <tr>
      <td><b>E Agglomération</b></td>
      <td><b>121</b></td>
    </tr>
    <tr>
      <td>E.1 Article décrivant une approche fine pour l'agglomération de particules avec un critère d'efficacité de collision</td>
      <td>121</td>
    </tr>
    <tr>
      <td>E.2 Article décrivant une approche fine pour l'agglomération de particules prenant en compte le résultat de la collision</td>
      <td>122</td>
    </tr>
    <tr>
      <td>E.3 Article décrivant une approche macroscopique pour l'agglomération de particules</td>
      <td>123</td>
    </tr>
    <tr>
      <td>E.4 Article sur l'accrétion de petites particules sur une particule de grande taille</td>
      <td>124</td>
    </tr>
    <tr>
      <td>E.5 Article décrivant un algorithme pour une meilleure précision dans la simulation de l'agglomération</td>
      <td>125</td>
    </tr>
    <tr>
      <td>E.6 Article décrivant un algorithme pour une meilleure précision dans la simulation de l'agglomération</td>
      <td>126</td>
    </tr>
    <tr>
      <td>E.7 Article sur la morphologie d'aggrégats obtenus par des simulations à l'échelle microscopique</td>
      <td>127</td>
    </tr>
    <tr>
      <td><b>F Fragmentation</b></td>
      <td><b>129</b></td>
    </tr>
    <tr>
      <td>F.1 Article sur la fragmentation de fibres en turbulence</td>
      <td>129</td>
    </tr>
  </tbody>
</table>
# Chapitre 1

## Introduction aux écoulements à phase dispersée : de la physique à la modélisation

**Définitions** : Un écoulement est dit *multiphasique* lorsque plusieurs phases immiscibles s'écoulent en même temps (donc en présence d'interfaces entre ces phases). Comme illustré sur la Figure 1.1, il peut s'agir de l'écoulement mélangeant des phases liquides (par exemple des liquides immiscibles comme l'eau ou l'huile), des phases gazeuses (comme l'air) et/ou des phases solides (par exemple du sable ou de la poussière dans l'air).

<figure><img src="image_3.png" /><figcaption>FIGURE 1.1 – Schéma d'un écoulement multiphasique (liquide + gaz + solide).</figcaption></figure>

Dans mes travaux, je m'intéresse notamment aux écoulements à *phase dispersée*, c'est-à-dire lorsque l'une des phases est présente sous forme d'inclusion (occupant des zones délimitées et discontinues de l'espace). Par exemple, c'est le cas de bulles d'air ou de feuilles de menthe contenues dans un verre d'eau. Par opposition, la seconde phase est appelée *phase continue*, puisqu'elle occupe continuement une région de l'espace contigüe (l'eau entre les bulles/feuilles dans l'exemple précédent).

Plus spécifiquement, mes recherches portent sur un cas particulier des écoulements à phase dispersée, à savoir les **écoulements chargés en particules solides** où des objets solides sont en suspension dans un fluide (liquide ou gazeux). Je me suis focalisé sur les **particules colloïdales**, qui ont une taille suffisamment faible pour que le mouvement Brownien ait un effet significatif (typiquement pour des tailles allant de quelques nanomètres à quelques micromètres).
## 1.1 Les phénomènes et mécanismes physiques en jeu

Il existe plusieurs manières d'aborder les problématiques liées aux écoulements à phase dissépée. En effet, plusieurs points de vue peuvent être adoptés, par exemple en s'intéressant à (A) au contexte de l'étude et au domaine d'application, (B) aux enjeux de compréhension et (C) aux phénomènes en jeu. Dans cette section, je donne quelques exemples concrets sur ces manières d'aborder ces problèmes, illustrant ainsi leur complémentarité.

(A) **Contexte et applications** : Les écoulements chargés en particules sont omniprésents dans l'environnement, dans l'industrie et dans la vie quotidienne [52]. En voici quelques exemples :

* La dispersion d'aérosols [32] ou de polluants dans l'atmosphère [41] (comme les pollens, la poussière, le sable, la suie ou les micro- et nano-plastiques). Ces petites particules en suspension dans l'air, abrégées PM pour 'Particulate Matter' en anglais, sont classées en fonction de leur diamètre aérodynamique. Dans le domaine de la métrologie des aérosols, on distingue souvent 2 types de particules [2] : les PM10 ont un diamètre plus petit que 10 µm alors que les PM2.5 sont plus petites que 2.5 µm (appelées aussi particules fines). A noter que d'autres classifications font mention des PM1.0 (aussi appelées particules très fines) et des PM0.1 (aussi appelées particules ultrafines). Cette classification a été introduite pour différencier la façon dont ces particules se comportent dans un écoulement et leur destin (voir les enjeux de compréhension décrits ci-dessous).

* Le charriage d'éléments solides par les rivières et les fleuves jusque dans les océans. Il s'agit notamment d'éléments issus de l'érosion des sols et des montagnes ayant des tailles variables allant de quelques micromètres (comme les limons ou l'argile) à plusieurs millimètres (comme le sable ou les cailloux). Ces objets vont ensuite sédimenter à l'embouchure des rivières, participant ainsi à la formation des deltas de rivières. Les rivières peuvent aussi charrier des polluants (comme les plastiques), notamment lorsqu'elles sont alimentées par des eaux usées provenant des villes. Cette pollution plastique charriée par les rivières contribue à l'accumulation de plastiques dans les océans [81].

* Le transport par les océans de matières inorganiques (comme le sable) et organiques (comme le plancton, qui regroupe à la fois le phytoplancton et le zooplancton) [152].

* L'accumulation de dépôts calcaires, de résidus alimentaires ou de produits de corrosion dans les canalisations. Parfois, cette accumulation peut même aboutir à une obstruction complète de la section de passage du fluide. Cela peut avoir des répercussions sur la vie quotidienne (comme le calcaire dans les machines à laver) mais aussi sur des installations industrielles (comme les échangeurs de vapeur dans les centrales électriques) [62].

* Les processus de filtration, qui permettent de filtrer les résidus contenus dans de l'air ou de l'eau (par exemple, les filtres à café ou bien les filtres à polluants pour limiter l'émission de suies dans les moteurs de voiture).

* Les processus de flottation, qui sont couramment utilisés dans les usines de traitement des eaux usées afin d'éliminer les résidus non-solubles. Ces processus consistent à favoriser la formation de floc, c'est-à-dire de macro-particules composées d'un grand nombre d'éléments liés entre eux, puis d'introduire des micro-bulles d'air qui vont entraîner ces flocs vers la surface où les flocs sont récupérés [22]. Dans ce cas précis, il s'agit même d'écoulements triphasiques (avec trois phases).

* Les processus industriels de fabrication du papier, dans lesquels une suspension de fibres dans de l'eau (aussi appelée pâte à papier) est déposée sur une toile formant ainsi un
enchevêtrement de fibres sur une surface. L'empilement ainsi formé est ensuite pressé puis séché pour obtenir une feuille de papier [85].

Ces différents exemples révèlent déjà toute la complexité des phénomènes en jeu, qui font appel à des domaines scientifiques variés et jouent un rôle sur une large gamme d'échelles spatiales et temporelles, rendant sa modélisation d'autant plus complexe. Dans la suite de ce mémoire, je tâcherai de rendre compte de toute cette complexité en suivant une façon de présenter les choses qui met en regard cette diversité de sujets/phénomènes/domaines scientifiques à la manière dont Edgar Morin en parle, c'est-à-dire en concevant la complexité non pas comme une recette ou une réponse mais plutôt comme un défi et une incitation à penser a.

**(B) Enjeux de compréhension** : Comme le montrent les différents exemples précédents, l'étude des suspensions de particules nécessite d'avoir une bonne compréhension des phénomènes en jeu. En particulier, plusieurs questions se posent parmi lesquelles :

**Q-1 Quelle est l'origine des particules dans l'écoulement ?**

On distingue les sources naturelles et les sources anthropiques (c'est-à-dire issues des activités humaines). La caractérisation des sources de particules dans un écoulement est une étape clé. En effet, selon l'application concernée, il s'agit soit de limiter la quantité de particules en suspension (comme les polluants émis par l'industrie automobile) ou bien de la favoriser (comme pour la quantité de fibres présentes dans la pâte à papier).

Dans les écoulements atmosphériques, les sources naturelles sont nombreuses [30] : on peut citer par exemple l'émission de cendres par l'activité volcanique, l'émission de sable suite à l'action du vent sur des zones côtières ou désertiques, l'émission de pollens ou de spores suite à l'action du vent sur des zones agricoles ou forestières, ou encore l'émission de poussières minérales suite à l'érosion de régions montagneuses. En parallèle, les activités humaines sont responsables de l'émission d'une large gamme de particules dans l'atmosphère, comme les suies issues de processus de combustion (voitures, centrales à charbon, etc.), les poussières issues des freins de véhicules ainsi que les particules issues de l'abrasion et/ou de l'usure des pneus [114, 50], les particules radioactives issues d'incidents ou d'accidents nucléaires [146], les micro- et nano-plastiques émis suite à l'abrasion de textiles synthétiques (notamment lors des lavages et séchages) ou suite à la dégradation de matériaux plastiques (e.g., les pneus ou l'incinération incomplète de déchets) [23].

Dans les milieux marins, les sources naturelles sont souvent liées aux processus d'érosion qui sont responsables de la présence de sédiments (sable, calcaire ou limon) dans les rivières. Les milieux marins sont également peuplés par des composés organiques comme le plancton. Les activités humaines sont responsables de certaines pollutions marines (comme l'accumulation de matériaux plastiques dans les océans, rivières et lacs [7]).

Dans les milieux industriels, les particules sont soit issues de l'érosion des installations (comme l'oxyde de fer) soit issues de l'injection de composés organiques ou inorganiques dans le processus industriel (par exemple les déchets dans les stations d'épuration).

**Q-2 Quelle est la nature des particules ?**

La question de la nature des particules est généralement liée à leur source. Par exemple, la composition chimique des plastiques émis par les activités humaines dépend de la nature des matériaux initiaux, avec plusieurs types de polymères possibles [23] : le polypropylène

<footnote>a. Référence à l'ouvrage d'Edgar Morin, "Science avec conscience", Edition Seuil (2017).</footnote>
(PP), le polyéthylène (PE), le polystyrène (PS), le poly(téréphtalate d'éthylène) (PET), poly(chlorure de vinyl) (PVC). La connaissance de la nature des particules s'avère importante puisqu'elle a un impact direct sur leurs caractéristiques comme : leurs propriétés massiques (comme la densité), mécaniques (déformabilité, rigidité) ou chimiques (composition). Ces propriétés influencent à leur tour la dynamique de ces particules (voir Q-3 ci-dessous) ainsi que leur destin (voir Q-4 ci-dessous) : par exemple, la composition chimique est cruciale pour caractériser la toxicité des particules, surtout si elles sont inhalées/ingérées par des êtres vivants (comme les hommes et les autres animaux b) [122].

**Q-3 Quelle est la dynamique de ces particules ?**

Une fois les particules injectées dans un écoulement se pose la question de leur dynamique, c'est-à-dire de l'évolution temporelle des caractéristiques des particules (comme leur position, leur vitesse, leur taille, etc). Cette dynamique résulte des inter-actions entre les particules et leur environnement (fluide, paroi, etc), qui peuvent induire un mouvement des particules mais aussi des changements de taille ou de masse (e.g., via l'oxidation d'une surface métallique). Les particules peuvent aussi avoir une action sur le fluide lui-même (on parle alors de couplage-retour sur le fluide). Par ailleurs, l'étude des mouvements de particules fait souvent une distinction entre les particules passives et les particules actives, ces dernières étant capables de se mouvoir (comme le zooplancton dans l'océan). Dans les écoulements atmosphériques, la plupart des études se focalisent sur des particules passives et leur mouvement est communément appelé *dispersion*. Les exemples les plus courants sont les polluants émis par une cheminée puis dispersés dans l'air, ou encore la poussière et le sable arrachés aux zones arides puis déplacés par delà les océans (comme le sable du Sahara qui atteint le bassin Amazonien [77]). Le mouvement des particules peut aussi devenir chaotique lorsque l'écoulement est turbulent (on parle alors de dispersion turbulente). Lorsque les particules sont suffisamment massives, leur mouvement peut être lent et discontinu avec des déplacements par bond et/ou roulement successifs (comme les cailloux charriés par les rivières).

**Q-4 Quel est le destin de ces particules ?**

Suite à l'émission et la dispersion de particules, un autre sujet d'études consiste à caractériser ce qu'il advient de ces particules. En effet, la présence de particules en suspension dans un écoulement peut avoir un effet désirable ou indésirable selon les applications. Par exemple, dans le cas de l'encrassement de tuyaux par du calcaire, l'objectif est de savoir à quelle vitesse un tel dépôt se forme, dans quelles régions ce dépôt a lieu de manière privilégiée et si le dépôt est suffisamment robuste pour colmater le tuyau (c'est-à-dire le boucher durablement). Dans le cas de particules fines présentes dans l'atmosphère (comme le pollen), l'un des enjeux actuels est d'étudier l'impact sur la santé humaine des particules inhalées. Ces travaux nécessitent donc d'aborder à la fois les questions de dispersion mais aussi de pénétration dans les voies respiratoires. Or, la façon dont ces particules inhalées pénètrent le système respiratoire dépend de leur taille [1] : i) les particules plus grandes que 10 µm restent principalement dans les voies aériennes supérieures (le nez, la bouche, le pharynx et larynx) ; ii) les PM10 pénètrent dans les voies aériennes inférieures, notamment dans la trachée et jusqu'aux bronches ; iii) les PM2.5 peuvent pénétrer jusque dans les alvéoles pulmonaires. Certaines particules ultrafines (PM0.1) pourraient même franchir la barrière alvéolo-capillaire et se retrouver dans la circulation sanguine.

<footnote>b. Référence à l'ouvrage de Baptiste Morizot, "Manières d'être vivant", éditions Actes Sud (2020).</footnote>
Loin d'être exhaustives, ces questions illustrent déjà la richesse des problèmes mis en jeu dans
le cas des écoulements à phase dispersée.

(C) **Phénoménologie** : Une autre manière d'aborder les problématiques des particules en suspension consiste à caractériser les phénomènes en jeu. Complémentaire des questions citées précédemment, cette manière d'aborder le problème permet de mettre en lumière les mécanismes qui sous-tendent ces phénomènes mais aussi d'identifier et de lister les différentes forces en jeu. Ces mécanismes et forces en jeu sont présentés dans les paragraphes qui suivent.

(C-1) **Mécanismes** : L'intérêt de cette manière d'aborder le problème est de proposer une distinction entre différents mécanismes, chacun d'eux pouvant ensuite être étudié en détail en le prenant de manière indépendante des autres. Cette réduction du problème général à un ou plusieurs sous-problèmes a permis des avancées significatives dans la compréhension fine de chaque mécanisme détaillé ci-dessous. Néanmoins, elle ne doit pas faire oublier que l'enjeu est bien d'aborder le problème pris dans sa globalité et, donc, d'étudier ensuite des situations plus complexes et réalistes qui mettent en œuvre plusieurs mécanismes à la fois et leur inter-actions.

Dans le cas des suspensions de particules, les principaux mécanismes sont listés ci-dessous (voir aussi la figure 1.2) [62, 40]. Évidemment, d'autres mécanismes spécifiques peuvent avoir lieu selon les conditions de l'étude (tels que les réactions chimiques ou la combustion) mais ils ne sont pas détaillés ici par souci de concision et pour rester focalisé sur le mouvement des particules.

<figure><img src="image_4.png" /><figcaption>FIGURE 1.2 – Schéma des 4 principaux mécanismes de particules en suspension.</figcaption></figure>

○ Transport

Le transport de particules correspond aux mouvements de translation et de rotation des particules en suspension sous l'action d'un écoulement. Ces mouvements résultent des différentes forces qu'un fluide exerce sur une particule en suspension. En fait, tout corps placé dans un écoulement subit une force hydrodynamique (aussi appelée force aérodynamique dans les écoulements d'air). Cette force est issue des frottements entre les molécules du fluide (liquide ou gazeux) et du corps en suspension. Elle résulte donc d'une vitesse relative non nulle, c'est-à-dire d'une différence entre la vitesse de l'écoulement et la vitesse de déplacement du corps en suspension. Cela donne donc lieu à un mouvement des particules qui tend à s'accorder à la vitesse du fluide autour de chacune d'elles. Néanmoins, la façon dont les particules s'adaptent à des fluctuations de la vitesse du fluide dépend aussi de leur inertie : ainsi, la vitesse de particules sans inertie s'ajuste instantanément à celle du fluide (comme un grain de poussière ou une plume) alors que les particules très inertielles sont peu affectées par les fluctuations rapides de la vitesse du fluide (comme une boule de pétanque). La dynamique des particules en suspension est également influencée par
la gravité : les particules plus denses que le fluide ont tendance à sédimenter (comme les graviers dans l'eau ou le sable dans l'air) alors que les particules plus légères que le fluide peuvent remonter (comme un bouchon de liège dans l'eau).

### Déposition

Le dépôt de particules a lieu lorsque ces dernières sont transportées par le fluide dans une région suffisamment proche d'une paroi pour que les interactions entre les particules et la paroi deviennent prépondérantes, donnant lieu parfois à leur adhésion. Ce mécanisme a donc lieu dans des milieux confinés uniquement, c'est-à-dire en présence de parois ou de bords qui délimitent l'écoulement fluide. Il s'agit d'un mécanisme qui résulte d'un équilibre entre les forces hydrodynamiques (liées au transport) et les forces physico-chimiques entre les surfaces. Ces forces physico-chimiques sont le plus souvent à courte portée et leur effet est négligeable dès que les particules s'éloignent des parois. Lorsque le dépôt s'accumule de manière trop importante, cela peut donner lieu au colmatage de la section de l'écoulement (par exemple, le bouchage des tuyaux par un dépôt calcaire).

### Resuspension

La resuspension correspond à la remise en suspension de particules précédemment déposées sur une surface sous l'action du fluide. Il s'agit donc du mécanisme inverse à celui de déposition. La resuspension résulte d'une rupture d'équilibre entre, d'un côté, les forces qui maintiennent une particule en adhésion avec la paroi (forces physico-chimiques) et, de l'autre, les forces qui tendent à mettre cette particule en mouvement (forces hydrodynamiques, etc.). Cette remise en suspension est l'une des sources d'émission de particules dans l'atmosphère (comme le sable issu de dunes exposées au vent ou le pollen). Ce mécanisme est très proche de l'érosion/abrasion, où des particules primaires sont directement émises dans l'atmosphère par des sources naturelles ou anthropiques (comme les poussières ou le sable). Il peut également être combiné à d'autres phénomènes donnant lieu à la mise en suspension de particules potentiellement nocives pour les êtres vivants (comme les suies issues des incendies et feux de forêts, donc en présence d'écoulements à forte variation de température).

### Agglomération et fragmentation

L'agglomération de particules a lieu lorsque plusieurs particules transportées par un écoulement se rencontrent et adhèrent entre elles pour former des objets de plus grande taille, appelés **agrégats** ou **agglomérats**. C'est un processus impliquant deux étapes : d'abord, les particules sont transportées par le fluide dans le voisinage l'une de l'autre, puis, lorsqu'elles sont suffisamment proches, elles interagissent via des forces physico-chimiques. Le résultat de ces interactions est soit la formation d'agrégats (lorsque les particules adhèrent entre elles) soit des collisions suivies de rebond (avec ou sans dissipation d'énergie selon que l'on soit dans le cas de collisions élastiques ou inélastiques). Cette agglomération joue un rôle important dans la stabilité des suspensions de colloïdes ou dans la formation de delta de rivières (où les limons charriés par les fleuves tendent à s'agglomérer, facilitant ainsi leur sédimentation).

Le mécanisme inverse correspond à la fragmentation d'agrégats déjà formés sous l'action de forces extérieures. La fragmentation est souvent induite soit par les collisions entre des particules (comme dans les formations de proto-planètes) soit par des forces hydrodynamiques (comme un cisaillement fluide à l'échelle de la particule). Dans certains cas spécifiques, la fragmentation résulte de phénomènes plus complexes (par exemple la fragmentation issue de vibrations et/ou de la cavitation de bulles par ultrasons).
(C-2) Forces en jeu : Une autre lecture possible des phénomènes de mouvement de particules consiste à s'intéresser aux forces qui jouent un rôle dans la dynamique de ces particules plutôt qu'aux mécanismes sous-jacents. En pratique, ces deux points de vue sont souvent associés puisqu'ils sont intimement liés. Dans le cas des écoulements chargés en particules, on distingue typiquement quatre types de forces [62, 40] : (i) les forces hydrodynamiques qui s'exercent entre les particules et le fluide, (ii) les forces externes agissant sur les particules (comme la gravité ou l'électro-magnétisme), (iii) les forces particule-surface et (iv) les forces particule-particule. En général, les deux premières sont responsables du transport et de la dispersion de particules dans un écoulement (sur des distances allant de quelques millimètres à plusieurs kilomètres dans le cas de dispersion atmosphériques). A l'inverse, les deux dernières forces régissent la façon dont ces particules vont soit adhérer à une paroi soit se lier entre elles (avec des forces agissant plutôt sur des distances inférieures au micromètre). C'est cette séparation d'échelle qui justifie parfois le recours à un autre point de vue, où ces phénomènes sont décrits comme résultant d'un processus impliquant deux étapes : l'étape de transport (c'est-à-dire liée aux interactions entre une particule et le fluide ainsi que le champ environnant) et l'étape d'attachement (faisant intervenir les forces physico-chimiques entre les objets, que ce soit des particules et/ou des parois).

Comme le laisse transparaître le paragraphe précédent, l'intérêt de cette décomposition du problème est de lister les différentes forces régissant la dynamique de particules en suspension. Cela s'avère très utile pour classer chacune des forces en fonction de leur domaine scientifique et pour s'inspirer des avancées réalisées dans chacun de ces domaines. Ainsi, les forces hydrodynamiques ont largement été étudiées dans le domaine de l'aérodynamique et de l'hydrodynamique (par exemple pour améliorer le profil aérodynamique des avions ou le profil hydrodynamique des bateaux). De la même manière, les forces entre particules et/ou parois ont été largement étudiées dans le contexte des forces inter-moléculaires (c'est-à-dire entre des molécules) et des forces inter-surfaciques (donc entre des objets composés d'un large nombre de molécules). Grâce aux développements de nouvelles techniques de mesure à l'échelle nanométrique (notamment par AFM pour "Atomic Force Microscope" ou SEM pour "Surface Electron Microscope" en anglais), des progrès significatifs ont récemment permis l'expression de nouvelles formulations pour les forces à courte portée entre deux surfaces (voir le livre [68] pour plus de détails).

<figure><img src="image_5.png" /><figcaption>FIGURE 1.3 – Résumé des types de particules en suspension (classées par taille et par nature).</figcaption></figure>
**Résumé : un sujet trans-disciplinaire et multi-échelle.** A travers ces exemples et cette introduction à la phénoménologie des écoulements diphasiques, il est évident que la compréhension de ces systèmes requiert des compétences sur des sujets trans-disciplinaires comme : la mécanique des fluides, la turbulence, les écoulements diphasiques, la physico-chimie des surfaces ou encore la mécanique des matériaux.

Par ailleurs, ces écoulements chargés en particules font souvent intervenir toute une gamme d'échelles, allant de l'échelle nanométrique (avec les forces inter-surfaciques) jusqu'à l'échelle kilométrique (avec le transport d'aérosols sur des distances inter-continentales). La figure 1.3 illustre cette pluri-disciplinarité et ce large éventail en résumant la nature et les tailles de particules typiquement rencontrées dans de tels écoulements diphasiques. De la même manière, la figure 1.4 résume la terminologie rencontrée dans le domaine des écoulements à phase dispersée en fonction de la concentration volumique de particules en suspension dans l'écoulement turbulent ainsi que la gamme de phénomènes ayant lieu (surtout pour la dynamique des particules).

<figure><img src="image_6.png" /><figcaption>FIGURE 1.4 – Résumé des phénomènes ayant lieu en fonction de la concentration de particules en suspension (pour une taille et une densité des particules données).</figcaption></figure>

Mes travaux de recherches s'inscrivent dans ce contexte multi-dimensionnel (au sens à la fois scientifique, spatial et temporel). En particulier, mes travaux consistent à développer des modèles permettant de rendre compte de ces phénomènes multi-échelles et multi-physiques. Ces dernières années, je me suis surtout intéressé aux suspensions diluées (cf. zone rouge dans la figure 1.4) et notamment à la dynamique de particules colloïdales (zone en rouge dans la figure 1.3) en suspension dans des écoulements turbulents.
## 1.2 Modélisation de ces phénomènes

L'objectif de cette section est d'apporter un éclairage sur la question : "Qu'est-ce qu'un modèle ?". Pour ce faire, je propose d'abord une introduction générale à ce que j'appelle l'art de la modélisation (cf. Section 1.2.1), en identifiant les grandes étapes dans la formulation d'un modèle. Une fois ces étapes identifiées, je donne plus de détails sur chacune d'elles dans le contexte des écoulements à phase dispersée (cf. Sections 1.2.2-1.2.4). Enfin, dans la Section 1.2.5, j'illustre le caractère évolutif des modèles en mettant en lumière les relations entre ces étapes.

### 1.2.1 Introduction à l'art de la modélisation

**Définitions** : Avant tout, il est important de rappeler quelques-unes des définitions d'un modèle dans un contexte scientifique.

* D'un point de vue étymologique, le mot *modèle* provient de l'italien *modello*, lui-même issu du latin vulgaire *modellus* qui est dérivé du latin *modulus* et signifie mesure. Un modèle est ainsi lié à la mesure d'un phénomène.

* Selon le dictionnaire Larousse c, un modèle scientifique est une représentation abstraite des relations entre les paramètres caractéristiques d'un phénomène ou d'un processus. Ainsi, il vise à simplifier le phénomène étudié [...], à unifier des phénomènes a priori sans rapport [...] et à prédire des phénomènes nouveaux [...].

* Selon l'Académie Française d, le terme *modèle* fait référence à une représentation physique, graphique ou, plus généralement, mathématique qui formalise les relations unissant les différents éléments d'un système, d'un processus, d'une structure, en vue de faciliter la compréhension de certains mécanismes ou de permettre la validation d'une hypothèse.

Ces définitions illustrent le fait qu'un modèle englobe plusieurs notions formant un triptyque :

* La formulation d'un modèle physique :
  Le modèle physique correspond à une représentation simplifiée utilisée pour expliquer le comportement d'un système physique ou le fonctionnement d'un processus physique. En pratique, l'élaboration d'un tel modèle implique de prendre en compte un certain nombre de choix/hypothèses de travail. Certains de ces choix/hypothèses découlent de considérations théoriques (comme le choix des mécanismes et des forces en jeu), d'autres des applications étudiées (avec le domaine de validité attendu), tandis que d'autres renvoient à des questions pratiques (par exemple, l'applicabilité d'un modèle dans des situations réelles ou bien la précision requise).

* Sa mise en oeuvre (implémentation numérique) :
  La résolution des équations d'un modèle peut se faire soit de manière analytique (avec ou sans approximation supplémentaire) soit en ayant recours à des méthodes numériques. C'est dans le cadre de la résolution numérique de telles équations que de nombreux logiciels ont vu le jour ces dernières décennies. Cela a donné lieu au développement de toute une gamme d'algorithmes avec divers niveaux de précision ainsi qu'au développement d'outils de calcul de plus en plus efficaces et rapides (par exemple en ayant recours à la parallélisation des tâches effectuées).

<footnote>c. Voir le dictionnaire Larousse (lien vers le site internet)</footnote>

<footnote>d. Voir le dictionnaire de l'Académie Française (lien vers le site internet)</footnote>
* L'analyse des résultats :
Une fois le modèle physique implémenté dans un logiciel, des simulations numériques peuvent être réalisées dans des cas idéalisés et/ou réalistes. Il s'agit alors d'analyser les résultats obtenus et notamment de confronter les données issues de simulations aux données théoriques/analytiques/expérimentales dans le but de valider/invalider un ou plusieurs aspect(s) du modèle.

### L'art de la modélisation.
La modélisation désigne l'action de formuler et développer un modèle. Cependant, sous cette définition faussement triviale se cache une tâche bien plus riche et complexe. Par exemple, Bonate mentionne dans son introduction au livre intitulé *The art of modelling* [16] : "La modélisation est plus que l'acte de modéliser. Pour être un modélisateur efficace, il faut comprendre les différents types de modèles et savoir quand un type de modèle est plus approprié qu'un autre, comment sélectionner un modèle dans une famille de modèles similaires, comment évaluer la qualité de l'ajustement d'un modèle, comment présenter les résultats de la modélisation à la fois verbalement et graphiquement à d'autres personnes, et accomplir toutes ces tâches dans un cadre éthique." e. Pour utiliser une autre analogie, il est aussi possible de voir la modélisation comme l'exploration d'un territoire partiellement connu, partiellement inconnu voire même partiellement inexploré. Le modèle devient alors similaire à une carte sur laquelle chacun peut s'appuyer pour atteindre une certaine destination ou suivre un chemin spécifique à travers ce territoire. Dans cette analogie, un modèle a donc peu de valeur intrinsèque lorsqu'il est pris en tant que tel. Son véritable intérêt/utilité ne peut s'évaluer que par rapport aux raisons mêmes pour lesquelles nous choisissons de nous tourner vers ce modèle.

Ainsi, la modélisation exige de combiner et d'équilibrer des aspects très différents, couvrant
à la fois des questions théoriques, numériques et algorithmiques. Cela donne lieu à un certain
nombre de paramètres/contraintes à respecter (comme un cahier des charges). La modélisation
consiste donc à faire des concessions/compromis entre ces différentes contraintes. En ce sens, le
développement d'un modèle s'apparente à un art, que j'appelle "l'art de la modélisation".

Dans les paragraphes qui suivent, je donne plus de détails sur cet art de la modélisation en
revenant sur chacune des trois étapes d'un modèle : sa formulation (cf. Section 1.2.2), sa mise en
oeuvre (cf. Section 1.2.3) et son analyse (cf. Section 1.2.4).

## 1.2.2 Formulation d'un modèle physique

La première étape dans la modélisation consiste à formuler/développer un modèle permettant
de reproduire un phénomène physique.

Néanmoins, cette tâche nécessite de faire des concessions/compromis entre les différentes
contraintes ayant trait aux enjeux de l'étude, aux objectifs, aux théories/hypothèses faites ainsi
qu'aux niveaux de description et de précisions requis. Cette large variété des éléments à prendre
en compte dans la formulation d'un modèle donne lieu à plusieurs grilles de lecture possibles.
Chacune d'elles varie en fonction de l'ordre utilisé pour introduire ces éléments/entrées.

Dans la suite, je détaille l'une des grilles de lecture possible, qui regroupe plusieurs éléments
ordonnés de la manière qui me paraît la plus en adéquation avec mes propres convictions f. Pour

<footnote>e. Traduit de l'anglais : "There is more to modeling than the act of modeling. To be an effective modeler requires understanding the different types of models and when one type of model is more appropriate than another, how to select a model from a family of similar models, how to evaluate a model's goodness of fit, how to present modeling results both verbally and graphically to others, and to do all these tasks within an ethical framework"</footnote>

<footnote>f. A noter que cette grille de lecture s'inspire de travaux récents sur un modèle de données commun, travaux</footnote>
le lecteur, il ne s'agit donc surtout pas de la prendre comme une grille de lecture universelle (c'est-à-dire qui devrait être suivie coûte que coûte), mais plutôt comme l'illustration d'une des façons de procéder (parmi tant d'autres possibilités). Ainsi, pour paraphraser Edgar Morin g,
cette grille de lecture ne doit pas être perçue comme une recette toute prête mais plutôt comme un défi et une incitation à penser.

## 1<sup>er</sup> élément : Domaine d'utilisation et d'application

La première chose à faire consiste à identifier clairement le type de problèmes que l'on souhaite résoudre grâce à ce modèle. Cela permet d'identifier et de lister les phénomènes que l'on souhaite prendre en compte et ceux qui sont négligeables. Par exemple, si l'on cherche à reproduire la vitesse à laquelle un dépôt de particules se forme dans un tuyau en présence d'échanges de chaleur, il est évident que les modèles choisis doivent décrire l'écoulement d'un fluide dans un tuyau, le transport de particules solides mais aussi l'effet des échanges de chaleur et variations de température sur cet écoulement et sur la dynamique des particules.

Le choix des applications visées peut avoir des conséquences sur la plage de validité du modèle, surtout en fonction des hypothèses qui vont être faites (voir élément 2). Par exemple, si l'on choisit de se focaliser sur les suspensions de particules dans des écoulements à température constante et qu'aucun modèle n'est inclus pour les échanges thermiques, le modèle ne sera applicable sans modification qu'à des écoulements isothermes. De la même manière, si l'on choisit de se focaliser sur l'étude de particules de très grande taille (pour lesquelles le mouvement Brownien devient négligeable), le modèle ne pourra pas être appliqué tel quel au cas de colloïdes (qui sont sensibles au mouvement Brownien).

## 2<sup>nd</sup> élément : Théories menant à une mise en équation.

Une fois le problème et le domaine d'application bien définis, la seconde étape nécessite de choisir les théories qui mettent en équation ces phénomènes, sans oublier les hypothèses qui ont permis de formuler ces équations. Cette étape est plus délicate pour les suspensions de particules puisqu'il y a des aspects théoriques à prendre en compte pour chaque phase :

* La phase fluide

<figure><img src="image_7.png" /><figcaption>FIGURE 1.5 – Schéma des propriétés d'un fluide et du type d'écoulement.</figcaption></figure>

<footnote>réalisés en partenariat avec Jean-Pierre Minier dans le cadre du projet Européen VIMMP (site web)</footnote>

<footnote>g. Référence à l'ouvrage d'Edgar Morin, "Science avec conscience", Edition Seuil (2017).</footnote>
Il existe plusieurs formulations théoriques pour décrire le comportement d'un fluide en écoulement. En particulier, ces formulations dépendent d'un certain nombre de propriétés du fluide et du type d'écoulement qui sont résumés ci dessous (voir Figure 1.5) :

- Le type d'écoulement :

A l'échelle macroscopique (par opposition à l'échelle atomique/moléculaire), l'écoulement d'un fluide est décrit grâce aux équations de Navier-Stokes. Ces équations résultent d'une description de type mécanique des milieux continus, c'est-à-dire qu'elles décrivent des champs (comme le champ de vitesse, le champ de pression). En fait, ces équations résultent de la combinaison de trois relations de conservation qui s'écrivent sous la forme [119] :

$$ \frac{\partial \rho_f}{\partial t} + \nabla \cdot (\rho_f \mathbf{U}_f) = 0 \qquad (1.1a) $$

$$ \frac{\partial (\rho_f \mathbf{U}_f)}{\partial t} + \nabla \cdot (\rho_f \mathbf{U}_f \mathbf{U}_f) = -\nabla P + \nabla \cdot \mathcal{T}_f + \mathbf{F} \qquad (1.1b) $$

L'équation (1.1a) découle de la conservation de la masse, qui permet de relier la variation temporelle de la densité du fluide $\rho_f$ au flux de masse (qui fait intervenir la vitesse du fluide $\mathbf{U}_f$). Elle est couramment appelée équation de bilan de masse ou bien équation de continuité. L'équation (1.1b) correspond au bilan de quantité de mouvement. Issue de la 2$^{nde}$ loi de Newton, cette équation relie l'accélération du fluide aux forces surfaciques et volumiques s'exerçant sur le fluide. Elle fait intervenir la pression $P$, le tenseur des contraintes visqueuses $\mathcal{T}_f$ et les forces extérieures $\mathbf{F}$ (comme la gravité ou parfois l'influence des particules sur le fluide).

Ces équations de Navier-Stokes permettent de décrire naturellement tout type d'écoulement, allant des écoulements laminaires aux écoulements turbulents en passant par tous les régimes transitoires. Néanmoins, la résolution de ces équations dans le cas turbulent n'est pas chose aisée, notamment à cause des termes non-linéaires dus au tenseur des contraintes visqueuses et de la nature chaotique des écoulements turbulents (caractérisés par des recirculations et tourbillons à plusieurs échelles) [119, 44]. Ainsi, la simulation d'écoulements turbulents introduit des questions de modélisation liées au niveau de description choisi pour simplifier la résolution de ces équations (cf. élément 4).

Néanmoins, les équations de Navier-Stokes ne sont pas nécessairement suffisantes pour décrire un écoulement. En effet, il existe des types d'écoulement faisant intervenir d'autres effets comme :

▷ Les effets thermodynamiques

Les écoulements sont souvent considérés comme des processus isothermes, c'est-à-dire ayant lieu à une température constante. Néanmoins, dans certaines situations, la température du fluide est variable en temps et/ou en espace. C'est le cas par exemple dans les procédés d'échange de chaleur, où l'écoulement d'un fluide froid va se réchauffer au contact d'une surface chauffée (jusqu'à entrer en ébullition). Dans de telles situations, les variations de température ont un effet sur l'écoulement fluide.

Il est alors nécessaire d'adapter les lois de comportement du fluide pour rendre compte de l'effet sur le fluide d'une variable thermodynamique, laquelle a été préalablement choisie parmi les potentiels thermodynamiques (e.g., enthalpie,
énergie interne) ou les variables naturelles (e.g., température, entropie). Dans le cas d'écoulements thermodynamiques, il est courant d'ajouter une équation sur le bilan de l'énergie. Par exemple, il est possible de décrire l'énergie totale par unité de masse E, considérée comme la somme de l'énergie interne par unité de masse e et de l'énergie cinétique (soit $E = e + \frac{|\mathbf{U}_f|^2}{2}$). Cela implique que le potentiel thermodynamique choisi est l'énergie interne et son évolution est donnée par :

$$ \frac{\partial (\rho_f E)}{\partial t} + \nabla \cdot (\rho_f E \mathbf{U}_f) = -\nabla \cdot (P \mathbf{U}_f) + \nabla \cdot (\mathcal{T}_f \cdot \mathbf{U}_f) + \mathbf{F} \cdot \mathbf{U}_f + \nabla \cdot \mathbf{q} \quad (1.2) $$

L'équation (1.2) fait intervenir le flux de chaleur $\mathbf{q}$, qui peut provenir soit de la conduction thermique (lorsqu'il y a des échanges avec des surfaces chauffées/refroidies) soit du rayonnement (e.g., dans les écoulements atmosphériques).

### ▷ Les effets de réactions chimiques

Les écoulements réactifs correspondent à un type particulier d'écoulements dans lesquels des espèces chimiques sont transportées et déclenchent des réactions chimiques. Ces réactions peuvent induire des variations de température. En effet, les réactions sont soit endothermiques (c'est-à-dire qu'elles captent de la chaleur) ou exothermiques (qui libèrent de la chaleur). C'est le cas dans les flammes de combustion, où l'air transporte de l'oxygène qui va réagir avec le gaz injecté au niveau de la flamme et libérer ainsi de la chaleur.

Pour rendre compte de l'effet de ces réactions chimiques sur l'écoulement, il est nécessaire d'inclure une équation sur le bilan d'énergie (comme dans les écoulements thermodynamiques) et de les coupler à des équations décrivant chacunes des réactions chimiques ayant lieu (incluant les informations sur les espèces chimiques et la cinétique des réactions).

### ▷ Les effets de magnéto-hydrodynamique

Lorsque le fluide a des propriétés électromagnétiques (comme dans un plasma ou dans un métal liquide), l'écoulement du fluide est régie par les phénomènes hydrodynamiques et par les phénomènes électro-magnétiques (c'est-à-dire liés à la présence dans l'écoulement d'éléments ayant des propriétés électriques et/ou magnétiques).

Dans ce cas, les équations de magnéto-hydrodynamique (MHD) se couplent aux équations de Navier-Stokes (pour la dynamique des fluides) aux équations de Maxwell (pour la dynamique des éléments/particules chargés et/ou magnétiques) [27].

### - Les lois de comportement :

Les lois de comportements visent à modéliser le comportement des fluides par des lois établies en considérant le fluide comme un milieu continu. En particulier, on distingue plusieurs aspects comme :

#### ▷ Compressibilité / incompressibilité :

Tout fluide est compressible, au sens où un changement de pression affecte la densité du fluide. Cependant, certains fluides sont souvent considérés comme incompressibles puisque les variations de densité en fonction de la pression et de la température sont négligeables. C'est le cas, par exemple, de l'eau. Certains
fluides sont également qualifiés de fluides dilatables, lorsque la densité varie en fonction de la température.

D'un point de vue théorique, la différence entre un fluide compressible et incompressible a un effet sur l'équation de bilan de masse (1.1a). En effet, dans le cas d'un écoulement incompressible, c'est-à-dire sous l'hypothèse que la densité $\rho_f$ est constante, l'équation (1.1a) se simplifie en :

$$ \nabla \cdot \mathbf{U}_f = 0 \qquad (1.3) $$

Un écoulement considéré comme incompressible est donc aussi appelé un écoulement à divergence nulle (sous entendu, la divergence de la vitesse du fluide).

### ▷ Fluide Newtonien ou non-Newtonien :

La distinction entre un fluide Newtonien et un fluide non-Newtonien a un impact sur l'expression de la viscosité du fluide. Dans un écoulement Newtonien (comme l'air ou l'eau), la contrainte de cisaillement évolue linéairement avec le gradient de vitesse (la viscosité étant la constante de proportionnalité). Il s'agit donc de fluides pour lesquels la viscosité est constante et indépendante de la contrainte de cisaillement. A l'inverse, les écoulements non-Newtoniens ne sont pas cette règle. Plus précisément, il y a toute une gamme d'écoulements non-Newtoniens, incluant les écoulements visco-élastiques (dont la viscosité et l'élasticité changent selon les contraintes, comme certains lubrifiants), les écoulements à viscosité variable en temps (où la viscosité apparente augmente ou diminue à mesure qu'une contrainte est appliquée, comme pour le cas du yaourt ou du beurre de cacahouète) et les écoulements à viscosité non-Newtonienne (où la viscosité augmente ou diminue lorsque la vitesse du fluide augmente, comme par exemple le ketchup ou le sang).

Dans un écoulement Newtonien, le tenseur des contraintes visqueuses est proportionnel au gradient de vitesse et s'écrit sous la forme $T_{f,ij} = \mu_f \left( \frac{\partial U_{f,i}}{\partial x_j} + \frac{\partial U_{f,i}}{\partial x_j} \right)$ (où $\mu_f$ est la viscosité dynamique du fluide). Cela permet de simplifier l'équation de bilan de quantité de mouvement (1.1b), qui s'écrit alors :

$$ \frac{\partial \mathbf{U}_f}{\partial t} + \mathbf{U}_f \cdot \nabla \mathbf{U}_f = -\frac{1}{\rho_f} \nabla P + \nu_f \nabla^2 \mathbf{U}_f + \mathbf{F} \qquad (1.4) $$

où $\nu_f$ est la viscosité cinématique du fluide (donnée par $\nu_f = \mu_f / \rho_f$).

### – Les contraintes spécifiques

L'écoulement d'un fluide dépend aussi d'un certain nombre de contraintes, notamment des contraintes spatiales (comme la géométrie du domaine dans lequel le fluide s'écoule et les conditions aux limites) ou bien des contraintes temporelles (e.g., les conditions initiales). Ces contraintes supplémentaires doivent être précisées afin d'être en mesure de reproduire certaines des variables (comme le champ de vitesse ou de pression).

Dans le cas des conditions aux limites, il existe plusieurs types de condition. Une première possibilité consiste à avoir recours à une condition aux limites de Dirichlet, qui consiste à imposer la valeur que la solution doit vérifier sur une frontière du domaine (par exemple, une vitesse nulle sur une paroi fixe). Une seconde option correspond aux conditions aux limites de Neumann, où l'on spécifie les valeurs que
prennent les dérivées de la solution sur une frontière du domaine (e.g., un flux nul en entrée d'un tuyau). Enfin, les conditions aux limites de Cauchy consistent à imposer à la fois une condition de Neumann et une condition de Dirichlet sur une frontière/limite du domaine.

o La phase dispersée

Tout comme pour la phase fluide, la dynamique de la phase dispersée dépend de plusieurs paramètres, notamment la nature des particules, leurs interactions avec le milieu ambiant et les interactions entre les particules elles-mêmes (voir la Figure 1.6).

<figure><img src="image_8.png" /><figcaption>FIGURE 1.6 – Schéma des propriétés de particules et de leurs interactions avec le milieu ambiant.</figcaption></figure>

Dans la suite, on s'intéresse à la dynamique de particules solides, considérées comme des objets dont on ne décrit pas la structure interne mais seulement les propriétés (comme par exemple sa masse, sa densité, sa forme). Pour ces raisons, les approches de type dynamique moléculaire (MD pour Molecular Dynamics en anglais) ne sont pas considérées puisqu'elles décrivent une particule comme un ensemble d'atomes/molécules reliés entre eux. À l'échelle macroscopique, la dynamique des particules est donnée directement par la $2^{nde}$ loi de Newton (c'est-à-dire le bilan de quantité de mouvement) :

$$ \frac{d\mathbf{X}_p}{dt} = \mathbf{U}_p \qquad (1.5) $$

$$ m_p \frac{d\mathbf{U}_p}{dt} = \sum \mathbf{F} \qquad (1.6) $$

$$ I_p \frac{d\mathbf{\Omega}_p}{dt} = \sum \mathbf{T}_F \qquad (1.7) $$

où $\mathbf{X}_p$ est la position des particules, $\mathbf{U}_p$ la vitesse de translation, $\mathbf{\Omega}_p$ la vitesse de rotation, $m_p$ sa masse, $I_p$ l'inertie de rotation et $\mathbf{F}$ (resp. $\mathbf{T}_F$) les forces (resp. moment des forces) s'exerçant sur la particule. La première équation décrit simplement l'évolution de la position des particules alors que la seconde équation (resp. la troisième) décrit la vitesse de translation de la particule (resp. la vitesse de rotation).

En fait, le choix des équations à résoudre découle du vecteur d'état retenu, c'est-à-dire de l'ensemble des caractéristiques associées aux particules que l'on souhaite suivre. Dans l'équation (1.7), le vecteur d'état $\mathbf{Z}_p$ se compose de la position, de la vitesse et de la rotation des particules $\mathbf{Z}_p = (\mathbf{X}_p, \mathbf{U}_p, \mathbf{\Omega}_p)$. Selon les situations, le vecteur d'état et les équations associées peuvent être simplifiées ou enrichies selon la nature des particules considérées. Par exemple, la rotation est souvent laissée de côté lorsque l'on
s'intéresse au cas de particules sphériques à la fois solides et non-déformables. Dans ce cas, le vecteur d'état se simplifie en $\mathbf{Z}_p = (\mathbf{X}_p, \mathbf{U}_p)$. À l'inverse, il est nécessaire de rajouter des équations supplémentaires pour rendre compte des déformations (comme l'extension, la flexion ou la torsion de particules allongées tels les polymères). Cela sera justement le thème du Chapitre 4, qui illustre comment ce vecteur d'état peut être complété pour rendre compte de l'orientation de sphéroïdes.

L'une des difficultés principales lorsque l'on établit un modèle pour la dynamique de ces particules est liée aux forces que l'on souhaite retenir dans l'expression de leurs mouvements. Comme indiqué sur la Figure 1.6, quatre types d'interactions ont lieu :

- Interactions particule-fluide :

Tout corps placé dans un écoulement fluide subit des forces hydrodynamiques liées aux interactions avec les molécules du fluide. Ainsi, il est possible d'estimer ces forces hydrodynamiques en intégrant, le long de toute la surface de l'objet, le tenseur des contraintes entre le fluide $\mathcal{T}_{hydro}$, ce qui donne :

$$ \mathbf{F}_{hydro} = \int_{S_p} \mathcal{T}_{hydro} \cdot \mathbf{n} \, dA \qquad (1.8) $$

où $\mathbf{n}$ est la normale sortante à la surface locale et $S_p$ la surface de l'objet. D'ailleurs, en suivant le principe d'action-réaction, c'est l'opposé de ces forces particules-fluide qui entrent dans l'équation de Navier-Stokes pour rendre compte de l'effet des particules sur l'écoulement (cf. Eq. (1.1b)). La prise en compte de ces effets est souvent appelé couplage retour sur le fluide ("two-way coupling" en anglais).

En pratique, il est fréquent d'avoir recours à un modèle simplifié pour les forces hydrodynamiques en approximant l'intégration du tenseur des contraintes par une formule analytique. Ce genre d'approche est courant dans les modèles plus macroscopiques (cf. la discussion sur le niveau de description dans l'élément 4) et également lorsque l'on considère le cas de particules ponctuelles [79], c'est-à-dire ayant une taille bien inférieure à la plus petite échelle de l'écoulement fluide (typiquement l'échelle de longueur de Kolmogorov $\tau_\eta$ dans les écoulements turbulents). L'intérêt des particules ponctuelles est que les forces hydrodynamiques peuvent se calculer en utilisant une formule simplifiée ne faisant intervenir que la vitesse du fluide à la position de la particule, appelée dans la suite la vitesse du fluide vue $\mathbf{U}_p = \mathbf{U}_f(\mathbf{X}_p(t), t)$ (comme dans [100, 103]).

Malgré ces hypothèses simplificatrices et des décennies de recherche sur le sujet (voir les travaux de Stokes au milieu du 19ème siècle puis ceux de Bassett-Boussinesq-Oseen entre 1888 et 1927), l'expression exacte de ces forces hydrodynamiques est encore aujourd'hui un sujet de recherche. Dans le cas de particules ponctuelles, les forces hydrodynamiques retenues sont généralement données par l'équation dite de Maxey-Riley-Gatignol [94, 47] (ou équation BBO). Cette équation prend en compte plusieurs contributions [100, 79, 3, 18] :

$$ \mathbf{F}_{hydro} = \mathbf{F}_{drag} + \mathbf{F}_{buoy} + \mathbf{F}_{PG} + \mathbf{F}_{AM} + \mathbf{F}_{hist} \qquad (1.9) $$

▷ La force de flottabilité $\mathbf{F}_{buoy}$ (buoyancy en anglais), provenant de la différence de densité entre l'objet en immersion et le fluide. Elle correspond à la différence
entre la force de gravité liée à la masse de la particule et la poussée d'Archimède (i.e., proportionnelle au volume d'eau déplacé par l'objet immergé). Elle s'exprime sous la forme :

$$ \mathbf{F}_{buoy} = \mathcal{V}_p (\rho_p - \rho_f) \mathbf{g} \qquad (1.10) $$

où $\mathcal{V}_p$ est le volume de la particule, $\rho_p$ sa densité, $\rho_f$ celle du fluide et $\mathbf{g}$ l'accélération de la pesanteur à la surface de la Terre. Cette force tend donc à faire remonter à la surface les particules plus légères que le fluide (comme un bouchon en liège dans l'eau) ou à faire couler les particules plus lourdes (comme un caillou dans l'eau).

▷ La force de traînée $\mathbf{F}_{drag}$ (drag en anglais), qui s'oppose au mouvement d'un corps dans un fluide et agit comme un frottement. Cette force s'exprime en général par rapport à la vitesse relative entre la particule et le fluide $\mathbf{U}_r = \mathbf{U}_s - \mathbf{U}_p$ :

$$ \mathbf{F}_{drag} = \frac{1}{2} \rho_f C_D A_p |\mathbf{U}_r| \mathbf{U}_r \qquad (1.11) $$

où $A_p$ est la surface de la particule exposée à l'accoulement et $C_D$ est le coefficient de traînée. Grâce aux avancées dans le domaine de l'aéro- et hydro-dynamique, il existe des formules empiriques approximant ces coefficients pour des objets de forme variée dans des écoulements isotropes [24, 19]. Par exemple, dans le cas de particules sphériques, l'une des formules empiriques couramment utilisée est :

$$ C_D = \begin{cases} \frac{24}{\mathcal{R}e_p} (1 + 0.15 \mathcal{R}e_p^{0.687}) & \text{si } \mathcal{R}e_p \le 1000 \\ 0.44 & \text{si } \mathcal{R}e_p \ge 1000 \end{cases} \qquad (1.12) $$

où $\mathcal{R}e_p = d_p|\mathbf{U}_r|/\nu_f$ est le nombre de Reynolds particulière lié à la vitesse relative $\mathbf{U}_r$. D'ailleurs, quand $\mathcal{R}e_p \ll 1$, on retrouve l'expression simplifiée de la force de traînée dans un écoulement de Stokes : $\mathbf{F}_{drag} = 3\pi d_p \rho_f \nu_f \mathbf{U}_r$. Néanmoins, il n'existe pas de consensus sur une formule approximant ces forces lorsque les particules sont proches d'une paroi.

▷ La force de gradient de pression $\mathbf{F}_{PG}$, qui correspond à la force exercée par le champ de vitesse non perturbé par la présence de la particule [79]. Elle s'exprime sous la forme :

$$ \mathbf{F}_{PG} = \mathcal{V}_p \rho_f \frac{D\mathbf{U}_s}{Dt} \qquad (1.13) $$

où $D/Dt = \partial/\partial t + \mathbf{U} \cdot \nabla$ est la dérivée particulière.

▷ La force de masse ajoutée $\mathbf{F}_{AM}$, issue de l'accélération/décélération du volume de fluide qui est déplacé par la particule en mouvement [83]. Elle correspond donc à l'une des contributions liées à la vitesse perturbée par la présence d'un objet en mouvement dans le fluide. Elle s'exprime sous la forme :

$$ \mathbf{F}_{AM} = \frac{1}{2} \rho_f \mathcal{V}_p \left( \frac{D\mathbf{U}_s}{Dt} - \frac{d\mathbf{U}_p}{dt} \right) \qquad (1.14) $$

▷ La force de Basset $\mathbf{F}_{hist}$, qui résulte de l'effet visqueux de l'accélération relative de l'accoulement non perturbé à l'emplacement de la particule par rapport à l'accélération de la particule. Cette force correspond ainsi à la force apparente exercée sur une particule dont la trajectoire est accélérée dans un écoulement. Elle
est également appelée force d'histoire de Basset, puisque l'écoulement autour de l'objet dépend de l'histoire de la particule (notamment de sa trajectoire). Elle s'exprime sous la forme d'une intégration le long de chaque trajectoire :

$$ \mathbf{F}_{hist} = \frac{3}{2} d_p^2 \rho_f \sqrt{\pi \nu_f} \int_{-\infty}^{t} \frac{ds}{ds} (\mathbf{U}_s(s) - \mathbf{U}_p(s)) \frac{ds}{\sqrt{t-s}} \quad (1.15) $$

En réalité, l'équation BBO n'est valable que pour $Re_p = d_p|\mathbf{U}_r|/\nu_f \ll 1$. Néanmoins, elle est communément utilisée même lorsque $Re_p \gg 1$ (mais en appliquant des expressions généralisées pour le terme de traînée).

D'autres forces peuvent également être incluses en plus de ces contributions. Un premier exemple est celui des forces de portance, correspondant à des forces agissant dans la direction perpendiculaire à la vitesse relative entre le fluide et la particule. Elles peuvent soit être issues d'un cisaillement de la vitesse fluide à l'échelle de la particule (force de portance dite de Saffman) soit être issues de la rotation de la particule elle-même (force de portance dite de Magnus) [79]. Un second exemple est celui du mouvement Brownien, qui est issu de l'interaction entre les molécules du fluide et une particule. Les effets du mouvement Brownien sont surtout visibles pour des particules plus petites que quelques micromètres (notamment les colloïdes). Enfin, une troisième catégorie de forces concerne les forces électro-magnétiques (par exemple avec des particules chargées dans un plasma) ou les échanges thermiques (lorsque les particules sont en suspension dans un fluide non-isotherme).

## – Interactions particule-particule

Les interactions entre deux particules jouent un rôle dans la formation d'amas de particules, qu'il s'agisse d'agrégats formés dans le fluide (cf. Chapitre 3) ou de tas de particules s'accumulant sur une surface (cf. Chapitre 2).

On distingue typiquement trois types de forces :

* Les forces de contact sont responsables de l'adhésion entre deux particules (composante normale aux surfaces en contact) et des frottements entre ces objets (composante tangentielle de résistance au mouvement). Les forces d'adhésion ont été largement étudiées dans le cadre des théories sur la mécanique du contact [68, 123, 147]. Ces théories prennent notamment en compte la déformation des surfaces lorsqu'elles entrent en contact. Plus précisément, il existe différentes théories selon la dûreté des matériaux, c'est-à-dire leur capacité à se déformer. D'un côté, certaines théories prennent en compte les déformations élastiques : par exemple, les théories JKR (du nom des auteurs Johnson, Kendall et Roberts [71]) et DMT (nommée d'après les travaux de Derjaguin, Muller et Toporov [31]) se fondent sur la théorie de Hertz pour les déformations élastiques. Elles diffèrent par le fait que, en plus des effets de déformation, la théorie DMT prend en compte les forces de non-contact agissant dans la zone de contact (liées aux forces de van der Waals). De l'autre côté, les théories M-P (d'après les travaux de Maugis et Pollock [92]) incluent les déformations plastiques.

Néanmoins, ces théories sur la mécanique du contact ont été construites en supposant les surfaces parfaitement lisses. En réalité, les surface sont rugueuses (surtout à l'échelle nanoscopique) et les forces de contact s'en trouvent profondément modifiées. La prise en compte de l'effet de la rugosité de surface est possible en ayant recours à divers types de modèles statistiques pour décrire une surface rugueuse (cf. [123] pour plus de détails ainsi que le Chapitre 2).
▷ Les forces à courte portée correspondent aux interactions ayant lieu sur des distances plus courtes ou du même ordre de grandeur que la taille caractéristique des particules [38, 66].

C'est le cas des forces inter-surfaciques, qui résultent des interactions microsco-piques entre chacune des molécules composant les objets en interaction. Il s'agit notamment des forces de van der Waals [115], qui sont dues à l'interaction entre des dipôles (permanents ou induits), ou bien des forces de type Lennard-Jones entre les molécules [68]. Ces forces ont typiquement lieu lorsque la distance entre les deux surfaces est de l'ordre de quelques nanomètres ou micromètres. Des forces dites de double-couche électrique (EDL pour "Electrical Double-Layer" en anglais) ont également lieu lorsque des particules chargées sont mises en suspension dans un liquide [68]. Dans ce cas, il se forme près de la surface une couche de liquide avec une sur-concentration d'ions ayant une charge opposée à celle de la surface. Cette couche de fluide est appelée double-couche électrique et s'étend sur une distance de l'ordre de la longueur de Debye (qui dépend de la concentration d'ions dans la solution). Lorsque deux surface chargées se rapprochent, le recouvrement de ces doubles couches donne lieu à des forces attractives (resp. répulsives) lorsque les charges sont opposées (resp. de même signe). Enfin, certaines forces physico-chimiques voir biologiques peuvent avoir lieu dans des cas particuliers (par exemple, les forces dues à la présence de ré-gions hydrophiles et hydrophobes sur les surfaces ou bien des forces de nature biologiques entre des objets organiques) [68].

D'autres forces à courte portée peuvent également avoir lieu à des échelles com-parables à la taille des particules. Par exemple, c'est le cas des forces capillaires qui ont lieu lorsqu'un pont liquide se forme entre deux surfaces en suspension dans l'air. Un autre exemple est celui des forces liées à la présence de polymères dans le fluide, qui peuvent donner lieu à des forces attractives ou répulsives selon les structures polymériques qui se forment entre les surfaces des objets [68].

▷ Les forces à longue portée correspondent aux interactions ayant lieu sur des distances beaucoup plus grandes que la taille caractéristique des particules. C'est le cas des forces électro-magnétiques, qui existent lorsque les particules sont chargées et/ou présentent des propriétés magnétiques.

- Interactions particule-surface

Les interactions entre une particule et une surface ont lieu lorsqu'une particule se rapproche suffisamment d'une paroi/surface pour que ces forces ne soient plus négligeables. La nature de ces forces est la même que les interactions particule-particule. Ainsi, on distingue les mêmes types de forces (c'est-à-dire entre les forces de contact, les forces à courte portée et les forces à longue portée).

- Effets internes à chaque particule

Dans le cas où les particules sont des organismes biologiques (comme les poissons, les oiseaux ou bien le plancton), elles peuvent avoir une motilité propre [87], c'est-à-dire qu'elles se déplacent dans une direction donnée et avec une certaine vitesse. On parle alors de particules actives, par opposition aux particules passives qui n'ont pas de motilité. Pour reproduire ces phénomènes, des modèles simplifiés proposent d'imposer la magnitude de la vitesse liée à la motilité propre de chaque particule, orientée selon l'orientation de l'objet à chaque instant [87].
Nous venons de voir que le choix des théories utilisées pour établir les équations décrivant les phénomènes sous-jacents s'avère une tâche délicate à cause du large éventail des possibilités s'offrant à nous. Néanmoins, une fois les équations établies, il est impératif de garder en mémoire les hypothèses associées à l'établissement de chacune d'elles, puisqu'elles peuvent s'avérer être un facteur limitant l'application d'un modèle à certaines situations.

### 3ème élément : Approche de modélisation

L'étape suivante consiste à définir l'approche de modélisation choisie pour résoudre les équations qui viennent d'être établies précédemment. Pour ce faire, il existe deux grandes approches [119, 124] (voir aussi la figure 1.7) :

<figure><img src="image_9.png" /><figcaption>FIGURE 1.7 – Résumé des types principales approches pour simuler les particules en suspension.</figcaption></figure>

* Les approches par champ (appelées aussi approches Eulériennes) décrivent les valeurs de champs et leur évolution en fonction du temps, et ce dans chacune des zones de l'espace ainsi pré définies.

Pour la phase fluide, les chercheurs ont souvent recours à ces approches par champ pour résoudre les équations de Navier-Stokes. Ces approches donnent une information sur le champ fluide (à savoir le champ de vitesse, le champ de pression, ainsi que le champ d'énergie turbulente lorsqu'il s'agit d'écoulements turbulents). Ces approches reposent donc sur une discrétisation de l'espace pour résoudre ces équations. Elles regroupent notamment les méthodes par volume fini (abrégé FVM pour "Finite Volume Method" en anglais), les méthodes par élément fini (abrégé en FEM pour "Finite Element Method") ou encore les méthodes spectrales.

Pour la phase dispersée (ici des particules solides), les champs ainsi décrits sont généralement la concentration de particules (scalaire), la vitesse moyenne des particules (vecteur) ou bien la taille moyenne des particules (scalaire). Les équations d'évolution de ces champs sont données par des équations de convection-diffusion, qui rendent compte du transport de ces champs via les processus physiques sous-jacents.

* Les approches trajectorielles (appelées aussi approches Lagrangiennes) proposent de suivre la trajectoire de plusieurs objets à mesure qu'ils se déplacent dans le temps et l'espace. La dynamique de chaque élément suivi est donnée par la $2^{nde}$ loi de Newton comme dans l'équation (1.7).

Pour la phase fluide, les approches Lagrangiennes consistent à suivre les trajectoires d'un grand nombre de particules fluides. Les équations pour la dynamique de ces particules fluides sont dérivées des équations de Navier Stokes (considérant donc le bilan de masse, de quantité de mouvement et d'énergie). L'intérêt de ces approches est qu'elles
ne requièrent pas de discrétisation spatiale (de type maillage) et sont donc plus adaptées aux écoulements dans des géométries complexes. Ces approches Lagrangiennes ont eu un essor plus tardif, notamment à cause du coût numérique plus important que les approches Eulériennes. Néanmoins, avec les progrès significatifs sur les algorithmes parallèles ainsi que sur la puissance de calcul, les approches Lagrangiennes sont de plus en plus fréquemment utilisées pour résoudre des écoulements fluides. C'est le cas notamment des approches SPH (pour "Smooth-Particle Hydrodynamics").

Pour la phase dispersée, les approches trajectorielles sont plus intuitives puisqu'elles reproduisent naturellement la dynamique de chaque particule en suspension dans l'écoulement. Comme détaillé dans l'élément suivant, il existe en fait toute une gamme d'approches Lagrangiennes avec des niveaux de description différents.

Il existe aussi des approches alternatives [93], comme par exemple les méthodes de Boltzmann sur réseau (LBM pour "Lattice-Boltzmann Method" en anglais). Les méthodes LBM consistent à résoudre l'équation discrète de Boltzmann pour simuler le comportement de fluides à une échelle mésoscopique. Ces méthodes se fondent sur des schémas de collision-propagation, décrivant les processus de collision entre les éléments ainsi que leur transport/propagation sur le réseau.

À ce stade, il est important de noter que d'autres méthodes numériques sont utilisées spécifiquement pour l'étude des particules de taille finie (c'est-à-dire comparable ou plus grande que l'échelle de longueur de Kolmogorov). Dans ce cas, il est parfois nécessaire de reconstruire l'interface exacte entre la particule et le fluide. Il existe en fait toute une littérature sur la reconstruction d'interface, en ayant recours à des approches comme IBM (pour "Immersed Boundary Method") [124]. Néanmoins, je n'aborderai pas ces méthodes dans le reste du mémoire puisque, à ce jour, je n'ai pas encore réalisé de simulations DNS avec des particules de taille finie.

**4<sup>eme</sup> élément : Niveau de description et précision**

Une fois les équations établies et l'approche de modélisation choisie, l'étape suivante consiste à définir le niveau de description de chacune de ces approches. Il existe en effet différents niveaux de description, allant de la résolution la plus exacte des équations à des résolutions approchées (c'est-à-dire moins précises mais plus rapides). Cette étape résulte ainsi d'un savant compromis entre la précision requise pour les résultats attendus ainsi que la rapidité des calculs effectués (voir aussi la Section 1.2.3 sur les aspects numériques).

* La simulation de l'écoulement fluide turbulent :

Les dernières décennies ont vu l'essor de modèles pour la simulation d'un écoulement fluide turbulent [119]. En effet, un écoulement turbulent est caractérisé par une large gamme d'échelles spatio-temporelles. En conséquence, il est nécessaire de choisir si l'on a besoin d'avoir accès à des données entièrement résolues ou si l'on peut se contenter d'un nombre plus réduit d'informations (et si oui, lequel). En pratique, plusieurs niveaux de description existent en fonction de la précision requise, du niveau d'information et des coûts numériques. Je résume ci-dessous les principaux outils disponibles en fonction du niveau de description, allant du microscopique vers le macroscopique (cf. la figure 1.8 et [119] pour plus de détails) :
<figure><img src="image_10.png" /><figcaption>FIGURE 1.8 – Résumé des niveaux de description pour simuler un écoulement fluide.</figcaption></figure>

- Les simulations numériques directes (ou DNS pour “Direct Numerical Simulation” en anglais) se fondent sur une résolution numérique des équations de Navier-Stokes sans modèle de turbulence. En pratique, cela veut dire que toutes les échelles (spatiales et temporelles) du fluide turbulent sont explicitement résolues. Pour ce faire, les approches se fondent souvent sur des méthodes spectrales (utilisant des décompositions de Fourier et de Chebyshev) et sur des discrétisations spatiales très fines (la distance séparant deux points successifs devant être plus petite que l’échelle spatiale de Kolmogorov). Ainsi, les simulations DNS donnent accès à la vitesse instantanée à n’importe quel point de l’espace considéré. Néanmoins, ces approches ne sont applicables qu’à des géométries simples (turbulence homogène isotrope, canal plan) et engendrent des coûts numériques élevés les rendant encore inadaptées pour traiter des situations réalistes (typiquement avec des domaines à géométrie complexe).

- Les simulations des grandes structures de la turbulence (ou LES pour “Large Eddy Simulation” en anglais) consistent à résoudre exactement les grandes échelles et à filtrer les petites échelles de la turbulence. En conséquence, des modèles dits de sous-maille sont introduits pour rendre compte des effets de la turbulence à une échelle plus petite que celle résolue numériquement. Grâce aux avancées en matière de performance algorithmiques et de puissance de calcul, les LES sont de plus en plus utilisées dans l’étude des écoulements autour de géométries complexes (écoulement autour d’obstacle, d’aile d’avion etc...) puisqu’elles donnent accès à une information assez fine pour des coûts de calcul acceptables.

- Les simulations SPH (pour “Smoothed Particle Hydrodynamics”) correspondent à des approches entièrement Lagrangiennes. Développées initialement dans le domaine de l’astrophysique [105, 93], l’approche consiste à suivre les trajectoires d’un grand nombre de particules fluides et, pour limiter les coûts de calcul, de ne résoudre les interactions entre ces particules fluides que dans le voisinage de chacune des particules (via des fonctions dites de lissage qui pondèrent la contribution de chacun des voisins en fonction de leur distance à la particule considérée).

- Les modèles par moyenne de Reynolds (ou RANS pour “Reynolds-Averaged Navier-Stokes”) proposent de résoudre une équation de Navier-Stokes moyennée. Plus précisément, la décomposition de Reynolds consiste à écrire la vitesse du fluide $U_f$ comme la somme d’une contribution moyenne $\langle U_f \rangle$ et d’une fluctuation $u'_f$. En
injectant dans les équations de Navier-Stokes (1.1), on obtient :

$$ \frac{\partial \langle \mathbf{U}_f \rangle}{\partial t} + \langle \mathbf{U}_f \rangle \cdot \nabla \langle \mathbf{U}_f \rangle = -\nabla \cdot (\langle \mathbf{u}'_f \otimes \mathbf{u}'_f \rangle) - \frac{1}{\rho_f} \nabla \langle P \rangle + \nu_f \nabla^2 \langle \mathbf{U}_f \rangle + \langle \mathbf{F} \rangle. \quad (1.16) $$

L'équation (1.16) fait apparaître un terme non-linéaire additionnel, appelé tenseur de Reynolds (le 1er terme à droite de l'équation). Pour résoudre cette équation, une loi de fermeture supplémentaire est donc nécessaire puisque le tenseur de Reynolds n'est pas connu. Cette loi de fermeture est souvent inspirée des travaux de Boussinesq qui relient la viscosité turbulente au tenseur de Reynolds. Deux méthodes existent dans la littérature : (i) les modèles locaux, où le tenseur de Reynolds est approximé directement avec les valeurs locales connues du fluide (comme les modèles $k - \epsilon$ et $k - \omega$); (ii) les modèles non-locaux, c'est-à-dire qui résolvent des équations supplémentaires pour approximer ce tenseur.

## o La simulation du transport de particules

Il existe là aussi plusieurs niveaux de description pour la simulation de la dynamique de particules [62, 98]. Dans la suite du mémoire, j'aborde notamment les approches trajectorielles (voir Figure 1.9).

<figure><img src="image_11.png" /><figcaption>FIGURE 1.9 – Résumé des niveaux de description pour simuler la dynamique de particules.</figcaption></figure>

- Les approches dites de suivi à N-particules (aussi appelées DEM pour "Discrete Element Method") consistent à suivre en même temps la trajectoire de N-particules en considérant toutes les forces qui s'exercent entre chaque couple de particules. Ces approches ont initialement été développées pour les milieux granulaires, c'est-à-dire les milieux formés par un ensemble de grains se déplaçant de façon conjointe (comme le sable) [26]. Cela signifie que, pour calculer la dynamique de N-particules, il est nécessaire de calculer $O(N^2)$ interactions. En termes de niveau d'information, ces approches sont donc comparables à des DNS puisqu'elles consistent à résoudre entièrement les équations. Ainsi, elles fournissent des informations sur la PDF à N-particules (y compris sur les moments d'ordre N). Par exemple, ces approches ont été utilisées pour obtenir la dynamique de particules et la morphologie exacte des agrégats formés. Néanmoins, tout comme la DNS, elles sont limitées à des situations simplifiées à cause de leur coût de calcul important (typiquement quelques milliers ou dizaine de milliers de particules).

- Des approches à gros grains ont été développées en s'inspirant des approches DPD (pour "Dissipative Particle Dynamics"). Ces approches consistent à suivre numériquement des méso-particules, aussi appelées parcelles, qui regroupent un grand
nombre de particules réelles. L'intérêt de ce genre d'approche est que cela permet de réduire le nombre de particules à suivre dans un écoulement, rendant les simulations dans des conditions réalistes plus abordables numériquement. Néanmoins, cela implique d'adapter les formulations des forces s'exerçant sur ces méso-particules et au sein de ces particules. En effet, l'information sur la position exacte de chacune des particules réelles n'est plus disponible puisque seule l'information sur la dynamique du centre de masse de chaque méso-particule est connue (position, vitesse, etc.). Ainsi, un peu comme dans une simulation LES pour la phase fluide, une partie de l'information est filtrée. Ainsi, si l'on s'intéresse à la formation d'agrégats, il est nécessaire d'avoir un modèle rendant compte de la manière dont les particules réelles présentes dans une même méso-particule interagissent.

- Des approches encore plus réduites consistent à faire un suivi Lagrangien de particules qui interagissent, non pas directement entre elles, mais plutôt avec un champ moyen. Ces approches reposent d'abord sur une évaluation des statistiques liées aux particules dans différentes régions de l'espace via des méthodes de type Monte-Carlo (par exemple pour calculer la vitesse moyenne ou la concentration des particules). La dynamique des particules est ensuite décrite par des modèles utilisant ces informations de type champ moyen. Ces approches sont appelées approches PDF à une particule ("One-particle PDF approach" en anglais). En effet, ce genre d'approche ne contient pas d'informations sur les statistiques à deux points, puisque chaque particule numérique interagit uniquement avec le champ moyen généré par l'ensemble des autres particules. Le recours à une discrétisation de l'espace (ou maillage) pour le calcul des statistiques particulaires rend ces approches naturellement compatibles avec les approches de type RANS pour la résolution de l'écoulement fluide.

À ce stade, il est important de noter que la liste précédente n'est pas exhaustive. L'objectif est simplement d'illustrer certains des choix possibles en termes de niveaux de modélisation. D'autres alternatives existent : par exemple, certaines s'inspirent des méthodes SPH avec la résolution des forces uniquement entre les voisins les plus proches via des fonctions de lissage (comme par exemple les méthodes SDPD pour "Smoothed DPD"). Il est également possible de raffiner certaines des approches décrites précédemment, en simplifiant certaines des expressions des forces. Ainsi, certaines approches proposent de ne résoudre exactement que les forces ayant une variation lente et de modéliser les forces ayant une variation rapide (comme les méthodes de type Langevin dynamics).

## o La simulation des forces inter-surfaciques

Dans le cas des écoulements chargés en particules se pose également la question du niveau de description des forces inter-surfaciques jouant un rôle dans les processus d'adhésion. Néanmoins, la difficulté est que ces forces ont une portée de l'ordre du nano- ou micro-mètre, donc à des distances bien plus faibles que les échelles de l'écoulement fluide.

- Les approches de type Dynamique Moléculaire (MD) consistent à simuler le mouvement de l'ensemble des molécules formant un objet au cours du temps [129]. De la même manière que les approches trajectorielles, la dynamique moléculaire classique résout donc les équations de mouvement de chaque molécule dans le système mais au niveau atomique, c'est-à-dire en calculant les interactions entre chaque molécule
(via des forces potentielles étant décrites par des modèles de type Lennard-Jones par exemple). Leur niveau de description à l'échelle de l'atome/la molécule permet donc d'obtenir un niveau d'information très fin sur les forces inter-surfaciques ou bien sur la déformation de surfaces en contact. Néanmoins, les coûts de calcul de ces méthodes limitent leur application à des situations simplifiées (par exemple l'adhésion entre deux particules micrométriques).

- Tout comme les approches trajectorielles, des approches à gros grains ont été développées dans le cadre des approches DPD, qui proposent de suivre des mésoparticules composées d'un ensemble de molécules.

- Enfin, il existe des formulations macroscopiques de ces forces. Dans l'élément 2, j'ai mentionné l'existence des théories de mécanique du contact (par exemple les théories JKR, DMT ou M-P) [68, 123]. Ces théories proposent justement une formulation macroscopique des forces d'adhésion entre des particules en prenant en compte des informations sur les propriétés des matériaux (dureté, déformabilité, etc.) et sur les formes des surfaces en contact. De la même manière, des formules macroscopiques des forces à courte portée entre deux particules ont été développées en prenant parfois en compte l'effet de la rugosité de surface [62, 61].

### 5ème élément : Cohérence et consistance

L'une des difficultés majeures lors du développement de tels modèles est de veiller à la cohérence et à la consistance en termes de niveau de description et d'hypothèses pour chacun des modèles utilisés. Cet aspect est d'autant plus important que, comme indiqué dans la section 1.1, les écoulements chargés en particules impliquent une multitude d'échelles spatio-temporelles ainsi que plusieurs disciplines (voir aussi la Figure 1.10).

<figure><img src="image_12.png" /><figcaption>FIGURE 1.10 – Illustration des aspects multi-échelles et multi-physiques en jeu.</figcaption></figure>
Par exemple, il est tout à fait possible numériquement de coupler un suivi Lagrangien à N-particules avec une approche de type RANS pour la simulation d'un écoulement turbulent. Néanmoins, il y a une incohérence profonde dans les modèles physiques choisis. En effet, l'approche RANS pour l'écoulement fluide ne donne pas accès à une information sur la vitesse instantanée à la position d'une particule, mais à une vitesse moyenne et à des fluctuations turbulentes. En conséquence, le niveau de description de l'approche globale RANS+DEM se retrouve limité à celui de l'approche RANS. Pire, les modèles RANS ne donnent qu'une information en un point, les rendant ainsi incompatibles avec les suivis Lagrangiens à N-points puisque le niveau de description n'est pas le même. Ainsi, par souci de consistance, les approches RANS sont compatibles avec des suivis Lagrangiens à 1 point. À l'inverse, si l'on tient absolument à avoir recours à une approche DEM pour la phase dispersée, il est nécessaire d'avoir une information à N-points sur le fluide : cela implique, pour des particules colloïdales, l'usage de DNS alors que, pour des particules plus larges, des approches de type LES peuvent être suffisantes (à condition que la méthode résolve explicitement les échelles de la turbulence auxquelles les particules sont sensibles).

Ces questions de cohérence et de consistance ne sont pas toujours aussi évidentes. Par exemple, si l'enjeu de l'étude est d'évaluer finement la formation de dépôts multicouches et son effet sur l'écoulement turbulent, il semble approprié de coupler des approches DNS pour la résolution du fluide turbulent et des approches DEM. Néanmoins, si les particules ont une taille plus grande que la plus petite échelle de l'écoulement turbulent, il est plus pertinent d'avoir recours à des simulations DNS avec des particules de taille finie PR-DNS (c'est-à-dire en résolvant l'écoulement autour des particules). A l'inverse, lorsque les particules sont beaucoup plus petites que l'échelle de Kolmogorov, l'approximation de particules ponctuelles est utilisée PP-DNS. De la même manière, toujours dans le cas d'études de dépôts multicouches, les simulations les plus fines ne correspondent pas forcément à des simulations DNS (avec ou sans taille finie) couplées à une approche DEM. En effet, la morphologie du dépôt formé sur les parois résulte non seulement du transport de particules par l'écoulement turbulent mais aussi des interactions particule-paroi et particule-particule. Dans la communauté diphasique, on trouve deux types de simulations microscopiques de tels dépôts : (i) celles qui supposent que les particules se déposent exactement à l'endroit où elles entrent en contact pour la première fois avec le dépôt déjà formé ; (ii) celles qui supposent que, une fois la particule ayant collisionné avec le dépôt déjà formé, cette dernière va rouler jusqu'à atteindre un équilibre (c'est-à-dire être en contact avec 3 autres particules). D'un point de vue des interactions inter-surfaciques, il apparaît bien que ces deux formulations correspondent à des niveaux de description macroscopiques et qu'elles constituent les deux extrêmes du phénomène d'adhésion. En réalité, une simulation microscopique à l'échelle des agrégats consisterait à coupler ces simulations à des simulations de type MD pour le calcul des forces inter-surfaciques. Cela permettrait ainsi de rendre compte de toute une gamme de phénomènes intermédiaires, où les particules peuvent rouler partiellement et s'arrêter dans des positions d'équilibre autre que ces deux extrêmes à cause de la présence de rugosité sur les surfaces des particules (c'est-à-dire de non-homogénéités dans la topologie des surfaces).

Ce dernier exemple révèle aussi l'importance que revêt la consistance entre les descriptions numériques lorsqu'il s'agit d'inclure des phénomènes ayant lieu à des échelles différentes. Ainsi, l'étape de transport de particules colloïdales dans des tuyaux peut s'effectuer sur
des distances allant jusqu'à plusieurs mètres alors que l'étape d'attachement de ces par-
ticules à une paroi a lieu sur des distances typiquement inférieures à quelques dizaines
de nanomètres. Dans ce cas, il apparaît inadéquat de coupler un calcul de type RANS
avec un suivi Lagrangien à 1 particule et avec une simulation de type MD pour les forces
inter-surfaciques. En effet, cela induirait une inconsistency sur le niveau de description
des phénomènes en proche paroi. D'un côté, le mouvement des particules transportés par
le fluide proche d'une paroi serait approximé (avec des erreurs à l'échelle du micromètre
voire du millimètre) alors que les forces physico-chimiques à plus petites échelles seraient
calculées de manière exacte (avec des erreurs à l'échelle du nanomètre). Ainsi, le facteur li-
mitant pour l'estimation des erreurs de calculs serait lié à l'approximation du transport en
proche paroi, et ce malgré une simulation exacte des forces inter-surfaciques. Néanmoins,
cet exemple montre qu'il est tout à fait possible de recourir à des approches multi-échelles
utilisant ce genre de modèle. Par exemple, il serait consistant d'utiliser une approche de
type MD pour estimer les valeurs moyennes des forces entre deux surfaces (voire des dis-
tributions de force), puis d'utiliser cette information filtrée dans une approche de type
RANS couplée à un suivi Lagrangien à une particule.

À travers ces différents éléments, il apparaît clairement que le choix d'une approche de modé-
lisation est une tâche ardue qui ne doit pas être sous-estimée. Elle exige en fait une compréhension
approfondie des enjeux, des objectifs spécifiques, ainsi que des limites des fondements/théories
actuels. Par ailleurs, d'autres éléments peuvent entrer en ligne de compte dans le développement
d'un modèle. Par exemple, certains modèles se veulent génériques et/ou flexibles, dans le sens
où ces modèles peuvent aisément être enrichis avec des notions initialement non abordées (par
exemple avec des aspects biologiques dans le cas de particules organiques).

Il apparaît aussi que les éléments exposés précédemment ne sont pas indépendants les uns
des autres. Il existe en fait une véritable inter-dépendance entre ces différents éléments. Par
exemple, il est évident que le choix d'une approche de type MD pour la résolution des forces inter-
surfaciques (relevant de l'élément 4) nécessite d'établir les équations de mouvement au niveau
moléculaire et non pas au niveau d'une particule prise comme un objet indivisible (notion relevant
de l'élément 2). De la même manière, un agrégat de particules peut être décrit numériquement
de deux manières : soit comme un ensemble de particules élémentaires collées les unes aux
autres (approche microscopique) soit comme une macro-particule ayant des propriétés spécifiques
(e.g., sa porosité, sa dimension fractale, son nombre de particules individuelle dans l'agrégat,
etc.). La première approche nécessite de décrire les interactions entre chaque particule prise
individuellement alors que la seconde nécessite de formuler des modèles pour l'interaction entre
des agrégats dont la morphologie n'est caractérisée que par des informations macroscopiques. Ces
inter-relations et inter-dépendances entre les éléments illustrent bien le fait que le développement
d'un modèle résulte d'un compromis entre ces différentes notions. Ainsi, l'ordre des éléments ne
doit pas être pris à la lettre, d'autant plus que j'ai choisi ici un certain ordre pour des raisons
purement narratives. Au contraire, cette inter-dépendance des éléments incite plutôt à cultiver
et entretenir ces allers-retours de pensée entre les différentes notions.

1.2.3 Implémentation numérique

La seconde étape dans la modélisation consiste à implémenter le modèle physique, c'est-à-dire
l'écriture d'un programme informatique (que j'appelle souvent code par la suite). En pratique,
cela a conduit au développement de logiciels de Mécanique des Fluides Numérique (CFD pour
"Computational Fluid Dynamics" en anglais). Ces logiciels permettent de simuler l'écoulement
fluide et de calculer la dynamique de particules en suspension. Ces logiciels CFD sont de plus en plus fréquemment utilisés par les chercheurs et les ingénieurs afin de réaliser des simulations numériques de problèmes donnés.

Le choix des algorithmes retenus pour implémenter un modèle physique et la façon de les programmer dépendent de deux éléments supplémentaires : les coûts numériques acceptables et la précision requise. En fait, ces contraintes sont intimement liées à l'élément 4 sur le niveau de description ainsi qu'au contexte de l'étude (élément 1). Le choix d'une approche numérique avec un certain niveau de description résulte toujours d'un subtil compromis entre la précision requise et le coût numérique admissible. Ainsi, les approches numériques avec un niveau de description très fin (comme DNS & DEM) sont préférables pour aborder des phénomènes locaux (comme par exemple la morphologie locale d'un agrégat). A l'inverse, les approches macroscopiques sont souvent préférées pour aborder des questions d'ingénierie, qui nécessitent d'obtenir rapidement une information à l'échelle macroscopique (par exemple la quantité de dépôt dans un échangeur de chaleur).

**Choix d'un algorithme** : La résolution des équations d'un modèle pose des questions supplémentaires liées au choix de l'algorithme. Les équations sont en effet résolues avec une discrétisation temporelle et spatiale (au sens de maillage dans les descriptions par champ ou de particules dans les descriptions trajectorielles). Cela aurait pu faire l'objet d'un critère supplémentaire dans la liste décrite précédemment, mais j'ai choisi ici de marquer la différence entre les aspects plus théoriques liés à l'établissement du modèle et les aspects plus algorithmiques liés à leur implémentation numérique (la frontière entre les deux étant poreuse).

Il existe en fait toute une littérature spécifique sur les méthodes numériques [55, 67]. Les schémas numériques sont classés en fonction de plusieurs catégories. On distingue notamment les méthodes implicites, les méthodes explicites et les méthodes hybrides (semi-implicite, etc.). Les schémas explicites évaluent l'état d'un système à un moment ultérieur uniquement à partir de l'état du système à l'instant présent alors que les schémas implicites prennent en compte à la fois l'état présent et l'état ultérieur. Le choix de l'une ou l'autre de ces méthodes est lié à un compromis entre la stabilité du calcul et son coût (les méthodes implicites étant plus stables mais lentes). Les schémas numériques peuvent aussi être répertoriés en fonction de leur précision, par exemple entre les schémas du 1<sup>er</sup> ordre, du 2<sup>nd</sup> ordre ou d'ordre plus élevé. Une autre catégorisation consiste à classer les schémas numériques en fonction du type d'équation à résoudre. C'est cette classification que j'ai choisi de rappeler brièvement ci-dessous puisqu'elle permet de faire la distinction entre :

* Les schémas numériques pour les équations aux dérivées ordinaires (EDO) :

Les EDO correspondent aux équations différentielles dont la ou les fonction(s) inconnue(s) ne dépendent que d'une seule variable. Par exemple, pour un problème de condition initiale sur un scalaire, l'EDO s'écrit sous la forme suivante :

$$y'(t) = f(t, y(t)), \quad y(t_0) = y_0 \tag{1.17}$$

où $f$ est une fonction dans $f : [t_0, +\infty[ \times \mathbb{R} \to \mathbb{R}$ et $y_0$ un scalaire donné.

Pour résoudre une telle EDO, il existe plusieurs types de schémas dont (liste non-exhaustive) :

- Les schémas du 1<sup>er</sup> ordre avec notamment les schémas d'Euler explicite (aussi appelée schéma d'Euler "forward") ou implicite (ou schéma d'Euler "backward"). Le schéma d'Euler explicite repose sur l'approximation de la dérivée par le développement limité
au 1er ordre (h étant le pas de l'intégration entre t_n et t_{n+1}) :

$$y_{n+1} \approx y_n + h f(t_n, y_n) \qquad (1.18)$$

Le schéma d'Euler implicite repose sur le même développement limité mais "backward" en temps :

$$y_{n+1} \approx y_n + h f(t_{n+1}, y_{n+1}) \qquad (1.19)$$

- Les schémas itératifs d'ordre plus élevé, comme les schémas de Runge-Kutta, qui proposent d'introduire des valeurs intermédiaires pour calculer la valeur y(t_{n+1}). Ainsi, le schéma de Runge-Kutta d'ordre 2 (RK2) consiste à estimer la dérivée au milieu du pas de temps d'intégration puis à refaire le pas d'intégration complet à partir de cette estimation :

$$y_{n+1} \approx y_n + h f\left(t_n + \frac{h}{2}, y_{n+\frac{1}{2}}\right) \qquad (1.20)$$
avec $y_{n+\frac{1}{2}} = y_n + \frac{h}{2} f(t_n, y_n)$.

Fréquemment utilisée, le schéma de Runge-Kutta d'ordre 4 (RK4) propose d'estimer la pente de la fonction via une moyenne pondérée des pentes au début, au milieu et à la fin de l'intervalle considéré :

$$\begin{align*}
y_{n+1} &\approx y_n + \frac{h}{6} (k_1 + 2k_2 + 2k_3 + k_4) \\
\text{où } &k_1 = f(t_n, y_n) \\
&k_2 = f(t_n + \frac{h}{2}, y_n + \frac{h}{2}k_1) \tag{1.21} \\
&k_3 = f(t_n + \frac{h}{2}, y_n + \frac{h}{2}k_2) \\
&k_4 = f(t_n + h, y_n + h k_3)
\end{align*}$$

- Les schémas à pas multiples qui, contrairement aux schémas précédents dits à pas unique (puisque les informations intermédiaires des schémas RK ne sont pas sauvegardées), sauvegardent et utilisent les informations sur plusieurs intervalles de temps. C'est le cas par exemple des schémas d'Adam-Bashforth qui s'écrivent à l'ordre 2 :

$$y_{n+2} \approx y_{n+1} + \frac{3}{2} h f(t_{n+1}, y_{n+1}) - \frac{1}{2} h f(t_n, y_n). \qquad (1.22)$$

À noter qu'il existe des formulations généralisées qui considèrent des combinaisons linéaires à ordre quelconque.

- Les schémas linéaires généralisés, qui proposent d'avoir recours à la fois à la notion de pas multiples et aux schémas itératifs (comme Runge-Kutta). Ce genre de schémas inclut par exemple les schémas dits de prédiction-correction, qui consistent à faire d'abord une première évaluation de la valeur à l'instant suivant (souvent via un schéma explicite) puis à corriger cette valeur (avec un schéma implicite). Par exemple, ces schémas de prédiction-correction sont utilisés dans les algorithmes de résolution des équations de Navier-Stokes de type SIMPLEC (le champ de vitesse est prédit en utilisant les informations sur le gradient de vitesse et pression au temps t puis corrigé en imposant une correction de pression imposant la divergence nulle).
- Les schémas exponentiels qui proposent d'intégrer de manière exacte les termes linéaires des EDO. Ainsi, en décomposant l'EDO sous la forme d'un terme linéaire $\mathcal{L}$ et d'un terme non-linéaire $\mathcal{N}$ :

$$y'(t) = \mathcal{L}y(t) + \mathcal{N}(y(t)), \quad y(t_0) = y_0, \qquad (1.23)$$

on peut écrire la solution exacte sous la forme :

$$y(t) = y_0 e^{\mathcal{L}t} + \int_{t_0}^{t} e^{\mathcal{L}(t-s)} \mathcal{N}(y(s))ds. \qquad (1.24)$$

L'intégration du $2^{nd}$ membre de l'équation peut se faire avec l'un des schémas numériques présenté précédemment (par exemple avec une méthode explicite où l'on considère comme constant la fonction $\mathcal{N}(y(s))$ au cours du pas de temps).

* Les méthodes numériques pour les équations aux dérivées partielles (EDP) :

Les EDP correspondent aux équations différentielles dont la ou le(s) fonction(s) inconue(s) dépendent(s) de plusieurs variables. Dans le cadre de mes travaux, les équations de dynamique des fluides et des particules en suspension font justement intervenir à la fois le temps et l'espace. Comme déjà mentionné dans le critère 2, il existe plusieurs types de schémas numériques pour résoudre ce type d'équation, avec notamment :

- les méthodes par différences finies (FDM), qui consistent à approximer les dérivées par des différences finies. C'est le cas par exemple du schéma de Cranck-Nicolson, qui correspond à un schéma implicite du $2^{nd}$ ordre en temps faisant une pondération entre le schéma d'Euler "backward" et "forward".

- les méthodes par volumes finis (FVM), qui se fondent sur des approximations d'intégrales directement. Plus précisément, dans les méthodes FVM, l'EDP est résolue sur une discrétisation spatiale (appelée maillage), qui divise le domaine en un ensemble de volumes élémentaires finis. Les EDPs sont ensuite discrétisées et approximées en évaluant les terme de flux aux interfaces entre les volumes finis. C'est le cas par exemple des schémas de Godunov ou des solveurs de Riemann.

- les méthodes par éléments finis (FEM), qui sont similaires aux méthodes FVM puisqu'elles reposent sur la discrétisation spatiale du domaine en éléments de taille finie. La différence par rapport aux FVM est que ces méthodes calculent ensuite une solution de l'EDP à l'aide des principes hérités de la formulation variationnelle. C'est le cas des méthodes spectrales qui reposent sur l'utilisation de polynômes orthogonaux (e.g., les polynômes de Legendre ou de Tchebychev) comme fonctions de base.

Il existe aussi des schémas hybrides, comme les méthodes de Galerkine discontinue qui réunissent des propriétés de la méthode des éléments finis (approximation polynomiale de la solution par cellule) et de la méthode des volumes finis (définition locale de l'approximation et calcul des flux aux interfaces des cellules du maillage).

* Les méthodes numériques pour les équations aux dérivées stochastiques (EDS) :

Les EDS correspondent à des équations différentielles qui comportent un terme de bruit blanc (c'est-à-dire la réalisation d'un processus aléatoire). Dans le domaine des écoulements multiphasiques, les EDS sont fréquemment utilisées pour les processus de diffusion ou pour modéliser des particules ayant des trajectoires aléatoires [46]. Dans le cas multidimensionnel (les résultats sont facilement transposables au cas multidimensionnel), on
cherche à résoudre l'EDS qui s'écrit avec un bruit blanc Gaussien $\xi(t)$ :

$$ \frac{dy(t)}{dt} = \mu(y(t), t) + \sigma(y(t), t) \xi(t) \qquad (1.25) $$

où $\mu(y(t), t)$ et $\sigma(y(t), t)$ sont deux fonctions.

Les schémas numériques pour résoudre ces EDS sont similaires aux schémas pour les EDO, à l'exception près qu'ils doivent respecter les contraintes du calcul stochastique (voir [100, 98, 46]). L'objectif est de donner du sens à la version intégrée de l'EDS :

$$ y(t) = y(t_0) + \int_{t_0}^{t} \mu(y(s), s)ds + \int_{t_0}^{t} \sigma(y(s), s)dW(s) \qquad (1.26) $$

où le 3ème membre de droite correspond à une intégrale stochastique faisant apparaître l'incrément $dW(t) = \xi(t)dt$ d'un processus de Wiener $W$. Pierre angulaire des EDS, le processus de Wiener $W$ est défini comme un processus Gaussien, c'est-à-dire entièrement caractérisé par sa moyenne et sa variance. Le processus de Wiener $W$ possède les propriétés suivantes : les trajectoires sont continues mais différentiables nulle part, ses incréments $dW_t$ sont stationnaires et indépendants (avec une moyenne nulle et une variance égale au pas de temps $dt$). Si l'on revient à la définition de l'intégrale stochastique, on peut l'écrire en divisant l'intervalle de temps en un certains nombres de petits pas de temps $\tau_i = [t_i; t_{i+1}]$ puis en prenant la limite

$$ \int_{t_0}^{t} \sigma(y(s), s)dW(s) = \lim_{N \to \infty} \sum_{i=0}^{N} \mu(y(s_i), s_i) (W_{t_{i+1}} - W_{t_i}) \qquad (1.27) $$

La difficulté supplémentaire est que cette intégrale n'est pas indépendante du choix des temps intermédiaires $\tau_i$ et que ce choix doit être précisé a priori. Une intégrale stochastique est définie soit avec la prescription d'Ito, qui suppose $\tau_i = t_i$, soit avec la prescription de Stratonovich, faisant l'hypothèse $\tau_i = 0.5 (t_{i+1} + t_i)$.

Avec cette règle du calcul stochastique en mémoire, il est ainsi possible d'écrire des schémas numériques pour l'EDS (1.25). Par exemple, avec un schéma d'Euler explicite, on a :

$$ y(t) = y(t_0) + h \mu(y(t_0), t_0) + \sigma(y(t_0), t_0) \Delta W_{t_0} \qquad (1.28) $$

où le terme aléatoire $\Delta W_{t_0}$ est simplement donné par $\sqrt{h} \zeta$ ($\zeta$ étant une valeur aléatoire tirée dans une loi Gaussienne).

À travers cet aperçu rapide de l'implémentation numérique pour résoudre des équations, il apparaît clairement que les choix algorithmiques résultent d'un subtil compromis entre la stabilité, la précision, la convergence et le coût/l'efficacité des algorithmes.

**Programmation informatique** : En lien avec l'efficacité et la rapidité d'exécution des algorithmes, ces algorithmes peuvent être parallélisés. Le principe du parallélisme est de réaliser plusieurs opérations de manière simultanée, donc de réaliser le plus grand nombre d'opérations en un minimum de temps. Cette parallélisation peut être faite soit en utilisant les ressources disponibles sur un ou plusieurs processeurs CPU (pour "Central Processing Unit" en anglais) soit sur des processeurs graphiques (GPU pour "Graphics Processing Unit"). Les progrès informatiques ont permis l'essor des calculs massivement parallèles, qui reposent sur les super-calculateurs (des machines regroupant un grand nombre de processeurs et donc une grande puissance de calcul). L'intérêt de ces calculs parallèles est donc de résoudre le même problème dans un temps plus court que si la résolution se fait de manière séquentielle (c'est-à-dire sans parallélisation).
### 1.2.4 Simulations et analyse des résultats

La troisième et dernière étape dans la modélisation consiste à réaliser des simulations numériques, c'est-à-dire de lancer le programme informatique (qui va résoudre les équations établies dans le modèle). Les résultats fournis par cette simulation numérique doivent ensuite être analysés. Cela nécessite donc de paramétrer les simulations numériques avec un ensemble de données d'entrée (tirée d'expériences, de simulations précédentes ou de bases de données) et ensuite d'analyser les données de sortie fournies par la simulation numérique.

**Paramétrisation des simulations numériques :** Pour réaliser une simulation numérique d'un phénomène physique, il est nécessaire de renseigner un certain nombre de paramètres d'entrée. Afin d'illustrer cette paramétrisation, je considère ici le cas particulier de la dispersion turbulente de particules illustré sur la Figure 1.11. La simulation consiste en un couplage d'une approche RANS pour la phase fluide et d'une approche à 1-particule (modèle Lagrangien stochastique) pour la phase dispersée. Les paramètres d'entrée peuvent être classés en fonction de leur nature :

<figure><img src="image_13.png" /><figcaption>FIGURE 1.11 – Illustration des paramètres d'entrée pour une simulation de dispersion turbulente.</figcaption></figure>

* Les paramètres physiques :

Dans cet exemple, il s'agit de fournir des informations sur la géométrie du cas étudié (la forme du domaine, les surfaces de bord ainsi que les conditions de bord), sur les propriétés du fluide (densité, viscosité, vitesse/débit) et sur les propriétés des particules (forme, taille, densité, zone d'injection, débit).

* Les paramètres de modélisation :

Les paramètres de modélisation se réfèrent aux choix d'une approche de modélisation (voir section 1.2.2) et d'un modèle numérique (voir section 1.2.3). En fait, il n'est pas rare que les logiciels CFD contiennent plusieurs modèles numériques et que l'utilisateur a donc la possibilité de choisir parmi ces modèles disponibles. Dans cet exemple, l'approche choisie est le couplage entre une approche RANS et un modèle Lagrangien stochastique.

* Les paramètres numériques :

Les paramètres numériques correspondent principalement aux informations sur la discrétisation spatiale (maillage) et temporelle (pas de temps). Ces paramètres ne sont pas considérés comme des paramètres physiques puisqu'ils dépendent intrinsèquement de la méthode numérique choisie. Dans le cas présent, la résolution de l'écoulement fluide par une approche Eulérienne nécessite un maillage, qui correspond à une discrétisation possible de la géométrie globale (qui est elle un paramètre physique).
**Analyse des résultats** Une fois la simulation numérique réalisée, les résultats numériques fournis par cette simulation sont analysés. Comme illustré sur la figure 1.12 dans le cas de la dispersion turbulente, il s'agit d'analyser les données de sortie de ces simulations numériques. Dans le cas d'une dispersion turbulente de particules dans un tuyau, on peut récolter des informations sur l'écoulement (comme le champ de vitesse, champ de pression, etc.) et des informations sur les particules (e.g., concentration dans certaines régions, flux de sortie ou flux déposé).

<figure><img src="image_14.png" /><figcaption>FIGURE 1.12 – Illustration des paramètres de sortie (sur le fluide et sur les particules) pour une simulation de dispersion turbulente.</figcaption></figure>

L'analyse des simulations numériques peut prendre plusieurs formes, incluant notamment :

* La convergence des résultats :

Les simulations étant réalisées avec une certaine discrétisation spatiale (de type maillage ou particules) et une discrétisation temporelle (pas de temps), il est nécessaire de vérifier que les résultats obtenus sont convergés. Par exemple, la convergence en maillage indique que les résultats obtenus deviennent indépendants de la taille du maillage, à condition que la taille des mailles soit plus petite qu'une certaine valeur. De la même manière, la convergence en temps peut être étudiée. Il est important de noter que les résultats convergent plus ou moins vite avec le pas de temps en fonction de la stabilité des modèles numériques choisis.

* La validation des résultats :

Les résultats obtenus sont validés en comparant ces derniers à des résultats théoriques, analytiques et/ou expérimentaux. Cela permet d'évaluer l'erreur commise par le modèle et son implémentation numérique (et parfois l'évolution de cette erreur en fonction de la variation des paramètres d'entrée).

* L'analyse de sensibilité et la quantification d'incertitudes

<figure><img src="image_15.png" /><figcaption>FIGURE 1.13 – Illustration des principes d'analyse de sensibilité & quantification d'incertitudes.</figcaption></figure>

L'analyse de sensibilité (SA) correspond à l'étude de la façon dont l'incertitude sur des données d'entrée d'un code se propage et donne lieu à une incertitude de la sortie (voir aussi la figure 1.13). Il s'agit alors d'estimer des indices de sensibilité qui quantifient
l'influence d'une entrée ou d'un groupe d'entrées sur chacune des sorties. La quantification des incertitudes (UQ) consiste à déterminer la probabilité de certains résultats si certains aspects du système considéré ne sont pas exactement connus.

Lors de récents travaux financés par une mission Covid-Inria, j'ai eu l'occasion de mettre en place une méthodologie reposant sur l'utilisation de ces outils d'analyse de sensibilité et de quantification d'incertitudes pour analyser les résultats de simulations sur la dispersion de gouttelettes dans l'air. Bien qu'encore novice sur l'utilisation de tous ces outils, cette expérience m'a convaincu de l'utilité de ces outils à la fois pour identifier les données d'entrée les plus importantes (ce qui permet de focaliser ensuite les efforts d'amélioration sur ces points) et pour aider au développement de modèles macroscopiques (inspirés de ce qui peut être obtenu via des méta-modèles).

### 1.2.5 Le cycle modèle / implémentation / analyse

J'ai choisi de présenter ici l'art de la modélisation sous la forme d'un triptyque modèle physique, implémentation numérique et analyse. Ce triptyque, bien que très pratique, ne doit pas occulter le fait qu'un véritable dialogue s'opère entre ces trois étapes (voir Figure 1.14). Ainsi, plutôt que de parler d'un triptyque, je devrais plutôt employer le mot "cycle de vie" de la modélisation. En effet, la comparaison entre les résultats numériques fournis par le modèle et les résultats théoriques/expérimentaux peut soulever de nouvelles questions théoriques et donc induire une nouvelle modélisation lorsqu'il n'y a pas un accord parfait entre les deux. C'est justement cet échange entre les différentes étapes qui participe à l'évolution des modèles, qui ne sont donc pas des notions figées mais plutôt des entités vivantes en constante évolution.

<figure><img src="image_16.png" /><figcaption>FIGURE 1.14 – Illustration du cycle de vie des modèles scientifiques : les modèles évoluent constamment au fur et à mesure qu'ils sont testés dans de nouvelles configurations/applications et, en de rares occasions, un changement soudain de paradigme peut conduire au développement d'un tout nouveau domaine/catégorie de modèles.</figcaption></figure>

Par exemple, dans la communauté des écoulements chargés en particules, les simulations numériques ont commencé d'abord par prendre en compte le transport des particules par un écoulement (potentiellement turbulent) en simplifiant voire en négligeant les interactions particulières. Dans le cas des approches hybrides Euler/Lagrange, les particules atteignant la paroi étaient considérées soit comme détruites (informatiquement parlant) soit comme déposées irréversiblement à la surface [62]. Cela a permis d'obtenir des informations pratiques comme la quantité de particules déposées dans un tuyau en fonction du temps. Plus tard, de nouvelles mesures expérimentales ont fourni des informations plus précises sur la morphologie des dépôts.
multicouches formés, ainsi que sur leur effet sur l'écoulement (pouvant aller jusqu'au colmatage complet de la section de passage). Pour reproduire ces phénomènes, les questions liées à la modélisation des interactions particule-paroi sont très vite apparues puisqu'elles jouent un rôle dans la morphologie des dépôts formés. Les modèles ont donc été enrichis avec de nouvelles descriptions prenant en compte ces interactions particule-paroi (et de la même manière particuleparticule). Plus récemment, les questions liées à l'effet de la rugosité des surfaces sur les forces inter-surfaciques et à leur modélisation ont fait leur apparition dans cette communauté (c'était l'objet de mes travaux de thèse).

Un point important est que cette évolution des modèles ne se manifeste pas nécessairement
par un progrès scientifique constant. Ce point a été joliment analysé par Arthur Koestler dans son
livre *Les Somnanbules : Une histoire de la vision changeante de l'homme sur l'univers* [75], où
il a montré que les conceptions scientifiques décrivant la physique de l'Univers ont alterné entre
des phases successives de progrès et parfois de régression. Il a également souligné que les progrès
étaient loin d'être continus, avec des changements soudains dans les théories précédemment bien
acceptées (par exemple, avec la théorie de la gravitation de Newton). Ceci est illustré dans la
Figure 1.14, où le cycle entre modèle/implémentation/analyse peut parfois être rompu après
validation pour conduire à un changement de paradigme. Pour rejoindre ce que j'ai dit dans la
préface sur les connaissances scientifiques, ces nouveaux paradigmes correspondent à de nouvelles
manières de penser les faits connus (comme l'a dit Lawrence Bragg).

Cette évolution des modèles, de leur implémentation et de leur validation est liée à la grande
variété des phénomènes en jeu et à la gamme d'échelles spatiales et temporelles impliquées. Cela
donne lieu à des modèles encore incomplets, qui ne rendent compte que de la partie actuellement
connue de ces phénomènes. Cette incomplétude ne doit cependant pas être vue comme une imper-
fection / défaillance des modèles mais plutôt comme une incitation à poursuivre des travaux de
recherche permettant de repousser ces limites et d'améliorer la compréhension et modélisation de
ces phénomènes. La modélisation ouvre ainsi des questions plus méta-physiques (qui se reflètent
dans les réflexions que j'ai introduites dans la préface). Par exemple, Karl Popper a introduit en
1934 la notion de falsifiabilité : un modèle est falsifiable (ou réfutable) s'il peut être logiquement
contredit par un test empirique potentiellement exécutable avec les technologies existantes. C'est
devenu l'un des piliers de la modélisation scientifique, où les modèles sont considérés comme
valides tant qu'ils n'ont pas été réfutés. Cela souligne encore le fait que la modélisation implique
également des questions éthiques et métaphysiques, donnant ainsi plus de sens à la terminologie
de "l'art de la modélisation".
# Chapitre 2

## Suspensions en milieux confinés : resuspension et phénomènes en proche paroi

### 2.1 Introduction

**Contexte** : Ce chapitre est consacré à mes travaux sur le développement de modèles pour les suspensions de particules sphériques dans des écoulements turbulents en milieux confinés, c'est-à-dire en présence de paroi. En particulier, je me focalise sur les phénomènes ayant lieu en proche paroi, là où le fluide est fortement anisotrope et là où les particules interagissent avec la paroi pour s'y déposer ou être remises en suspension.

Ces situations se rencontrent très fréquemment dans la vie de tous les jours. Un exemple bien connu est celui du sable, de la poussière ou des feuilles dans l'atmosphère. Ces particules ont tendance à sédimenter puis à se déposer sur le sol. Elles peuvent alors être remises en suspension lorsque l'écoulement d'air local est suffisamment intense (que ce soit en raison du vent, d'une personne marchant à proximité ou même d'une porte que l'on ouvre/ferme) [76]. Un autre exemple notable est lié à la qualité de l'eau [84] : celle-ci peut être compromise par la présence de minuscules éléments qui s'accumulent à l'intérieur des tuyaux des systèmes de distribution d'eau (comme les bactéries, les métaux ou même le micro-plastique). Ces phénomènes se produisent aussi dans bon nombre d'applications environnementales. Par exemple, c'est le cas en hydrologie (avec le transport de sédiments par les rivières) [113] ou en océanographie (avec les sédiments charriés par les rivières ou par les vagues) [15]. Enfin, de nombreux processus industriels sont aussi concernés par ces phénomènes. Pour ne citer que quelques exemples, ces phénomènes sont jugés indésirables lorsque les particules sont potentiellement dangereuses pour la santé comme dans les applications médicales (e.g., les particules contaminées provenant des sols d'hôpitaux qui peuvent être resuspendues en raison d'activités humaines ou de systèmes de ventilation [17]), dans l'industrie énergétique (avec les particules radioactives emportées par le vent à la suite d'incidents nucléaires [73]) ou encore dans l'industrie automobile (avec les particules de suie issues des pots d'échappement ou de l'abrasion des pneus ou des freins [6]). À l'inverse, dans d'autres situations, l'accumulation de particules sur des surfaces est un phénomène souhaité, notamment dans les systèmes de filtration [54] qui sont utilisés dans les usines de traitement des eaux usées ou dans l'industrie alimentaire/viticole.
**Verrous** : L'un des enjeux est donc d'élaborer des modèles capables de reproduire ces phénomènes tout en restant applicables dans des situations réalistes (afin d'aider à la conception de nouvelles technologies). Néanmoins, comme indiqué dans l'introduction (voir Chapitre 1), la dynamique des particules dans des écoulements turbulents pariétaux est régie par plusieurs forces : les forces hydrodynamiques (agissant sur le transport de particules par l'écoulement), et les forces physico-chimiques à courte portée responsables de l'adhésion particule/paroi ou particule/particule. La difficulté est ainsi liée au caractère pluri-disciplinaire et multi-échelle de ces phénomènes. Cela rend les simulations directes quasi-impossibles d'un point de vue numérique, puisqu'il faudrait par exemple suivre des milliards de particules en suspension dans un fluide s'écoulant dans une géométrie complexe. L'enjeu est donc de réduire le nombre de degrés de liberté à résoudre numériquement jusqu'à obtenir une description simplifiée qui rend compte des phénomènes les plus importants et reste applicable dans des situations pratiques à grande échelle.

**Objectifs** : L'objectif de ce chapitre est donc d'illustrer l'art de la modélisation (cf. Chapitre 1) dans le cas de modèles macroscopiques pour ces écoulements chargés en particules. Pour ce faire, les modèles (ainsi que leur implémentation numérique dans un code de calcul) ont été développés en respectant un certain nombre de contraintes, parmi lesquelles :

i. les phénomènes importants doivent être traités avec le minimum d'approximations ;

ii. les informations présentes dans le modèle doivent être suffisantes pour traiter correctement certains aspects physiques clés (comme la dispersion en taille des particules ou la formation de dépôts de forme complexe) ;

iii. l'implémentation numérique doit permettre de réaliser des simulations dans des cas complexes et réalistes (e.g. la dispersion d'aérosols dans l'atmosphère).

iv. le modèle doit être flexible, c'est-à-dire qu'il doit pouvoir être facilement enrichi avec de nouvelles formulations/approches par la suite (par exemple avec des aspects biologiques pour simuler la déformation et le dépôt d'organismes vivants).

**Hypothèses** : Dans ce chapitre, je décris mes travaux réalisés en considérant le cas simplifié de particules sphériques dans un écoulement turbulent. Je reviendrai sur le cas de particules non-sphériques ultérieurement dans le Chapitre 4. Pour des raisons de simplicité, je me concentre ici sur la dynamique des particules en laissant de côté les effets thermiques, les changements de rayon ou de masse des particules. Par ailleurs, je me suis focalisé sur le cas de particules colloïdales, c'est-à-dire ayant une taille inférieure à quelques micromètres. Cette taille étant bien plus faible que l'échelle spatiale de Kolmogorov $r_p \ll \eta$, cela permet de considérer ces objets comme des particules ponctuelles et d'utiliser des formules approchées pour les forces hydrodynamiques s'exerçant sur ces particules (voir Section 1.2.2). Les particules sont supposées être beaucoup plus denses que le fluide (ce qui permet de négliger les forces de masse ajoutée) et la force de Basset est également négligée (pour des raisons de simplicité). Enfin, mes travaux ont été réalisés dans le cas de suspensions diluées, c'est-à-dire avec très peu de particules en suspension. Cela permet de négliger les effets de collision entre particules, qui seront abordées dans le Chapitre 3.

**Méthodologie** : La méthodologie suivie dans ce chapitre illustre la manière dont un modèle est développé en respecant les objectifs, les hypothèses et autres contraintes imposées. Dans le cas présent, un cadre de modélisation approprié consiste à adopter un point de vue Lagrangien stochastique pour le suivi de particules. Ces approches proposent de simuler la PDF (fonction de densité de probabilité) d'un certain nombre de variables qui ont été sélectionnées pour décrire des
informations pertinentes sur la phase dispersée et qui forment le vecteur d'état $\mathbf{Z}_p$. Ainsi, elles sont aussi appelées approches Lagrangienes PDF. Ces approches ont été initialement développées dans les années 1980 et 1990 (cf. travaux de Pope et ses collaborateurs [57, 120]) puis adaptées pour traiter le cas d'écoulements à phase dispersée (cf. travaux de Minier [97, 98, 100]).

L'intérêt de ces approches est multiple. D'abord, le choix d'une approche Lagrangienne permet de respecter les contraintes (ii.) et (iv.) qui nécessitent de traiter certains aspects au niveau local pour chacune des particules. Ensuite, le choix d'une formulation stochastique est lié aux critères (i.) et (iii.) qui requièrent des méthodes applicables dans des situations réelles. Ces critères rendent les approches microscopiques (typiquement via des simulations DNS pour le fluide) hors de propos. En conséquence, il est nécessaire d'avoir recours à des modèles de turbulence, puisqu'ils permettent de réaliser des simulations à grande échelle. Cela donne lieu à des approches hybrides où la phase fluide est résolue avec une approche par champ alors que la phase dispersée est simulée par un suivi Lagrangien. De plus, comme les particules étudiées ici sont de très petites tailles, elles sont sensibles aux petites échelles de la turbulence, celles-ci même qui ne sont pas explicitement résolues par ces modèles de turbulence. Il est donc nécessaire de développer des modèles rendant compte de l'effet de ces fluctuations à petite échelle sur la dynamique de particules. C'est dans ce cadre que les approches Lagrangienes stochastiques prennent tout leur sens puisqu'elles distinguent les variables rapides (modélisées) et les variables lentes (explicitement calculées). Je donne dans le paragraphe suivant, plus de détails sur ces approches.

Un système d'équations différentielles stochastiques (EDS) est utilisé pour décrire l'évolution dans le temps des variables attachées à la particule. Ces variables composent le vecteur d'état de la particule $\mathbf{Z}_p$, qui inclut ici la position, la vitesse et la vitesse du fluide vu par la particule : $\mathbf{Z}_p = (\mathbf{X}_p, \mathbf{U}_p, \mathbf{U}_s)$. L'équation d'évolution est alors donnée par :

$$dX_{p,i} = U_{p,i} dt, \qquad (2.1a)$$

$$dU_{p,i} = \frac{U_{s,i} - U_{p,i}}{\tau_p} dt + A_{p,i} dt, \qquad (2.1b)$$

$$dU_{s,i} = A_{s,i}(t, U_{s,i}) dt + B_{s,i}(t, U_{s,i}) dW_i. \qquad (2.1c)$$

La première équation (2.1a) donne simplement la dynamique de la position de la particule en intégrant sa vitesse. La seconde équation (2.1b) ne fait intervenir que deux contributions : le 1er terme correspond à la force de traînée alors que le 2nd terme regroupe les autres forces $\mathbf{A}_p$ (ici la gravité). Ces deux contributions sont les termes dominants dans l'équation BBO puisque $\rho_p \gg \rho_f$. La force de traînée est exprimée en fonction du temps de relaxation de la particule $\tau_p$ défini comme suit :

$$\tau_p = \frac{\rho_p}{\rho_f} \frac{8 r_p}{3 C_D \| \mathbf{U}_s - \mathbf{U}_p \|}, \qquad (2.2)$$

avec $\rho_f$ la densité du fluide, $\rho_p$ celle des particules, $r_p$ le rayon des particules et $C_D$ le coefficient de traînée. À noter que l'on néglige la force de masse ajoutée (puisque $\rho_p \gg \rho_f$) et la force de Basset (par simplicité). La force de traînée fait aussi apparaître la vitesse du fluide vu par la particule $\mathbf{U}_s$, qui est définie comme la vitesse instantanée du fluide échantillonnée à l'emplacement de la particule $\mathbf{U}_s(t) = \mathbf{U}_f(t, \mathbf{X}_p(t))$. Comme indiqué dans le paragraphe précédent, la vitesse instantanée à la position d'une particule ne peut pas être obtenue ici à partir de l'information réduite disponible sur le champ de vitesse du fluide (typiquement, ses premier et deuxième moments). C'est pourquoi, une équation supplémentaire (2.1c) est ajoutée pour la variable $\mathbf{U}_s$. À noter que cette équation supplémentaire n'est pas nécessaire lorsque la simulation fluide donne déjà cette information (comme dans le cas d'une DNS). Dans l'équation (2.1c), cette vitesse du
fluide vu $U_s$ est représentée par un modèle de Langevin : le terme de dérive $A_{s,i}$ regroupe les effets ayant des variations lentes (comme le terme de gradient de pression) alors que le terme de diffusion $B_{s,i}$ modélise les effets des fluctuations rapides (comme la turbulence à petite échelle). Chacun des termes $dW_{i,s}$ correspond aux incréments d'un processus de Wiener. Le détail des expressions retenues pour ces termes de dérive et de diffusion est disponible dans [100, 97, 101, 99].

Pour la suite de ce chapitre, il est important de réaliser que le recours à une approche Lagrangienne PDF implique que les objets numériques suivis dans chaque simulation correspondent à des particules stochastiques. Chacune de ces particules, modélisée par des EDS telles que l'équation (2.1), doit être considérée comme un échantillon de la PDF correspondante [100, 98, 99]. Dans le cadre actuel, je traite donc essentiellement une PDF lagrangienne à deux temps et à une particule. Pour le dire différemment, cela correspond à simuler la dynamique de particules indépendantes les unes des autres, puisque les nombres aléatoires générés pour chaque particule sont indépendants des autres particules. Malgré tout, cela permet d'extraire des informations statistiques, comme par exemple une PDF eulérienne à un temps et à un point (en intégrant les variables d'intérêt sur chaque élément de volume défini dans le maillage utilisé).

Dans la suite de ce chapitre, j'illustre comment ce genre de modèles a été enrichi pour rendre compte des phénomènes de remise en suspension des particules (cf. Section 2.2) et des phénomènes liés à la dynamique des particules en proche paroi (cf. Section 2.3).

## 2.2 Phénomènes de resuspension de particules

Depuis ma thèse, j'ai continuellement travaillé sur le sujet de la remise en suspension de particules déposées sur une surface et exposée à un écoulement turbulent.

**Méthodologie** : La modélisation de la remise en suspension de particules est rendue complexe par la gamme d'échelles temporelles et spatiales impliquée dans ce phénomène, par les diverses façons dont une particule peut être remise en suspension et aussi par la variété des forces en jeu (voir les revues sur le sujet [137, 109, 157, 49, 59, 113, 130]).

En matière de phénoménologie, la resuspension d'une particule isolée peut se produire de plusieurs manières (voir Figure 2.1) : par arrachement direct (c'est-à-dire qu'une particule est immédiatement détachée de la surface), par glissement (où une particule glisse sur la paroi avant d'être détachée), par roulement (lorsqu'une particule roule sur la paroi avant d'être détachée) ou par une combinaison de ces mécanismes. Une fois la particule détachée de la surface, elle est ré-entraînée par le fluide et peut soit se re-déposer sur la surface plus en aval soit s'éloigner de la surface et rester en suspension. L'ensemble de ces phénomènes correspond à la remise en suspension d'une particule.

<figure><img src="image_17.png" /><figcaption>FIGURE 2.1 – Illustration des 3 mécanismes de resuspension d'une particule isolée.</figcaption></figure>
En termes de forces, le détachement d'une particule isolée a lieu lors d'une rupture de l'équilibre entre les différentes forces ou des moments de ces forces. En écrivant le bilan des forces/moments, on obtient les conditions suivantes :

$$ \text{Lift-off} \quad \mathbf{F}_{\text{lift}} > \mathbf{F}_{\text{grav}} + \mathbf{F}_{\text{adh}}, \qquad (2.3a) $$

$$ \text{Sliding} \quad \mathbf{F}_{\text{drag}} > \mu_{\text{fr}} (\mathbf{F}_{\text{adh}} + \mathbf{F}_{\text{grav}} - \mathbf{F}_{\text{lift}}), \qquad (2.3b) $$

$$ \text{Rolling} \quad \mathbf{M}_{\text{drag}} + \mathbf{M}_{\text{lift}} > \mathbf{M}_{\text{adh}} + \mathbf{M}_{\text{grav}}, \qquad (2.3c) $$

où $\mathbf{F}_{\text{lift}}$ (resp. $\mathbf{M}_{\text{lift}}$) est la force hydrodynamique de portance (resp. le moment des forces), $\mathbf{F}_{\text{drag}}$ la force hydrodynamique de traînée, $\mathbf{F}_{\text{grav}}$ la force liée à la gravité (le poids) et $\mathbf{F}_{\text{adh}}$ la force d'adhésion entre la particule et la paroi. L'équation (2.3) fait aussi apparaître le coefficient de friction statique $\mu_{\text{fr}}$ (qui mesure la résistance au déplacement lorsque deux surfaces sont en contact). La résolution de telles équations nécessite d'utiliser des formules pour exprimer chacune de ces forces (empiriques ou théoriques).

En matière de modèle, la description retenue ici est fondée sur un scénario de resuspension en trois étapes (voir aussi la Figure 2.2) [53, 60] : (i) une particule est d'abord mise en mouvement par l'action du fluide lorsqu'il y a rupture d'équilibre des moments, (ii) elle se déplace ensuite sur le substrat rugueux par roulement, et (iii) elle se détache finalement de la surface pour réintégrer l'écoulement en entrant en collision avec une aspérité présente sur la surface.

<figure><img src="image_18.png" /><figcaption>FIGURE 2.2 – Illustration du scénario en trois étapes utilisé dans le modèle de resuspension : (a) la mise en mouvement (lorsque l'équilibre est rompu) ; (b) le roulement sur la surface et (c) le détachement par impact sur des rugosités de grande taille.</figcaption></figure>

Le choix de ne retenir que le mécanisme de roulement dans la resuspension de particules est motivé par le fait que les particules étudiées ici sont de taille colloïdale, donc entièrement immergées dans la sous-couche visqueuse. Dans ce cas, les expériences semblent indiquer que le roulement est le mécanisme prépondérant [59]. Cela permet aussi de négliger la contribution des moments des forces de portance par rapport à celles de la traînée. A partir de ces différentes hypothèses, on peut obtenir directement la vitesse de translation de la particule (supposée proportionnelle à la vitesse angulaire $U_{p,//} \approx r_p \Omega_p$) :

$$ I_p \frac{dU_{p,//}}{dt} \approx r_p M_{\text{drag},//} - r_p M_{\text{adh}} \qquad (2.4) $$

où $I_p = 7m_p r_p^2/5$ est l'inertie de rotation de la particule (avec $m_p$ sa masse). Cette équation nécessite de calculer le moment des forces de traînée et celui des forces d'adhésion. Le moment des forces de traînée est obtenu en exprimant les forces dans l'approximation de Stokes et en appliquant la formule de O'Neill [111], qui donne :

$$ M_{\text{drag}} = 1.4 r_p F_{\text{drag}} = 1.4 r_p \times 6 f \pi \rho_f \nu_f r_p \mathbf{U}_r \qquad (2.5) $$
où *f* est un facteur correctif pour l'anisotropie de l'écoulement en proche paroi [111].

De manière similaire, le moment de forces d'adhésion est explicitement calculé via un mo-
dèle stochastique. La stochasticité du modèle d'adhésion vient du fait que les surfaces sont ru-
gueuses (au moins à l'échelle nanométrique et parfois à des échelles plus micrométriques voire plus
grandes). Ces imperfections sur la surface induisent une grande distribution des forces d'adhé-
sion, puisque les zones réellement en contact dépendent de l'état de surface autour du point
de contact. Dans mes travaux, les surfaces rugueuses sont décrites par un plan lisse recouvert
d'aspérités hémisphériques : leur taille varie selon les rayons mesurés expérimentalement et leur
nombre étant donné par la couverture de surface de ces aspérités. La distribution de la force
d'adhésion avec de telles surfaces rugueuses est obtenue grâce à des simulations Monte-Carlo où
les surfaces rugueuses sont générées aléatoirement. Ce modèle sur le moment des forces d'adhé-
sion fait intervenir plusieurs informations : (a) le nombre d'aspérités en contact avec la particule,
et (b) la distance au point de pivot (le point de contact le plus en aval de la particule).

**Implémentation numérique** : Ce modèle a été implémenté dans le code CFD Code_Saturne.
Le détail de l’algorithme est donné dans l’article D.1. Dans cet article, un effort particulier a
été réalisé pour permettre de simuler la dynamique de particules en roulement sur la paroi.
En effet, cela nécessite de bien calculer les deux paramètres intervenant dans le modèle des
moments d’adhésion (à savoir le nombre d’aspérités en contact et la distance au point de pivot).
Or, lorsqu’une particule roule sur une surface donnée, les statistiques obtenues pour ces deux
paramètres doivent être indépendantes du déplacement de la particule pendant un laps de temps
(donc du pas de temps utilisé et de la vitesse de roulement de la particule). Cela a nécessité de
développer des algorithmes prenant en compte les zones de recouvrement des surface de contact
entre deux pas de temps (lorsqu’elles existent, voir Figure 2.3). Le détail de ces algorithmes est
donné dans l’article D.1.

<figure><img src="image_19.png" /><figcaption>(a) Cas sans recouvrement (grands déplacements) (b) Cas avec recouvrement (petits déplacements)</figcaption></figure>

FIGURE 2.3 – Illustration de la façon dont le calcul du nombre d'aspérités en contact avec une particule dépend du déplacement d'une particule pendant un pas de temps.

**Résultats sur la resuspension** : Le modèle des moments et forces d'adhésion a été validé durant mes travaux de thèse [61]. Cette validation a montré l'intérêt de ce genre de modèles, qui permet de reproduire diverses formes de la distribution des forces d'adhésion (Gaussiennes, log-normales ou bien à plusieurs modes) sans avoir à se donner a-priori une forme de la distribution. Suite au développement d'algorithmes dédiés dans l'article D.1, les résultats numériques obtenus avec le modèle complet ont été comparés à des données expérimentales sur la resuspension de particules. Par exemple, dans l'article D.2, la quantité de particules restantes sur la paroi après
avoir été exposées à un écoulement (pendant une durée de 1 s) a été caractérisée en fonction de la vitesse de frottement du fluide. Comme on peut le voir sur la Figure 2.4a, les résultats obtenus pour des particules de verre de 20, 30 ou 40 µm sont en accord avec les données expérimentales lorsque la rugosité de surface comprend des petites aspérités de 9 nm couvrant 3.14 % de la surface et des grandes aspérités de 0.5 µ couvrant 4 % de la surface. En particulier, cela valide la capacité du modèle à prédire le seuil médian de la vitesse fluide $u_{50}$%, qui correspond à la vitesse fluide à laquelle 50 % des particules ont été remises en suspension.

<figure><img src="image_20.png" /><figcaption>FIGURE 2.4 – Fraction des particules de GBM restant sur un substrat de verre lisse en fonction de la vitesse de friction : données expérimentales (symboles et barres d'erreur), résultats numériques (lignes) obtenus pour deux surfaces rugueuses ayant la même rugosité moyenne $R_a$.</figcaption></figure>

Néanmoins, ces travaux ont aussi révélé les limites du modèle. En effet, les résultats obtenus sont beaucoup moins bons lorsque l'on impose des petites aspérités de 25 nm couvrant 2 % de la surface (voir Figure 2.4b). Or, ces deux mises en données de la rugosité de surface correspondent toutes les deux aux valeurs de rugosité moyenne $R_a$ mesurée expérimentalement. Cette rugosité moyenne correspond à la moyenne arithmétique de la hauteur de rugosité (absolue) au-dessus d'un plan de référence. Ces résultats confirment que des informations moyennes sur la rugosité de surface ne sont pas suffisantes pour caractériser entièrement une surface rugueuse. En plus, comme le montre la Figure 2.5, cette information moyenne sur la rugosité de surface ne permet pas d'étudier les phénomènes extrêmes de resuspension puisqu'ils sont liés aux queues de distribution de la force d'adhésion (elle-même en lien avec les valeurs extrêmes de la rugosité de surface). Pour mieux mettre en lumière ces phénomènes, des simulations ont été réalisées en modifiant les 10% des valeurs maximales obtenus par le modèle dans le 1er cas (aspérités de 9 nm couvrant 3.14 % de la surface) et en les multipliant par un facteur aléatoire entre 1 et 3. Cela permet d'obtenir une queue de distribution des forces/moments d'adhésion qui s'étale sur une plage plus importante et, ainsi, de mieux reproduire les observations sur la resuspension aux grandes valeurs de vitesse fluide. Cela peut se comprendre simplement en réalisant que, par cette manipulation arbitraire, les 10 % des particules les plus adhérentes à la paroi se sont retrouvées encore plus collées à la paroi et donc plus difficiles à resuspendre. D'où une queue de distribution sur la resuspension plus étalée.
<figure><img src="image_21.png" /></figure>

FIGURE 2.5 – Illustration des valeurs moyennes de rugosité extraites d'un profil de surface (à gauche en haut) comparées à une mesure complète des distributions (à gauche en bas). Analyse de l'effet d'amplification des forces d'adhésion les plus importantes sur la fraction de particules restantes sur la surface (zoom sur les queues de distributions).

Ainsi, si l'on cherche maintenant à avoir des approches qui permettent de prédire non seulement le seuil médian de la vitesse fluide $u_{50\%}$ mais plus généralement n'importe quelle seuil de la vitesse $u_{k\%}$ (avec $k \in [0, 1]$), il devient absolument nécessaire d'extraire des profils de surface non pas des valeurs moyennes (comme $R_a$) mais bien toute la distribution des tailles de rugosité et de leur couverture de surface.

**Effets de collisions entre particules** : Plus récemment article D.3, des travaux réalisés en partenariat avec Grégory Lécrivain au HZDR à Dresden ont permis de révéler le rôle de collisions entre des particules micrométriques lorsque la surface est recouverte d'un grand nombre de particules [10]. En effet, lorsque la surface couverte par les particules initialement déposées sur la surface atteint environ 10 %, des collisions successives entre particules ont été observées expérimentalement (voir Figure 2.6).

<figure><img src="image_22.png" /><figcaption>(a) Schéma des collisions entre particules induisant une resuspension.</figcaption></figure>

<figure><img src="image_23.png" /><figcaption>(b) Clichés consécutifs de billes de verre de 40 µm montrant une suite d'événements de resuspension dûs aux collisions.</figcaption></figure>

FIGURE 2.6 – Illustration du rôle des collisions dans la resuspension.
Ces observations viennent corroborer des indices récents indiquant le rôle des collisions dans la resuspension [134] : dès lors qu'une particule micrométrique se met en roulement sur la surface, elle peut entrer en collision avec d'autres particules déposées si la concentration surfacique de particules est suffisamment élevé. Cela induit même des phénomènes en cascade (comme avec des boules de billard). Cela signifie que la resuspension par collision de particules n'est pas une suite d'événements indépendants (comme dans le cas des dépôts monocouche éparses où les particules sont isolées les unes des autres) mais est lié à des événements corrélés. Cela ouvre la voie à de nouveaux modèles qui rendent compte de ces phénomènes.

Une première version grossière d'un modèle incluant l'effet de collisions entre particules a été réalisée dans l'article D.2. Pour tenir compte de la propagation des collisions, le modèle a été étendu comme suit : lorsqu'une particule commence à se déplacer (c'est-à-dire lors d'une rupture d'équilibre donné par les équations (2.3)), la probabilité de collision avec une autre particule déposée est évaluée avec des critères géométrique simples (utilisant le temps de séjour en proche paroi de la particule en mouvement et la surface couverte par les particules encore déposées). Lorsqu'un candidat est trouvé, une collision a lieu. Ce premier modèle simplifie aussi le résultat de la collision : la particule en mouvement est immédiatement resuspendue et celle qui vient d'être impactée se met en mouvement. Malgré sa simplicité, ce modèle a permis de confirmer le rôle des collisions successives dans la resuspension de particules mais aussi le rôle de la répartition initiale des particules sur la surface (voir Figure 2.7) : en effet, les résultats observés expérimentalement ne sont reproduisibles numériquement qu'en forçant les particules à former des clusters denses sur la surface initialement (et non en partant d'un dépôt parfaitement homogène).

<figure><img src="image_24.png" /><figcaption>FIGURE 2.7 – Évolution de la fraction restante en fonction de la vitesse de l'écoulement : comparaison entre les mesures et les simulations numériques avec/sans modèle de collisions successives montrant le rôle du clustering initial des particules.</figcaption></figure>

## 2.3 Dynamique en proche paroi : conditions de bord

Une autre de mes contributions a porté sur le développement de modèles pour la dynamique de particules en proche paroi. En effet, une fois les particules remises en suspension, ces dernières vont alors se mouvoir dans la région proche de la surface. Ces mouvements complexes sont responsables de phénomènes comme la saltation ou la reptation des particules sur les surfaces. Ainsi, il est nécessaire d'avoir des modèles rendant correctement compte de cette dynamique.
<figure><img src="image_25.png" /><figcaption>FIGURE 2.8 – Schéma de l’interaction paroi/particule montrant comment les particules arrivant à $z = z_{mirror}$, choisi dans la couche logarithmique, sont réfléchies dans le domaine de manière à respecter un flux de masse nul dans la direction normale à la paroi ainsi que l’échange de quantité de mouvement dû à la friction dans la direction longitudinale. Les valeurs $U_p^{in}$ (orientées vers le haut) et $U_p^{out}$ (orientées vers le bas) représentent respectivement les vitesses d’entrée et de sortie.</figcaption></figure>

**Méthodologie** : L’un des principaux enjeux dans les approches hybrides est d’utiliser des formulations pour les conditions aux limites qui soient consistantes entre le modèle pour la phase fluide et le modèle pour la phase dispersée.

Dans le cas présent, il s’agit d’appliquer une condition de bord adéquate lorsqu’une particule rebondit sur une paroi. En pratique, dans une simulation numérique avec une approche hybride, le fait que le pas de temps soit souvent choisi par rapport à la dynamique du fluide implique que les particules en proche paroi peuvent se déplacer sur des distances plus grandes que plusieurs fois leur diamètre. Néanmoins, il est tout à fait possible de détecter la collision avec une paroi, simplement en repérant les instants où une particule se retrouve en dehors du domaine physique. Lorsqu’une collision est détectée, il faut alors appliquer une condition de bord pour rendre compte du rebond de la particule sur la paroi. En supposant des collisions purement élastiques, il est naturel de penser en premier lieu à appliquer une condition de rebond spéculaire (cf. Fig. 2.8). Cela consiste simplement à réfléchir la position de la particule par rapport au plan formé par la surface et à ne modifier que la composante normale de la vitesse de la particule $U_p$. En pratique, en supposant que la normale de la surface est orientée selon la direction $z$ comme dans la figure 2.8, cela donne :

$$U_{p,x}^{in} = U_{p,x}^{out}, \qquad (2.6a)$$

$$U_{p,y}^{in} = U_{p,y}^{out}, \qquad (2.6b)$$

$$U_{p,z}^{in} = -U_{p,z}^{out}. \qquad (2.6c)$$

où les exposants $^{in}$ (resp. $^{out}$) correspondent aux vitesses des particules entrantes (resp. sortantes) du domaine.

Cependant, cette condition de rebond spéculaire n’est pas compatible avec les conditions de bord qui sont appliquées dans les modèles de turbulence (comme RANS). En effet, cette condition revient à écrire une condition de flux nul à la paroi pour la vitesse parallèle à la paroi, négligeant ainsi les échanges de quantité de mouvement dans les régions proches de la paroi. En suivant les idées utilisées pour formuler les conditions de bord dans les modèles de turbulence (où la condition est appliquée dans la couche logarithmique et non pas à la paroi directement) [119],
il est possible de dériver une formulation qui rend compte du cisaillement en proche paroi. Cela donne (cf. [35, 102] pour plus de détails) :

$$ U_{p,x}^{in} = U_{p,x}^{out} - 2 \frac{\langle U'_{f,x} U'_{f,z} \rangle}{\langle U'^2_{f,z} \rangle} U_{p,z}^{out}, \qquad (2.7a) $$

$$ U_{p,y}^{in} = U_{p,x}^{out} - 2 \frac{\langle U'_{f,y} U'_{f,z} \rangle}{\langle U'^2_{f,z} \rangle} U_{p,z}^{out}, \qquad (2.7b) $$

$$ U_{p,z}^{in} = -U_{p,z}^{out}, \qquad (2.7c) $$

À noter que les équations données précédemment peuvent aisément se généraliser au cas de collisions inélastiques. Pour ce faire, il suffit d'appliquer un coefficient de restitution $e$, qui donne le ratio entre la quantité d'énergie avant le rebond et après le rebond.

**Résultats sur les tenseurs de Reynolds** : Pour illustrer l'importance d'avoir recours à des conditions de bord consistantes dans les approches hybrides, des simulations numériques ont été réalisées dans le cas d'un canal plan turbulent dans lequel des traceurs (particules fluides) ont été injectées. Les détails sont décrits dans l'article C.1. L'intérêt de ce cas est de vérifier que les résultats obtenus pour la phase fluide et pour la phase dispersée sont les mêmes, surtout pour les tenseurs de Reynolds puisqu'ils mesurent la quantité d'énergie dans en proche paroi (et sont donc directement impactés par les conditions de bord). Les profils obtenus sont visibles sur la Figure 2.9, où l'on voit clairement que les résultats sont meilleurs avec une condition de bord consistance avec celle du fluide plutôt qu'avec la condition de rebond spéculaire. En particulier, l'amélioration des résultats est significative pour les composantes du tenseur de Reynolds $R_{11}$ et $R_{13}$. On peut aussi remarquer que les résultats pour la composante $R_{33}$ ne sont pas affectés par les conditions aux limites puisque la composante de la vitesse normale à la paroi reste inchangée entre les deux formulations.

<figure><img src="image_26.png" /><figcaption>FIGURE 2.9 – Profils des tenseurs de Reynolds dans un canal plan : comparaison entre les résultats obtenus avec la simulation CFD pour la phase fluide (trait plein) et avec un modèle Lagrangien stochastique pour la phase dispersée (symboles).</figcaption></figure>
## 2.4 Perspectives

Ces travaux ouvrent un certain nombre de pistes pour des travaux futurs dont je liste ci-dessous les plus importants selon moi :

**Dynamique corrélée** : Le développement de nouveaux modèles pour rendre compte des collisions entre les particules nécessite d'inclure des informations sur les corrélations temporelles et spatiales entre les événements de resuspension. Cela peut être réalisé grâce à une analyse approfondie des données expérimentales, notamment pour extraire des informations sur la fréquence de ces collisions et sur les résultats de ces collisions (bilan de quantité de mouvement, bilan d'énergie). L'un des enjeux serait de développer un modèle de collision similaire à ceux utilisés pour les collisions de particules dans l'écoulement (avec un coefficient de restitution pour les collisions inélastiques). En parallèle, des informations précises sur les corrélations entre les collisions et les corrélations liés aux structures turbulentes en proche paroi peuvent être extraites de simulations numériques directes de type PR-DNS (collaborations futures envisagées avec Bernhard Vowinckel sur ce sujet).

Ces travaux ouvriront aussi la voie au développement de modèles Lagrangiens stochastiques incluant des informations à deux points en espace.

<figure><img src="image_27.png" /><figcaption>FIGURE 2.10 – Illustration des différents mécanismes en jeu dans la resuspension des particules en fonction de la structure du dépôt et de l'intensité du flux (pour une taille, une forme et une densité de particules fixés).</figcaption></figure>

**Modèle unifié de resuspension** : Plus généralement, un modèle rendant compte des collisions entre particules dans le cas d'un dépôt monocouche dense ouvrira la voie à des modèles unifiés pour la remise en suspension de particules depuis des dépôts monocouches éparses aux dépôts multicouches. En effet, la resuspension depuis des dépôts multicouches est bien plus complexe, notamment à cause du rôle de la morphologie du dépôt (dans les
forces de cohésion entre particules), du rôle joué par les corrélations spatiales ou temporelles sur la fragilisation de tels dépôts mais aussi de l'influence des effets collectifs (les particules pouvant se déplacer par paquet suite aux collisions très fréquentes). Il sera alors possible de voir émerger une description allant continûment des dépôts monocouches éparses aux dépôts multicouches (comme celui dessiné dans la Figure 2.10).

**Simulations de dépôts multicouches** Enfin, de telles études approfondies sur la resuspension de particules permettront de mieux rendre compte des effets de la morphologie du dépôt et de son evolution au cours du temps sous l'effet des interactions avec le fluide ou des impacts par d'autres particules en suspension. Cela ouvrira la voie à des modèles plus fins pour la formation de tels dépôts, notamment ceux qui permettent de déduire la morphologie du dépôt formé comme un résultat des interactions hydrodynamiques sur le dépôt et des forces de cohésion entre les particules. En effet, les simulations actuelles ont tendance à considérer deux cas extrêmes : soit la particule se dépose à l'endroit où elle entre en contact avec le dépôt déjà formé, soit la particule continue à se déplacer (en roulant/glissant) jusqu'à atteindre une situation d'équilibre (avec 3 points de contact avec des particules formant le dépôt). Le premier cas donne lieu à des dépôts très poreux alors que le second donne lieu à des dépôts très denses. En réalité, la présence de rugosité sur la particule (qui induit des forces de friction localement fortes) peut donner à des situations d'équilibre décrivant toute la gamme de dépôts entre ces deux cas extrêmes.
# Chapitre 3

## Agglomération et fragmentation de particules

### 3.1 Introduction

**Contexte** Jusqu’à présent, j’ai décrit les modèles pour la dynamique de particules sphériques en considérant le cas de suspension suffisamment diluées pour que les collisions entre particules soient négligeables. En réalité, les particules solides en suspension dans un écoulement peuvent se rencontrer et, lorsque les interactions à courte portée le permettent, se coller les unes aux autres pour former des agrégats [148, 95, 21]. La fréquence de collisions entre des particules augmente avec leur concentration dans l’écoulement (sans forcément avoir besoin d’aller jusqu’à des suspensions très denses) et elle dépend également de l’inertie des particules en suspension dans un écoulement turbulent.

Ce processus est communément appelé agglomération ou agrégation de particules solides (le terme coalescence étant réservé pour les gouttes/bulles et le terme flocculation étant employé pour les polymères). La fragmentation correspond au phénomène opposé, où un agrégat déjà formé se brise en plusieurs éléments de taille plus petite sous l’effet de contraintes (liées aux forces hydrodynamiques ou à des collisions).

Les processus d’agglomération et de fragmentation sont très fréquents dans la nature. Par exemple, la formation d’agrégats joue un rôle dans la géo-morphologie (avec la sédimentation de flocs formés rapidement par les sédiments en suspension dans les rivières lorsqu’ils atteignent les deltas [96]), dans la météorologie (e.g., la formation des nuages [48] ou la croissance de gouttelettes [42, 139]) ou bien dans l’astrophysique (avec la formation des planètes par croissance des agrégats de poussière [14]). L’agglomération de particules se produit également dans un large éventail de situations industrielles, comme dans les usines de traitement des eaux usées (séparation des solides par décantation des agrégats [121] ou procédés de flottation [136]), dans les systèmes de combustion (croissance des particules de suie [45]), ou dans l’industrie alimentaire (poudres de lait sec [58]). Enfin, l’agglomération a des implications dans le domaine médical (avec l’agrégation des protéines étant liées à certaines maladies dégénératives [13]).
**Verrous** : L'une des principales difficultés des phénomènes d'agglomération et de fragmentation est qu'ils nécessitent de rendre compte de plusieurs effets à la fois, parmi lesquels :

○ Les collision entre des particules, dont la fréquence dépend de l'écoulement (cisaillement pur, turbulence homogène isotrope), de la dynamique des particules (transport turbulent, sédimentation, mouvement Brownien) [38]. D'un point de vue numérique, la détection d'événements de collision nécessite aussi de développer des algorithmes spécifiques [142].

○ Le résultat de chaque événement de collision, chacune pouvant donner lieu soit à la formation d'un agrégat (lorsque les deux particules se collent entre elles), soit à une simple collision (rebond élastique ou inélastique) voire même à la fragmentation (lorsque des agrégats entrent en collision). Le résultat de chaque collision dépend de plusieurs paramètres, dont la vitesse relative, l'angle de collision, les propriétés des particules (composition chimique, charge, forme et taille) [21] ainsi que des forces à courte ou longue portée agissant entre les particules (comme les forces de capillarité ou les forces de contact) [39].

○ La morphologie des agrégats, qui dépend évidemment de la manière dont les agrégats ont été formés mais aussi de leur âge et de leur histoire. En effet, les agrégats peuvent se restructurer sous l'effet de contraintes extérieures (collisions, cisaillement du fluide) tandis que les liaisons entre les particules individuelles peuvent se consolider avec le temps (à cause des effets liés à la croissance de surface ou au sintering) [37].

○ La dynamique des agrégats, notamment leur sédimentation dans un écoulement (qui dépend de la masse/taille totale de chaque agrégat mais aussi de sa morphologie) [38].

Ces phénomènes sont étroitement couplés les uns aux autres, puisque la fréquence de collision entre des agrégats est pilotée par leur dynamique, ... qui dépend de leur morphologie, ... elle-même résultante de la façon dont chaque agrégat s'est formé, ... etc. C'est ce qui rend si ardu la modélisation de tels phénomènes.

**Objectifs** : L'objectif de ce chapitre est de présenter différents types de modèles pour l'agglomération et la fragmentation de particules, chacun d'eux permettant de reproduire un certain nombre des effets identifiés précédemment. Ainsi, en lien avec l'art de la modélisation (cf. Chapitre 1), ce chapitre illustre le fait que chaque modèle prend tout son sens en fonction des choix qui ont été faits (ici par rapport au niveau de description et aux paramètres de sortie).

**Hypothèses** : Dans ce chapitre, je me concentre sur le cas de particules sphériques, ayant une taille colloïdale et dont la concentration reste faible (typiquement dans le régime dilué). Cela permet de simplifier la dynamique des particules individuelles et de les considérer comme des particules ponctuelles. Par ailleurs, le choix de s'intéresser à des particules colloïdales implique de rendre compte des phénomènes de physico-chimie des interfaces qui sont responsables de la cohésion entre chacun des grains individuels.

**Méthodologie** : Au cours de mes recherches, j'ai développé divers modèles d'agglomération répondant à un cahier des charges précis. Chaque modèle a donc ses propres spécificités, notamment selon le niveau de description qui a été retenu et surtout selon le type d'information qui est supposé connu (les données d'entrée du modèle) et le type d'information qui est prédit par le modèle (les données de sortie). Dans ce chapitre, j'ai donc choisi de faire une distinction entre les modèles qui permettent de prédire le noyau d'agglomération (fréquence + résultat des collisions) et les modèles qui requièrent une information sur ces noyaux d'agglomération. Pour chacun de ces modèles que j'ai développés, je les classe ensuite en fonction de leur niveau de description.
## 3.2 Simulations directes de l'agglomération

La première catégorie de modèles auxquels j'ai contribué permet de calculer explicitement la fréquence des collisions entre les particules ainsi que le résultat des collisions (formation d'un agrégat ou bien simple rebond). Certains d'entre eux (ceux avec le niveau de description le plus complet) permettent également de calculer la morphologie des agrégats ainsi obtenus.

### 3.2.1 Méthodes Monte-Carlo

Entre 2020 et 2021, j'ai eu l'opportunité de collaborer avec Jérôme Yon et José Moran (du département "Optiques et Laser" du CORIA à Rouen) sur le développement d'algorithmes utilisant une approche de type Monte-Carlo pour la simulation de la morphologie d'agrégats.

**Méthodologie** : Les approches DEM, où les N-particules sont suivies simultanément, permettent de calculer la morphologie des agrégats puisque la position de chacun des grains individuels est connue au cours de la simulation. Lorsque ces approches DEM sont couplées à des simulations DNS pour un écoulement turbulent, cela permet également de calculer la fréquence de collision entre les particules, puisque les collisions entre particules sont pilotées par leur dynamique dans l'écoulement. Dans le cas présent, j'ai participé à des travaux utilisant une approche de type Monte-Carlo, c'est-à-dire où la dynamique des particules est simplifiée pour faciliter les simulations. Dans l'approche MCAC développée au CORIA (cf. Fig. 3.1a), l'algorithme consiste à commencer chaque simulation avec une répartition aléatoire et homogène des particules dans le domaine, puis à déplacer tour à tour chaque particule (tirée au hasard parmi toutes les particules) dans une direction aléatoire et sur une distance donnée $\lambda_p$ [106]. Lorsqu'une collision est détectée (par des critères géométriques simples), les deux particules se collent exactement au premier point de contact. Les particules au sein d'un même agrégat sont ensuite considérées comme se déplaçant de manière concomitante.

<figure><img src="image_28.png" /><figcaption>(a) Illustration de l'algorithme MCAC pour l'agglomération.</figcaption></figure>

<caption>(b) Illustration de la manière dont le potentiel d'interaction est pris en compte pour calculer le résultat de chaque collision détectée dans l'algorithme MCAC amélioré.</caption>

FIGURE 3.1 – Schéma de l'algorithme MCAC amélioré pour le calcul de l'agglomération de particules en prenant en compte l'énergie d'interaction.

Dans l'article E.7, nous avons couplé cette approche MCAC à un calcul des forces d'interactions entre les particules. Pour ce faire, l'énergie d'interaction a été calculée en prenant en
compte un potentiel de Lennard-Jones (qui comprend donc une partie attractive liée aux forces de van der Waals et une partie répulsive liée au principe d'exclusion de Pauli) ainsi que l'énergie liée aux interactions électrostatiques (les particules pouvant être chargées). Le calcul de l'énergie d'interaction entre des agrégats composés de plusieurs particules a été simplifié en ne prenant en compte que les particules qui sont en contact entre les deux agrégats. En effet, comme les forces d'origine moléculaire sont de très courte portée, la contribution liée aux autres particules qui ne sont pas en contact est négligeable [9]. Ce calcul permet d'extraire ensuite la valeur maximum de l'énergie d'interaction, appelée barrière d'énergie $E_{barr}$, ainsi que la valeur du puits de potentiel d'adhésion $E_{well}$ (correspondant au minimum de l'énergie d'interaction).

Comme illustré sur la Fig. 3.1b, l'introduction de cette énergie d'interaction dans l'algorithme a été réalisée en considérant le scénario suivant :

1. Aucune collision n'a lieu lorsque l'énergie cinétique relative lors de l'impact $E_{kin}$ est plus faible que la barrière d'énergie (c'est-à-dire quand $E_{kin} < E_{barr}$);

2. Une collision a lieu lorsque l'énergie cinétique permet de franchir la barrière d'énergie (lorsque $E_{kin} > E_{barr}$). Dans ce cas, deux résultats sont envisagés :

2.a. Les deux éléments se collent lorsqu'ils ne peuvent pas s'échapper du puits de potentiel (dès que $E_{kin} < |E|_{well} + E_{barr}$);

2.b. Les deux éléments rebondissent lorsque l'énergie cinétique est suffisante pour sortir du puits de potentiel (si $E_{kin} > |E|_{well} + E_{barr}$);

En suivant ce scénario, une probabilité de collision ainsi qu'une probabilité d'adhésion ont été évaluées puis inclus dans l'algorithme MCAC (en faisant un tirage aléatoire pour déterminer le résultat de chaque collision).

**Résultats** : À partir de cet algorithme MCAC modifié, nous avons pu suivre l'évolution de la taille des objets au cours du temps. Les résultats sur le rayon géométrique moyen $R_{v,geo}$ (cf. Figure 3.2 ci-contre) montrent bien que la taille des agrégats augmente de façon monotone au cours du temps. En particulier, l'accroissement de la taille des agrégats est plus important dans le cas de suies jeunes que dans le cas de suies matures (à cause de la densité apparente plus faible des suies jeunes). Cependant, ces résultats montrent aussi que la prise en compte de la répulsion électrostatique n'a pas d'effet sur la taille des agrégats obtenus. Ces résultats confirment que, pour les suies, il est possible de négliger l'effet des charges sur le processus d'agglomération.

<figure><img src="image_29.png" /><figcaption>FIGURE 3.2 – Évolution de la moyenne géométrique du rayon des agrégats/particules au cours du temps.</figcaption></figure>
### 3.2.2 Simulations DNS

Durant mon post-doctorat à l'OCA, j'ai eu l'occasion de participer à l'encadrement de la thèse de Robin Vallée dont une partie des travaux a porté sur l'accrétion de particules.

**Méthodologie** : L'accrétion de petites particules sur une grosse sphère en chute libre dans un écoulement a été étudiée en réalisant des simulations DNS de l'écoulement autour d'une sphère fixée. La vitesse verticale $U_z$ est imposée, et ces valeurs ont été choisies pour couvrir toute une plage de nombre de Reynolds particulaire $Re_p = Ud_{sphere}/\nu_f$ (allant de 1 à 400). Ensuite, des particules ponctuelles ont été injectées de manière uniforme en amont de la grosse sphère. Leur dynamique est résolue en considérant uniquement les forces de traînée (les forces de masse ajoutée et de Basset sont négligées puisque l'on considère des particules beaucoup plus denses que le fluide). Par ailleurs, plusieurs tailles de particules ont été considérées, de manière à avoir des nombres de Stokes $St = \tau_p U/d_{sphere}$ allant de 0.1 à 1 ($\tau_p$ correspondant au temps de relaxation des particules défini dans l'équation (2.2)). Enfin, lorsque les particules ponctuelles impactent la sphère, le traitement de la collision a été simplifié en imposant un coefficient de restitution $e$ qui mesure le ratio entre l'énergie cinétique avant et après impact : le cas $e = 1$ correspond au cas de collision purement élastique alors que le cas $e < 1$ correspond à des collisions inélastiques (où une certaine quantité d'énergie est perdue lors de la collision, pouvant aller jusqu'à $e = 0$ pour des collisions dites collantes).

**Résultats** : Dans l'article E.4, nous avons caractérisé comment l'efficacité de collision dépend des deux paramètres adimensionnels du système (à savoir le nombre de Stokes $St$ et le nombre de Reynolds $Re_p$). Dans le cas de l'accrétion sur une sphère, cette fréquence de collision est définie comme le ratio entre le nombre de particules ponctuelles qui entre en collision avec la sphère et le nombre de particules qui sont contenues dans le volume de fluide balayé par la sphère. Cela représente donc bien la quantité de particules qui vont effectivement entrer en collision avec la sphère parmi toutes celles qui le pourraient.

Les résultats obtenus pour la fréquence de collision sont présentés dans la Fig. 3.3a en fonction du nombre de Stokes des particules. Comme attendu, cette fréquence de collision tend vers 1 dans le cas de grand nombre de Stokes. On constate aussi que cette fréquence tend vers 0 en deçà d'un certain seuil critique $St < St_c$ ($St_c$ diminue en fonction du nombre de Reynolds). Enfin, pour une valeur fixe du Stokes, on remarque que la fréquence augmente en fonction du nombre de Reynolds (ce qui est une conséquence de la diminution de $St_c$ en fonction de $Re_p$).

Par ailleurs, nous avons aussi exploré si des collisions successives peuvent engendrer une adhésion des particules ponctuelles à la paroi en un temps limité. Pour ce faire, la dynamique des particules ponctuelles ayant rebondies sur la sphère a été caractérisée en fonction de leur Stokes $St$, du Reynolds $Re_p$ et du coefficient de restitution $e$. Comme on peut s'y attendre, une diminution du coefficient de restitution engendre un nombre de collisions successives plus important puisque l'énergie perdue lors de chacune des collisions est plus importante (cf. Fig. 3.3b). Néanmoins, cela ne donne pas lieu à un nombre infini de collisions en un temps fini (comme sur le cas d'un plan [12]) puisque les particules se déplacent dans la direction transverse à la paroi jusqu'à partir dans le sillage de la sphère. Malgré tout, les particules ponctuelles pourraient adhérer à la sphère au bout d'un certain nombre de collisions successives, puisque l'énergie cinétique décroît à mesure que le nombre de collision augmente. Ainsi, en prenant en compte le rôle des forces d'adhésion (avec un critère énergétique comme dans le modèle MCAC), on pourrait facilement introduire un seuil sur les énergies cinétiques en deçà duquel les particules adhèrent à la sphère.
<figure><img src="image_30.png" /></figure>

(a) Évolution de la fréquence de collision en fonction du nombre de Stokes pour plusieurs valeurs du nombre de Reynolds (avec les cas limite de l'écoulement de Stokes $Re_p = 0$ et d'Euler $Re_p = \infty$).

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

(b) Exemples de trajectoires de particules ($St = 10$) dans un écoulement de Stokes pour trois valeurs du coefficient de restitution : $e = 0.5, 0.75$ et $1$.

FIGURE 3.3 – Résultats obtenus sur la fréquence de collision et sur les collisions successives entre des particules ponctuelles et une sphère dans un écoulement.

### 3.2.3 Méthode à probabilité de collisions

Durant mon post-doctorat à l'IMP PAN, j'ai également travaillé sur une approche probabiliste pour évaluer la fréquence et l'efficacité de collision pour des particules dans le régime diffusif.

**Méthodologie** : Ces travaux de recherche ont été effectués en se focalisant sur le cas de particules ayant un mouvement purement diffusif dans un fluide au repos. Dans la limite sur-amortie ($\Delta t \gg \tau_p$), il est alors possible d'exprimer la dynamique de ces particules diffusives uniquement avec une équation sur la position des particules :

$$d\mathbf{X}_p(t) = B_x d\mathbf{W}(t) \qquad (3.1)$$

où $B_x$ est le coefficient de diffusion ($B_x = \sqrt{2k_B T \tau_p / m_p}$ dans le cas du mouvement Brownien à la température $T$, avec $k_B$ la constante de Boltzmann) et $d\mathbf{W}$ l'incrément d'un processus de Wiener. À noter que, dans le régime diffusif, la vitesse des particules existe bien mais qu'elle fluctue rapidement (dans une formulation continue, elle est représentée par un bruit blanc Gaussien) [116].

La difficulté qui se présente alors est de sélectionner un algorithme pertinent pour la détection des collisions entre particules. En effet, contrairement aux simulations PR-DNS ou PP-DNS où les particules ont des trajectoires rectilignes durant chaque petit pas de temps, la trajectoire des particules en régime diffusif n'est pas connue (cf. Fig. 3.4a). Ainsi, il n'est pas possible d'avoir recours aux techniques usuelles fondées sur le recouvrement des particules (que ce soit au début ou à la fin du pas de temps) ou sur des critères géométriques (croisement de trajectoires rectilignes entre les positions de départ et d'arrivée durant chaque pas de temps) [142]. Dans ce cas, la bonne façon de procéder est d'utiliser une approche stochastique pour déterminer la probabilité que deux particules se rencontrent durant le pas de temps. Pour ce faire, j'ai utilisé l'approche
initiée dans [104] : l'idée est de calculer la probabilité de rencontre entre deux particules en connaissant leurs positions initiales et finales et en utilisant le fait que leur mouvement relatif est une diffusion pure avec pour coefficient $B_{x,ij}^2 = B_{x,i}^2 + B_{x,j}^2$ [104, 63]. Dans ce cas, la probabilité de collision est donnée par la probabilité $P_a^b(R_{ij})$ que la distance entre cette paire atteigne une valeur inférieure ou égale à la somme des deux rayons des particules ($r_{ij} = r_i + r_j$) pendant un pas de temps. Cette probabilité est calculée à l'aide d'un pont diffusif [104] :

$$P_a^b(r_{ij}) = \frac{\exp\left(\frac{2r_{ij}(a+b-r_{ij})}{B_{x,ij}^2 \Delta t}\right) - 1}{\exp\left(\frac{2ab}{B_{x,ij}^2 \Delta t}\right) - 1} \quad (3.2)$$

où a (resp. b) est la distance entre les deux particules au début du pas de temps (resp. à la fin du pas de temps).

<figure><img src="image_32.png" /><figcaption>FIGURE 3.4 – Illustration de l'approche probabiliste utilisée pour détecter et traiter la collision entre des particules dans le régime diffusif.</figcaption></figure>

Dans l'article E.1, la méthode a été généralisée pour modifier la position des particules lorsqu'une collision a lieu. Pour ce faire, il est nécessaire d'introduire deux changements. Premièrement, le temps et la position de la collision doivent être estimés. Cela est possible grâce à une formule similaire à l'équation (3.2) qui donne la distribution des temps de première collision connaissant la distance initiale a et finale b entre une paire de particules purement diffusives. Deuxièmement, il faut modifier les trajectoires des particules suite à la collision. Dans ces travaux, on a supposé les collisions comme étant purement élastiques, ce qui implique que la dynamique du barycentre des deux particules n'est pas modifiée par la collision (puisqu'il y a conservation de l'énergie cinétique et de la quantité de mouvement). Ainsi, il suffit de modifier la trajectoire d'un des deux partenaires (en réitérant la dynamique sur le temps restant après la collision) pour en déduire la position du partenaire à la fin du pas de temps (cf. Fig. 3.4a).

Enfin, dans l'article E.2, la méthode a été étendue pour prendre en compte les phénomènes d'adhésion entre les particules. Ainsi, en utilisant le même genre de critère énergétique que dans l'algorithme MCAC (mais cette fois-ci basé sur la théorie DLVO pour les nanoparticules dans un liquide, c'est-à-dire incluant les forces de van der Waals et la double-couche électrostatique), la collision peut soit donner lieu à un rebond élastique entre les particules (cf. Fig 3.4a) soit
à la formation d'un agrégat (cf. Fig. 3.4b). Lorsqu'un agrégat est formé, la position du centre de masse de l'agrégat est calculée en partant de la position de la collision et en intégrant la dynamique de l'agrégat (c'est-à-dire d'un objet supposé sphérique ayant une masse égale à la somme des masses des deux particules impliquées dans la collision).

**Résultats** : Dans l'article E.1, nous avons caractérisé la fréquence de collision prédite par cette approche. La Figure 3.5a montre les résultats obtenus pour la fréquence de collision dans le cas d'une population bi-dispersée de particules (avec des tailles de 0.1 µm et de 1 µm) suivant un mouvement Brownien : les résultats convergent vers les valeurs théoriques pour de telles particules Browniennes à mesure que le nombre de collisions augmente (convergence statistique des approches Monte-Carlo). Néanmoins, pour obtenir de tels résultats, il est nécessaire de respecter deux conditions : le pas de temps doit être suffisamment grand pour que le régime de diffusion soit respecté ($\Delta t \gg \tau_p$) et pour que le déplacement moyen d'une particule pendant le pas de temps soit beaucoup plus grand que le diamètre des particules. Lorsque la $2^{nde}$ condition n'est pas respectée, des collisions multiples ont lieu entre une même paire de particule, puisque ces dernières se retrouvent très proches les unes des autres au début de chaque pas de temps, impliquant que leur probabilité de collision est très élevée.

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

(a) Fréquence de collisions en fonction du temps de simulation pour une population bi-disperse de particules (de 0.1 µm et 1 µm) ayant un mouvement Brownien à température ambiante. Compraison entre les résultats de l'approche probabiliste (symboles) et les valeurs théoriques (lignes) pour les trois types de collision : $1\mu m - 1\mu m (\times)$, $0.1\mu m - 0.1\mu m (+)$, $1\mu m - 0.1\mu m (o)$.

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

(b) Évolution de la concentration en particules primaires $n_1(\Delta)$, des doublets $n_2(\nabla)$, des triplets $n_3(\Box)$ et du nombre total d'agrégats/particules $n_T(o)$ en fonction du temps adimensionné $t/\tau_a$ dans le cas d'agglomération rapide (efficacité de collision $\alpha_{ij} = 1$) : comparaison entre les valeurs théoriques (lignes) et les simulations avec le modèle probabiliste (symboles).

**FIGURE 3.5** – Résultats obtenus avec l'approche probabiliste pour la fréquence de collision (à gauche) et l'évolution des tailles des agrégats (à droite).

Dans l'article E.2, l'évolution temporelle de la taille des agrégats formés a été étudiée. Pour ce faire, nous avons considéré le cas de particules ayant un mouvement Brownien à température ambiante et soumises à une agglomération rapide (c'est-à-dire que toutes les collisions donnent lieu à de l'adhésion, soit une efficacité de collision $\alpha_{ij} = 1$). La Figure 3.5b présente l'évolution du nombre de particules primaires, du nombre de doublets (agrégats formés de deux grains individuels), du nombre de triplets (agrégats formés par trois particules) et du nombre total
de particules/agrégats dans le domaine. Ces concentrations sont tracées en fonction du temps adimensionné par la demi-vie du processus d'agglomération $\tau_a = 3\mu/(4k_B T n_0)$ (qui correspond au temps au bout duquel le nombre de particules total a été diminué de moitié). Il est évident sur cette figure que les résultats obtenus avec cette approche probabiliste sont en très bon accord avec les valeurs théoriques, même en partant d'un nombre relativement limité de particules initialement dans le domaine (ici 1000).

Cela confirme l'intérêt de ce genre d'approches probabilistes pour permettre de calculer les fréquences et efficacités de collision pour des particules dans le régime diffusif. Néanmoins, pour rendre ces approches compatibles avec les modèles Lagrangiens stochastiques utilisés dans le Chapitre 2, il faudrait généraliser ces approches en sortant du cas sur-amorti. Cela n'est malheureusement pas chose facile (j'y reviendrai dans les perspectives 3.4) et cela nous a amené à formuler un autre type de modèle pour traiter l'agglomération dans le cadre de ces modèles Lagrangiens stochastiques (cf. Section 3.3.2).

## 3.3 Modèles requérant une fréquence de collision

La seconde catégorie de modèles auxquels j'ai contribué consiste à reproduire l'effet de l'agglomération sur la taille des agrégats formés en supposant la fréquence et l'efficacité des collisions apriori connues. Ce type d'approche peut dont directement bénéficier de couplages avec les approches décrites dans la Section 3.2 qui peuvent lui fournir des informations sur ces fréquences et efficacités de collision. En ce sens, cela illustre l'intérêt de réaliser des études avec différents niveaux de description, parfois en allant même jusqu'à les coupler étroitement dans ce que l'on appelle des "workflows" (où les simulations numériques sont couplées et lancées de manière automatique par un programme informatique).

### 3.3.1 Dynamique de population (PBE)

Strictement parlant, je n'ai pas développé de modèles spécifiques pour les équations de dynamique de population (abrégées en PBE pour "Population Balance Equation" en anglais). Néanmoins, je m'en suis inspiré pour développer une approche stochastique prenant en compte l'effet de l'agglomération sur l'évolution des tailles de particules. J'ai aussi utilisé ces approches PBE pour obtenir des formules analytiques afin de valider ce modèle (cf. article E.3). Pour ces raisons, je rappelle brièvement ci-dessous les grands principes de ces approches PBE ainsi que les limitations (les lecteurs intéressés par plus de détails sont invités à consulter les revues et livres sur ce sujet [126, 127, 131, 128, 70]).

**Méthodologie** : Les équations de dynamique de population sont dérivées de l'approche initiée par Smoluchowski [144], qui décrit l'évolution temporelle de la densité (en nombre) d'une population de particules dans un volume de contrôle $V_C$. Au fur et à mesure que les particules s'agglomèrent, la distribution des tailles change en suivant l'équation suivante [126] :

$$ \frac{\partial n(\mathcal{V}_i, t)}{\partial t} = \frac{1}{2} \int_0^{\mathcal{V}_i} \alpha(\mathcal{V}_i - \mathcal{V}_j, \mathcal{V}_j) \beta(\mathcal{V}_i - \mathcal{V}_j, \mathcal{V}_j) n(\mathcal{V}_i - \mathcal{V}_j, t) n(\mathcal{V}_j, t) d\mathcal{V}_i - \int_0^{\infty} \alpha(\mathcal{V}_i, \mathcal{V}_j) \beta(\mathcal{V}_i, \mathcal{V}_j) n(\mathcal{V}_i) n(\mathcal{V}_j) d\mathcal{V}_j, \quad (3.3) $$