# Modélisation numérique des matériaux granulaires immergés: applications aux avalanches sous-marines

Patrick Mutabaruka, Jean-Yves Delenne, Farhang Radjai, Soga Kenichi

► To cite this version:

Patrick Mutabaruka, Jean-Yves Delenne, Farhang Radjai, Soga Kenichi. Modélisation numérique des matériaux granulaires immergés: applications aux avalanches sous-marines. CFM 2011 - 20ème Congrès Français de Mécanique, Aug 2011, Besançon, France. hal-03422774

HAL Id: hal-03422774

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

Submitted on 9 Nov 2021

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 CC BY 4.0 - Attribution - International License
# Modélisation numérique des matériaux granulaires immergés : applications aux avalanches sous-marines

P Mutabaruka<sup>a</sup>, J-Y Delenne<sup>a</sup>, F Radjai<sup>a</sup>, S Kenichi<sup>b</sup>

a. Laboratoire de Mécanique et Génie Civil, Université de Montpellier 2, CNRS, Place Eugène Bataillon, 34095 MONPELLIER Cedex 5

b. University of Cambridge, Department of Engineering, Trumpington Street, Cambridge CB2 1PZ (UK)

## Résumé :

Les avalanches sous-marines apparaissant au voisinage des côtes constituent un risque majeur pour les structures offshores. Les mécanismes mis en jeu lors de leur initiation et propagation restent mal connus. Ils dépendent notamment de la pente, de la composition et de la quantité de matériau déstabilisé. Dans cette étude on aborde leur modélisation par une approche couplée entre la méthode de Dynamique des Contacts pour les interactions entre grains et Lattice Boltzmann pour les interactions fluide-grains. On s'intéresse notamment aux distances de propagation en fonction de la géométrie initiale.

## Abstract :

Underwater avalanches are a major risk to offshore structures. The mechanisms involved in their initiation and propagation need to be further investigated. They depend mainly on the slope, composition and quantity of material destabilized. In this study we approach their modeling by coupling the Contact Dynamics method for interactions between grains and the Lattice Boltzmann for fluid-grain interactions. We will mainly study the run out distances versus the initial geometry.

**Mots clefs :** Avalanches immergées, milieux granulaire, approches numériques discrètes

## 1 Introduction

Les avalanches sous-marines peuvent se produire sur de très faibles pentes, se propager sur de longues distances et mettre en jeu des volumes très importants de dépôts de sédiments [1]. Elles sont fréquentes au niveau des marges continentales, transportant les sédiments vers les bassins océaniques. Elles provoquent des courants de turbidité générant des dépôts dont la structure interne et la distance de propagation sont intimement liée aux effets hydrodynamiques et à la composition (granulométrie, cohésion...) des matériaux transportés [2].

Ces avalanches représentent un risque majeur pour les structures offshore et onshore (plate-forme pétrolières, pipelines, câbles sous-marins de télécommunications...). L'évaluation des risques passe par la détermination des distances de propagation en fonction de scénarios impliquant volumes, énergie gravitaire et topographie des fonts marins.

Dans cette étude on aborde la modélisation numérique de l'étalement de colonnes granulaires immergées par une approche couplée fluide-grains. Le comportement mécanique des grains est simulé par la méthode de Dynamique des Contacts [3, 4] alors que celui du fluide est décrit par la méthode Lattice Boltzmann [5, 6]. Une trentaine d'essais sont réalisés pour différentes largeurs et hauteur d'échantillons pour lesquels on s'intéresse à la détermination de la distance d'arrêt.
## 2 Méthode numérique

Différentes approches peuvent être utilisées pour étudier le couplage entre un milieu granulaire et un fluide. Concernant le fluide on peut distinguer les approches basées sur la loi de Darcy utilisant un maillage défini à une échelle supérieure au grains [7] et les approches qui s’appuient sur la résolution des équations de Navier-Stokes et qui nécessitent une discrétisation à une échelle subparticulaire. Ces dernières ont été mises en oeuvre en utilisant différentes méthodes pour le fluide telles que la méthode DNS [8, 9] ou la méthode Lattice Boltzmann (LBM) [10, 11] dans le cadre de simulations couplées avec une approche de type Discrete Element Method (DEM) pour la phase granulaire. Dans cette étude, on a choisi d’utiliser une approche LBM couplée à la méthode de Dynamique des Contacts (DC).

### 2.1 Dynamique des contacts

Dans la méthode DC, on considère que les grains sont parfaitement rigides. Les forces de contact sont déterminées à partir de lois ne faisant intervenir que les vitesses des particules en contact. La déformation du système n’étant observée qu’à partir des mouvements des particules, la discrétisation temporelle est adaptée à la dynamique des particules et ne fait pas intervenir l’échelle des contacts comme c’est le cas pour les approches en éléments discrets classiques telles que la dynamique moléculaire (DM). Ceci permet notamment d’avoir des pas de temps $\Delta t \simeq 10^{-4}$ bien supérieurs à ceux de la DM.

Dans cette étude, on utilise la librairie de calcul gdm-tk développée par V. Richefeu, C. Voivret et F. Radjaï (www.cgp-gateway.org). Des lois de contact de Signorini et de frottement de Coulomb sont implémentées avec un coefficient de frottement $\mu = 0,3$ et un coefficient de restitution faible.

### 2.2 Lattice Boltzmann

La méthode Lattice Boltzmann (LBM) a été développée initialement pour la résolution de problème de dynamique des gaz [12, 14]. Elle a été enrichie par de nombreux travaux qui ont permis d’étendre son domaine d’utilisation à des problèmes complexes d’hydrodynamique mais aussi de transports de solutés, d’écoulements biphasiques...

Cette méthode découle de l’équation de Boltzmann discrétisée sur réseau qui comporte deux termes principaux décrivant la collision et l’advection des molécules. La complexité de leur dynamique est prise en compte en considérant la densité de probabilité $f$ (fonction de distribution) d’avoir une molécule à une position $\vec{x}$ à l’instant $t$ avec une quantité de mouvement $\vec{p}$. En supposant notamment que les interactions moléculaires ont une portée faible, il est possible de considérer que le milieu est composé de particules de fluide mésoscopiques dont la dynamique est gouvernée par les équations classiques de Newton :

$$ \frac{d\vec{x}}{dt} = \frac{\vec{p}}{m} \qquad (1) $$

$$ \frac{d\vec{p}}{dt} = \vec{F} \qquad (2) $$

où $m$ est la masse d’une particule et $\vec{F}$ représente la force extérieure appliquée à la particule.

La description de la cinétique des particules à partir de la fonction de distribution peut s’écrire sous la forme :

$$ D_t f(\vec{x}, \vec{p}, t) = C \qquad (3) $$

avec $D_t = [\partial_t + \frac{\vec{p}}{m} \partial_{\vec{x}} + \vec{F} \partial_{\vec{p}}]$, et où le terme de gauche décrit l’advection de la particule et le terme de droite est un opérateur décrivant les collisions entre particules. Dans le cas des gaz, ce terme peut être établi à partir de l’échelle moléculaire en supposant que le système n’implique que des interactions binaires entre molécules. Différent travaux théoriques ont permis de généraliser aux cas de collisions simultanées permettant notamment d’étendre la méthode au cas des fluides denses. La résolution du terme de droite peut être complexe et conduire à un système impliquant plusieurs paramètres de relaxation. Un modèle simple, classiquement utilisé, est celui de Bhatnagar-Gross-Krook (BGK) [12]. Ce modèle est basé sur le fait que pour un fluide la présence de collisions est de nature à faire tendre le
système vers un équilibre local. Dans ce modèle un seul paramètre de relaxation $\tau$ est pris en compte et l'opérateur de collision peut se mettre sous la forme :

$$C = -\frac{1}{\tau}[f(\vec{x}, t) - f^{eq}(\vec{x}, t)] \quad (4)$$

où $f^{eq}(\vec{x}, t)$ est la fonction d'équilibre de Maxwell.

Ce modèle permet de retrouver les équations de Navier-Stokes avec une bonne approximation. Différents schémas de discrétisation sont possibles [13]. Le passage des densités locales aux variables mésoscopique (masse volumique, pression, vitesse) est obtenu par intégration des densités $f$.

Par la suite, le domaine fluide sera discrétisé selon un réseau régulier de pas $\Delta x$. Un schéma classique D2Q9 de discrétisation sur réseau bidimensionnel est utilisé [14]. Au niveau de chaque noeud du maillage, ce schéma permet de prendre en compte la probabilité que les particules restent sur place ou se propagent dans 8 directions de l'espace. La vitesse de propagation des ondes sur le réseau est donnée par $c_s = \frac{1}{\sqrt{3}} \frac{\Delta x}{\Delta t}$. Avec ces paramètres la viscosité cinématique est donnée par :

$$\nu = \frac{1}{\sqrt{3}} c_s \Delta x \left( \tau - \frac{1}{2} \right) \quad (5)$$

## 2.3 Couplage Dynamique des contacts – Lattice Boltzmann

Le couplage entre la méthode DC et la méthode LBM se fait par l'intermédiaire du calcul des efforts au niveau des frontières grain-liquide [11]. Une des difficultés majeures provient du fait qu'il est nécessaire de déterminer les intersections entre le réseau maillé et les grains. Différentes approches permettent de gérer ces conditions aux limites. Une approche intéressante, proposée par Bouzidi et al. [15] permet de limiter les effets de gradin dus à la régularité du maillage. En pratique, des résultats de bonne qualité sont obtenus pour une discrétisation de 10 mailles ($10 \Delta x$) selon le diamètre des grains. Une contrainte importante provient du fait qu'en 2D l'espace poral dans des échantillons granulaires denses est occlus. Ceci conduit à l'impossibilité d'avoir des écoulements entre grains. Pour remédier à ce problème on considère qu'une couche de liquide de 1 maille entoure les grains. Cette hypothèse qui a pour conséquence de modifier les forces de lubrification entre grains fait l'objet d'une étude en cours.

Le pas de temps pris en compte lors du couplage CD-LBM, est imposé par le plus petit pas de temps entre celui de la l'algorithmme CD et celui de LBM. Pour des viscosités prochent de celle de l'eau, la discrétisation temporelle de l'algorithmme LBM est inférieur à celle de l'algorithmme CD mais reste supérieur à celle de l'algorithmme DM. La durée de simulation totale étant gouvrnée par le temps de calcul de l'algorithmme LBM, le couplage CD-LBM et avantageux devant celui du DM-LBM. Cependant, il est possible de synchroniser le calcul DM-LBM en gardant deux discrétisations temporelles de sorte que pour n pas de temps DM on ait un seul pas de temps LBM ($n = \Delta t_{LBM} / \Delta t_{DM}$, avec n arrondi à l'entier supérieur). La contrainte de cette synchronisation, est de conduire à une évolution en éscalier des forces hydrodynamiques et donc une discontinuitée dans le calcul de l'algorithmme DM. Ceci peut engendrer une divergence ou des interpénétrations importantes au niveau des grains; ce qui rendrait le système instable. Ce problème de synchronisation donne l'avantage à la méthode CD devant la méthode DM où, la discrétisation temporelle peut prendre la même valeur sans générer une influence considéreble sur le temps de calcul.

La discrétisation temporelle de l'algorithmme LBM peut être obtenu en imposant un nombre de Mach suffisamment petit pour que les ondes se propagent à une vitesse au moins 10 fois supérieure à la vitesse des grains et en utilisant l'équation 5 on déduit la valeur du temps de relaxation $\tau$.

## 3 Etalement de colonnes immergées

Dans le cas sec, les instabilités gravitaires ont fait l'objet de nombreuses études expérimentales et numériques depuis plusieurs années. La statistique des tailles des avalanches est contrôlée par la différence entre les angles de déclenchement et d'arrêt d'un talus granulaire. Plus récemment, des
études ont été menées sur l'étalement d'une masse granulaire sous effet de la gravité. Des relations non trivialles ont été mises en évidence entre l'étalement de l'échantillon et sa géométrie initiale [16, 17, 18]. Cependant, la pertinence de ces résultats en présence d'un fluide reste très peu explorée. Dans le cadre de cet article on s'intéresse à des configurations géométriques similaires en situation immergée dans un cas bidimensionnel (Fig. 1).

<figure><img src="image_2.png" /><figcaption>FIGURE 1 – Etalement d'une colonne granulaire.</figcaption></figure>

## 3.1 Protocole

Des échantillons granulaires de différentes largeurs $L_0$ et hauteur $H_0$ sont réalisés en déposant les grains par gravité dans des boites. Les diamètres $d$ des grains sont choisis en utilisant un modèle basé sur la loi statistique beta cumulé $\beta(d_r, a, b)$ proposé par Voivret et al. [19] où $d_r = (d-d_{min})/(d_{max}-d_{min})$ et $a = 2$ et $b = 2$. L'étalement de la courbe granulométrique est donné par $s = (d_{max}-d_{min})/(d_{max}+d_{min})$ et est égal à $1/3$ et le rapport entre le diamètre maximamal et le diamètre minimal $r = d_{max}/d_{min}$ est de 2. Pour les différents échantillons réalisés la masse volumique des grains $\rho_s$ est choisie égale à $2600\text{kg/m}^3$ et le diamètre moyen est $<d> \simeq 2 \cdot 10^{-3}m$.

Trois séries $S^1$, $S^2$, $S^3$ sont réalisées avec des largeurs $L_0^1 = 6 <d>$, $L_0^2 = 2L_0^1$, $L_0^3 = 3L_0^1$. Pour chaque série différents échantillons sont obtenus en faisant varier les hauteurs $H_0$ de $0.33L_0$ jusqu'à $4.33L_0$. Par la suite on s'intéresse au rapport d'aspect de l'échantillon qui est donné par $a = H_0/L_0$.

## 3.2 Résultats

Une trentaine d'essais sont réalisés pour les échantillons des 3 séries. Les compacités initiales des différents échantillons sont égales à $\phi = 0.80$ ce qui est très proche de la compacité maximale accessible.

<figure><img src="image_3.png" /><figcaption>FIGURE 2 – Etalement d'une colonne de grains dans de l'eau. Les grandeurs sont croissantes du bleu au rouge. Sur la figure sont représentés : dans le fluide la norme des vitesses, dans les grains la pression exercée par les grains en contacts environnants.</figcaption></figure>

Pour chaque échantillon on note $L(t)$ la distance d'étalement de la colonne en fonction du temps. Cette distance est déterminée en considérant la position du dernier contact pour lequel les grains ne sont pas isolés du reste de l'échantillon. Un exemple de simulation d'étalement est donné sur la figure 2. La figure 3a représente les cinétiques d'étalement des échantillons de la série $S^2$. On constate une
première phase pour laquelle l'évolution du front est rapide puis une saturation qui apparaît pour
des temps légèrement croissants en fonction de a. L'énergie cinétique totale de l'échantillon E_c est la
somme de toutes les énergies cinétiques de translation et de rotation (figure 3b). La distance d'arrêt
est déterminée quand l'énergie devient inférieure à 1% de l'énergie au pic.

<figure><img src="image_4.png" /><figcaption>FIGURE 3 – (a) Evolution du front de propagation $L(t)$ pour différents rapport d'aspects $a$ dans le cas de la série S². (b) Evolution de l'énergie cinétique totale de l'échantillon $E_c$ en fonction du temps</figcaption></figure>

La figure 4 montre l'évolution de la distance d'étalement adimensionnée $(L - L_0)/L_0$ en fonction de la surface initiale de l'échantillon $L_0 \times H_0$ et du rapport d'aspect initial $a$. On constate que, comme dans le cas avalanches gravitaires aériennes, la distance d'étalement se met à l'échelle avec $a$. Une analyse dimensionnelle basée sur [20, 21] permet de montrer qu'avec les paramètres utilisés pour les grains et pour la viscosité de l'eau, les forces inertielles subies par les grains dominent les forces de viscosité. Pour des tailles de grains plus faibles, les forces de viscosité agissant sur les grains sont plus importantes. Des essais en cours permettront de préciser si l'évolution de la distance d'arrêt en fonction de $a$ présente ou non plusieurs régimes. Ce point est crucial pour mieux comprendre la transition entre les courants de turbidité et les coulées de débris et les phénomènes de ségrégation complexes que l'on rencontre au niveau du front de propagation des marines.

<figure><img src="image_5.png" /><figcaption>FIGURE 4 – Distance d'étalement adimensionnée en fonction de : (a) la surface de l'échantillon ; (b) le rapport d'aspect $a$. Séries S¹ (o), S² (+), S³ (□).</figcaption></figure>

4 Conclusions

Dans cette étude on s'est intéressé aux couplages entre un milieu granulaire 2D simulé par la méthode Dynamique des Contacts et un fluide modélisé par une approche Lattice Boltzmann. Des essais d'effondrement de colonnes immergées ont été réalisés pour différentes largeurs et hauteurs initiales d'échantillon. On montre que, comme dans le cas des avalanches aériennes, la distance d'étalement adimensionnée de ces colonnes se met à l'échelle selon une loi de puissance avec le rapport d'aspect initial de l'échantillon et non sa surface.
Une perspective directe de cette étude sera de préciser si pour des tailles de grain inférieures, les effets de viscosité du fluide s'expriment différemment. Enfin, des simulations 3D devront permettre de confirmer cette analyse.

Références

[1] Legros, F. 2002 The mobility of long-runout landslides. *Engineering Geology*. **63** 301-331

[2] Bouma, A. H. 1962 Sedimentology of some Flysch deposits : A graphic approach to facies interpretation. *Elsevier Pub. Co.*, New York.

[3] Moreau, J. J. 1994 Some numerical methods in multibody dynamics : application to granular materials. *European journal of mechanics A. Solids*. **13** 93-114

[4] Radjai, F., Richefeu, V. 2009 Contact dynamics as a nonsmooth discrete element method. *Mechanics of Materials*. **41** 715-728

[5] Mansouri, M., Delenne, J.-Y., El Youssoufi, M. S., Seridi, A. 2010 A 3D DEM-LBM approach for the assessment of the quick condition for sands. *C.R. Mecanique*. **337** 675-681

[6] Lallemand, P., Luo, L.-S. 2003 Lattice Boltzmann method for moving boundaries. *Journal of Computational Physics*. **184** 406-421

[7] McNamara, S., Flekkoy, E. G., Maloy, K. J. 1999 Grains and gas flow : Molecular dynamics with hydrodynamic interactions. *Physical Review E*. **62** 4

[8] Wachs, A. 2009 A DEM-DLM/FD method for direct numerical simulation of particulate flows : Sedimentation of polygonal isometric particles in a Newtonian fluid with collisions. *Computers and Fluids*. **38** 1608-1628

[9] Topin, V., Dubois, F., Monerie, Y., Perales, F., Wachs, A. 2011 Micro-rheology of dense particulate flow : application to immersed avalanches. *Journal of Non-Newtonian Fluid Mechanics*. **166** 63-72

[10] Feng, Z.-G., Michaelides, E.E 2004 The immersed boundary-lattice Boltzmann method for solving fluid-particles interaction problems. *J. Computational Physics*. **195** 602-628

[11] Ladd, A. J. C. 1994 Numerical simulation of particular suspensions via a discretized Boltzmann equation, Part 2, Numerical results. *J. Fluid Mech*. **271** 311-339

[12] Wolf-Gladrow, D. A. 2000 Lattice-Gas Cellular Automata and Lattice Boltzmann models : an introduction, *Lecture Notes in Mathematics*, Springer-Verlag Press, Berlin Heidelberg New York.

[13] Chen, S. and Doolen, G. D. 1998 Lattice Boltzmann Method for fluid flows. *Annu. Rev. Fluid Mech*. **30** 329-64

[14] Sukop, M. C., Thorne, D. T. 2006 Lattice Boltzmann Modeling. *An Introduction for Geoscientists and Engineers*. Springer-Verlag Press, Berlin Heidelberg New York.

[15] Bouzidi, M., Firdaouss, M., Lallemand, P. 2001 Momentum transfer of a Boltzmann-lattice fluid with boundaries. *Physics of Fluids*. **13** 3452-3459

[16] Lajeunesse, E., Monnier, J. B., Homsy, G. M. 2005 Granular slumping on a horizontal surface. *Physics of fluids* **17** 103302

[17] Mangeney-Castelnau, A., Bouchut, F., Vilotte, J. P., Lajeunesse, E., Aubertin, A., Pirulli, M. 2005 On the use of Saint Venant equations to simulate the spreading of a granular mass. *J. Geophys. Res.* **110**

[18] Staron, L., Hinch, E. J. 2007 The spreding of a granular mass : role of grain properties and intial condition. *Granular Matter*. **9** 205-217

[19] Voivret, C., Radjai, F., Delenne, J.-Y., El Youssoufi, M. S. 2007 Space-filling properties of poly-disperse granular media. *Physical Review E*. APS, **76** 021301

[20] Courrech du Pont, S., Gondret, P., Perrin, B., Rabaud, M. 2003 Granular Avalanches in Fluids. *Physical Review Letters*. **90**

[21] Cassar, C., Nicolas, M., Pouliquen, O. 2005 Submarine granular flows down inclined planes. *Physics of fluids*. **17** 103301