# IRN GeoMech workshop on hydromechanical instabilities. Booklet of abstracts

Antoine Wautier, Nadia Benahmed, Pierre Philippe

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

## To cite this version:

Antoine Wautier, Nadia Benahmed, Pierre Philippe. IRN GeoMech workshop on hydromechanical instabilities. Booklet of abstracts. IRN GeoMech/INRAE, 2023. hal-04457745

HAL Id: hal-04457745

https://hal.inrae.fr/hal-04457745

Submitted on 14 Feb 2024

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_2.png" /></figure>
IRN GeoMech workshop on hydromechanical
instabilities. Booklet of abstracts

Antoine Wautier, Nadia Benahmed, Pierre Philippe

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

To cite this version:

Antoine Wautier, Nadia Benahmed, Pierre Philippe. IRN GeoMech workshop on hydromechanical instabilities. Booklet of abstracts. IRN GeoMech/INRAE, 2023. hal-04457745

HAL Id: hal-04457745

https://hal.inrae.fr/hal-04457745v1

Submitted on 14 Feb 2024

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_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_3.png" /></figure>

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

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

# IRN GeoMech workshop on hydromechanical instabilities

INRAE Aix-en-Provence, May 2-3 2023

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

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

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

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

Organizing committee: Antoine Wautier, Nadia Benahmed, Pierre Philippe

<figure><img src="image_10.png" /></figure>
<figure><img src="image_3.png" /></figure>

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

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

# IRN GeoMech workshop on
hydromechanical instabilities

INRAE Aix-en-Provence, May 2-3 2023

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

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

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

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

Organizing committee: Antoine Wautier, Nadia Benahmed, Pierre Philippe

<figure><img src="image_10.png" /></figure>
# Welcome Address

Dear all,

It is our pleasure to welcome you in Aix-en-Provence for this international workshop co-organized by the international research network GeoMech and the research department AQUA of INRAE.

Granular materials are involved in many natural hazards, such as landslides, avalanches, dike or dam failures, etc. Due to their discrete nature, they have a complex mechanical behavior resulting from local interactions between grains and collective behaviors organized at different scales. Thus, the failure modes of geomaterials can take several forms: diffuse failure possibly resulting in liquefaction, localized failure with shear bands at several scales (from a few grains to the size of the complete system), mixed mode failure, etc. Moreover, granular materials being porous, they are very often subjected to water infiltration likely to induce important changes of microstructure that can affect their physical, hydraulic and mechanical properties.

The main objective of this workshop is to review recent advances in the understanding of the elementary mechanisms of destabilization of granular materials and their impact on failure modes (e.g. liquefaction, strain localization). The topics addressed during the workshop may be related, for instance, to the mechanical response and stability of geomaterials in the presence of capillary or solid bridges (sintering, bio-calcification, dissolution/precipitation, ...), internal erosion (suffusion and clogging), surface erosion, etc.

We hope that the list of forthcoming talks will lead to fruitful discussion and possibly to future collaborations.

Please sit back, relax and enjoy the workshop!

Antoine Wautier and Nadia Benahmed and Pierre Philippe
# Welcome Address

Dear all,

It is our pleasure to welcome you in Aix-en-Provence for this international workshop co-organized by the international research network GeoMech and the research department AQUA of INRAE.

Granular materials are involved in many natural hazards, such as landslides, avalanches, dike or dam failures, etc. Due to their discrete nature, they have a complex mechanical behavior resulting from local interactions between grains and collective behaviors organized at different scales. Thus, the failure modes of geomaterials can take several forms: diffuse failure possibly resulting in liquefaction, localized failure with shear bands at several scales (from a few grains to the size of the complete system), mixed mode failure, etc. Moreover, granular materials being porous, they are very often subjected to water infiltration likely to induce important changes of microstructure that can affect their physical, hydraulic and mechanical properties.

The main objective of this workshop is to review recent advances in the understanding of the elementary mechanisms of destabilization of granular materials and their impact on failure modes (e.g. liquefaction, strain localization). The topics addressed during the workshop may be related, for instance, to the mechanical response and stability of geomaterials in the presence of capillary or solid bridges (sintering, bio-calcification, dissolution/precipitation, ...), internal erosion (suffusion and clogging), surface erosion, etc.

We hope that the list of forthcoming talks will lead to fruitful discussion and possibly to future collaborations.

Please sit back, relax and enjoy the workshop!

Antoine Wautier and Nadia Benahmed and Pierre Philippe
Practical information

To reach Aix-en-Provence by train the most convenient is to stop at Aix-en-Provence TGV and then take the bus 40 to the city center. Participants coming by plane can land at Marseilles Airport and also reach Aix-en-Provence by bus 40. From Aix-en-Provence, to reach the conference premises, take bus 13 (or bus 110) to Le Tholonet (stop at Ferrageon).

During the two days of the workshop, lectures will be given in the Conference Room Cézanne (building D) at INRAE research center of Aix-en-Provence (address: 3275 route Cézanne, 13100 Le Tholonet). To access to internet, select the Wifi "Eduspot". It should automatically oppen a connection page (https://portail.invite.inrae.fr). Click on "invité à un séminaire", put your email address with the password "ZQHEZU".

Lunches will be served just outside the conference room in the form of a lunch cocktail by "La table de Cana". Vegetarian options will be available. Coffee will be served directly next to the conference room.

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

For those of you needing an accommodation during your stay in Aix-en-Provence, you will be hosted at the **Hotel Odalys City l'Atrium** close to the city center (15 Cr Gambetta, 13100 Aix-en-Provence).
# Practical information

To reach Aix-en-Provence by train the most convenient is to stop at Aix-en-Provence TGV and then take the bus 40 to the city center. Participants coming by plane can land at Marseilles Airport and also reach Aix-en-Provence by bus 40. From Aix-en-Provence, to reach the conference premises, take bus 13 (or bus 110) to Le Tholonet (stop at Ferrageon).

During the two days of the workshop, lectures will be given in the Conference Room Cézanne (building D) at INRAE research center of Aix-en-Provence (address: 3275 route Cézanne, 13100 Le Tholonet). To access to internet, select the Wifi "Eduspot". It should automatically oppen a connection page (https://portail.invite.inrae.fr). Click on "invité à un séminaire", put your email address with the password "ZQHEZU".

Lunches will be served just outside the conference room in the form of a lunch cocktail by "La table de Cana". Vegetarian options will be available. Coffee will be served directly next to the conference room.

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

For those of you needing an accommodation during your stay in Aix-en-Provence, you will be hosted at the **Hotel Odalys City l'Atrium** close to the city center (15 Cr Gambetta, 13100 Aix-en-Provence).
We are pleased to invite you to a nice dinner in town at the restaurant **Drôle d'endroit** (14 Rue Annonerie Vieille, 13100 Aix-en-Provence). We will wait for you at the restaurant on Tuesday night by 7.30 pm.

The restaurant is located within walking distance from the Hotel Odalys City l'Atrium. Please note that the restaurant is a little bit hidden away from the main streets of Aix city center. So you may better check the map before leaving the hotel ;)

In any case, if you need help you can give us a call :

+33 6 88 28 18 61 (Antoine)

+33 6 10 78 22 54 (Nadia)

+33 6 89 81 70 35 (Pierre).
We are pleased to invite you to a nice dinner in town at the restaurant **Drôle d'endroit** (14 Rue Annonerie Vieille, 13100 Aix-en-Provence). We will wait for you at the restaurant on Tuesday night by 7.30 pm.

The restaurant is located within walking distance from the Hotel Odalys City l'Atrium. Please note that the restaurant is a little bit hidden away from the main streets of Aix city center. So you may better check the map before leaving the hotel ;)

In any case, if you need help you can give us a call :

+33 6 88 28 18 61 (Antoine)

+33 6 10 78 22 54 (Nadia)

+33 6 89 81 70 35 (Pierre).
# Program

Tuesday 02/05/2023

<table><tr><td>10.30 - 11.00</td><td>Welcome coffee</td></tr><tr><td>11.00 - 11.10</td><td>Welcome address</td></tr><tr><td>11.10 - 11.30</td><td>Oscar POLANIA<br>Experimental study of column collapse: influence of the grain size poly-dispersity on the collapse sequence and runout</td></tr><tr><td>11.30 - 11.50</td><td>Jérôme DURIEZ<br>DEM capabilities with Polyhedral and Level Set shape descriptors</td></tr><tr><td>11.50 - 12.10</td><td>Na DENG<br>Multiscale modeling of oedometer test on crushable tubular particles</td></tr><tr><td>12.10 - 14.00</td><td>Lunch</td></tr><tr><td>14.00 - 14.40</td><td>Keynote: Pierre-Yves HICHER (online)<br>Instability in granular slopes subjected to internal erosion</td></tr><tr><td>14.40 - 15.00</td><td>Sylvie NICAISE<br>Field overflow tests on a Rhône levee</td></tr><tr><td>15.00 - 15.20</td><td>Ngoc Son NGUYEN (online)<br>Consequences of scalping and scalping/replacement procedures on strength properties of coarse-grained gap-graded soils</td></tr><tr><td>15.20 - 15.40</td><td>Abhijit HEGDE<br>An experimental study on the shear banding phenomenon in granular materials using the orthogonal cutting apparatus</td></tr><tr><td>15.40 - 16.10</td><td>Coffee break</td></tr><tr><td>16.10 - 16.30</td><td>Antoine WAUTIER<br>Extending the bifurcation domain concept at the structure scale</td></tr><tr><td>17.00 - 18.30</td><td>Free time or small walk in the Sainte-Victoire mountain area</td></tr><tr><td>19.30</td><td>Dinner in town</td></tr></table>
# Program

Tuesday 02/05/2023

<table><tr><td>10.30 - 11.00</td><td>Welcome coffee</td></tr><tr><td>11.00 - 11.10</td><td>Welcome address</td></tr><tr><td>11.10 - 11.30</td><td>Oscar POLANIA<br>Experimental study of column collapse: influence of the grain size poly-dispersity on the collapse sequence and runout</td></tr><tr><td>11.30 - 11.50</td><td>Jérôme DURIEZ<br>DEM capabilities with Polyhedral and Level Set shape descriptors</td></tr><tr><td>11.50 - 12.10</td><td>Na DENG<br>Multiscale modeling of oedometer test on crushable tubular particles</td></tr><tr><td>12.10 - 14.00</td><td>Lunch</td></tr><tr><td>14.00 - 14.40</td><td>Keynote: Pierre-Yves HICHER (online)<br>Instability in granular slopes subjected to internal erosion</td></tr><tr><td>14.40 - 15.00</td><td>Sylvie NICAISE<br>Field overflow tests on a Rhône levee</td></tr><tr><td>15.00 - 15.20</td><td>Ngoc Son NGUYEN (online)<br>Consequences of scalping and scalping/replacement procedures on strength properties of coarse-grained gap-graded soils</td></tr><tr><td>15.20 - 15.40</td><td>Abhijit HEGDE<br>An experimental study on the shear banding phenomenon in granular materials using the orthogonal cutting apparatus</td></tr><tr><td>15.40 - 16.10</td><td>Coffee break</td></tr><tr><td>16.10 - 16.30</td><td>Antoine WAUTIER<br>Extending the bifurcation domain concept at the structure scale</td></tr><tr><td>17.00 - 18.30</td><td>Free time or small walk in the Sainte-Victoire mountain area</td></tr><tr><td>19.30</td><td>Dinner in town</td></tr></table>
Wednesday 03/05/2023

<table><tbody><tr><td>09.10 - 09.30</td><td>Pierre PHILIPPE<br/>Submerged discharge of cemented granular subsoil as a sinkhole occurrence scenario in flood condition</td></tr><tr><td>09.30 - 09.50</td><td>Dat G. PHAN (online)<br/>Capturing onset of bifurcation and post-localisation behaviour in constitutive modelling of partially saturated soils</td></tr><tr><td>09.50 - 10.10</td><td>Takashi MATSUSHIMA (online)<br/>Observation of grain-scale erosion/sedimentation process using DEM/CFD coupling simulation</td></tr><tr><td>10.10 - 10.30</td><td>Marie MIOT<br/>Extending the 3D-H model to snow mechanics: model calibration using oedometer tests</td></tr><tr><td>10.30 - 11.00</td><td>Coffee break</td></tr><tr><td>11.00 - 11.40</td><td>Keynote: Didier MAROT<br/>Influence of stress state on the initiation and the whole development of suffusion</td></tr><tr><td>11.40 - 12.00</td><td>Shadi YOUSSEF<br/>Role of soil preparation in controlling the erodibility behavior of cohesive soils when subjected to surface erosion</td></tr><tr><td>12.00 - 12.20</td><td>Abbas FARHAT<br/>Experimental study of the hydromechanical instabilities of a cemented granular material</td></tr><tr><td>12.20 - 12.30</td><td>Closure and announcements</td></tr><tr><td>12.30 - 14.00</td><td>Lunch</td></tr><tr><td>14.00 - 16.00</td><td>Lab tour</td></tr></tbody></table>
Wednesday 03/05/2023

<table><tbody><tr><td>09.10 - 09.30</td><td>Pierre PHILIPPE<br>Submerged discharge of cemented granular subsoil as a sinkhole occurrence scenario in flood condition</td></tr><tr><td>09.30 - 09.50</td><td>Dat G. PHAN (online)<br>Capturing onset of bifurcation and post-localisation behaviour in constitutive modelling of partially saturated soils</td></tr><tr><td>09.50 - 10.10</td><td>Takashi MATSUSHIMA (online)<br>Observation of grain-scale erosion/sedimentation process using DEM/CFD coupling simulation</td></tr><tr><td>10.10 - 10.30</td><td>Marie MIOT<br>Extending the 3D-H model to snow mechanics: model calibration using oedometer tests</td></tr><tr><td>10.30 - 11.00</td><td>Coffee break</td></tr><tr><td>11.00 - 11.40</td><td>Keynote: Didier MAROT<br>Influence of stress state on the initiation and the whole development of suffusion</td></tr><tr><td>11.40 - 12.00</td><td>Shadi YOUSSEF<br>Role of soil preparation in controlling the erodibility behavior of cohesive soils when subjected to surface erosion</td></tr><tr><td>12.00 - 12.20</td><td>Abbas FARHAT<br>Experimental study of the hydromechanical instabilities of a cemented granular material</td></tr><tr><td>12.20 - 12.30</td><td>Closure and announcements</td></tr><tr><td>12.30 - 14.00</td><td>Lunch</td></tr><tr><td>14.00 - 16.00</td><td>Lab tour</td></tr></tbody></table>
List of extended
abstracts
List of extended
abstracts
Experimental study of column collapse: influence of the grain size polydispersity on the collapse sequence and runout

Oscar Polanía¹,², Nicolás Estrada², Mathieu Renouf¹, Emilien Azéma¹,³ and Miguel Cabrera⁴

¹ LMGC, Université de Montpellier, CNRS, Montpellier, France

² Dept. Ingeniería Civil y Ambiental, Universidad de Los Andes, Bogotá D.C., Colombia

³ Institut Universitaire de France (IUF), Paris

⁴ Department of Geoscience & Engineering, TU Delft, Delft, The Netherlands

oscar.polania@umontpellier.fr

os.polaniao@uniandes.edu.co

**Keywords:** Polydispersity, granular, flows, immersed

**Abstract**

Granular flows are found in geophysical processes like landslides, debris flows or pyroclastic flows. Those flows are often subaerial, but occasionally they can be submerged in a fluid. A characteristic among them is the presence of particles of different sizes, a property termed polydispersity. The flow of granular materials can be distinguished in three regimes, according to the motion of the grains with the Stokes number and the density ratio between the particles and fluid. The combination of these dimensionless numbers leads to three flow regimes, that are: the inertial regime, the free-fall regime, and the viscous regime [1].

A benchmark study case for the study of granular flows is the collapse of a granular column [2, 3, 4, 5, 6, 7, 8]. In this experiment, a column built up with an initial length $L_0$, an initial height $H_0$, hence an aspect ratio $A = H_0/L_0$, is let to collapse by self-weight in a horizontal plane, reaching a final runout Lf and a final height Hf. During collapse, the column front has an initial acceleration phase followed by a steady propagation stage with constant velocity, until it decelerates and stops. The final runout of granular columns is found to be longer in dry systems than in immersed systems [5, 6], and for dry cases, polydispersity does not influence the column mobility [9, 10]. However, for immersed granular columns, it has been observed that mobility increases in bidisperse columns due to a segregation layer between big and small grains, and this effect is enhanced when the ratio between big and small grains increases [11].

Although previous studies have done a great effort for understanding granular flows in the column collapse configuration, the influence that polydispersity has on it is still unclear. The role that polydispersity has on granular flows is important because polydispersity is intrinsic in many natural processes. Therefore, our research focuses on exploring the influence that polydispersity has on the collapse sequence and runout of dry and immersed granular columns in the free-fall and inertial regimes. We explore experimentally the three-dimensional nature of this process. We study short and tall columns and explore a wide size span, having systems with $\lambda = d_{max}/d_{min} = [1-20]$, where $d_{max}$ and $d_{min}$ stand for the maximum and minimum particle diameter in the system, respectively. Besides the novelty of using different GSDs, our setup is instrumented with a network of basal pore pressure transducers, allowing us to evaluate the polydispersity role on the pore pressure changes $\Delta P$ and the influence that it has on the collapse (see Figure 1).
Experimental study of column collapse: influence of the grain size polydispersity on the collapse sequence and runout

Oscar Polanía¹,², Nicolás Estrada², Mathieu Renouf¹, Emilien Azéma¹,³ and Miguel Cabrera⁴

¹ LMGC, Université de Montpellier, CNRS, Montpellier, France

² Dept. Ingeniería Civil y Ambiental, Universidad de Los Andes, Bogotá D.C., Colombia

³ Institut Universitaire de France (IUF), Paris

⁴ Department of Geoscience & Engineering, TU Delft, Delft, The Netherlands

oscar.polania@umontpellier.fr

os.polaniao@uniandes.edu.co

**Keywords:** Polydispersity, granular, flows, immersed

**Abstract**

Granular flows are found in geophysical processes like landslides, debris flows or pyroclastic flows. Those flows are often subaerial, but occasionally they can be submerged in a fluid. A characteristic among them is the presence of particles of different sizes, a property termed polydispersity. The flow of granular materials can be distinguished in three regimes, according to the motion of the grains with the Stokes number and the density ratio between the particles and fluid. The combination of these dimensionless numbers leads to three flow regimes, that are: the inertial regime, the free-fall regime, and the viscous regime [1].

A benchmark study case for the study of granular flows is the collapse of a granular column [2, 3, 4, 5, 6, 7, 8]. In this experiment, a column built up with an initial length $L_0$, an initial height $H_0$, hence an aspect ratio $A = H_0/L_0$, is let to collapse by self-weight in a horizontal plane, reaching a final runout Lf and a final height Hf. During collapse, the column front has an initial acceleration phase followed by a steady propagation stage with constant velocity, until it decelerates and stops. The final runout of granular columns is found to be longer in dry systems than in immersed systems [5, 6], and for dry cases, polydispersity does not influence the column mobility [9, 10]. However, for immersed granular columns, it has been observed that mobility increases in bidisperse columns due to a segregation layer between big and small grains, and this effect is enhanced when the ratio between big and small grains increases [11].

Although previous studies have done a great effort for understanding granular flows in the column collapse configuration, the influence that polydispersity has on it is still unclear. The role that polydispersity has on granular flows is important because polydispersity is intrinsic in many natural processes. Therefore, our research focuses on exploring the influence that polydispersity has on the collapse sequence and runout of dry and immersed granular columns in the free-fall and inertial regimes. We explore experimentally the three-dimensional nature of this process. We study short and tall columns and explore a wide size span, having systems with $\lambda = d_{max}/d_{min} = [1-20]$, where $d_{max}$ and $d_{min}$ stand for the maximum and minimum particle diameter in the system, respectively. Besides the novelty of using different GSDs, our setup is instrumented with a network of basal pore pressure transducers, allowing us to evaluate the polydispersity role on the pore pressure changes $\Delta P$ and the influence that it has on the collapse (see Figure 1).
We observe that the collapse process of polydisperse dry columns is independent of the polydispersity level. For dry monodisperse and polydisperse columns, the collapse initiates at a similar instant and evolves with a similar front velocity. When the collapse comes to an end, the runout slightly varies, being greater for polydisperse columns. On the contrary, for immersed columns, the interaction with a fluid changes the collapse of monodisperse and polydisperse columns. Moreover, the collapse initiation is delayed in polydisperse columns. Then, during the steady propagation stage, the velocity of both monodisperse and polydisperse columns is similar, but polydisperse systems keep this stage for a longer period and their runout considerably increases, compared with monodisperse (see Figure 1).

We show that the differences between monodisperse and polydisperse immersed columns arise from the pore pressure changes $\Delta P$ that columns have beneath them. For the collapse process, we identify the pore pressure evolution in three stages. First, there is a negative change of the pore pressure due to the in-flow of fluid into the column. The fluid comes in into the granular structure when the column dilates to start the motion [4]. Then, during the spreading of the granular flow, the changes of the pore pressure are positive. Finally, pore pressure changes stabilize $\Delta P \approx 0$ (see Figure 2).

Our research provides evidence that the collapse of immersed granular columns strongly depends on the polydispersity level. Our results could be of great use in processes that involve polydisperse granular materials that account for extreme grain size variations, e.g., debris flows, pyroclastic flows, submarine landslides. Our experimental work provides a novel insight of the basal pore pressure variations and the influence that it has on the collapse sequence of immersed flows.

## References

1. S. Courrech du Pont, P. Gondret, B. Perrin, M. Rabaud, Phys. Rev. Lett., 90, 044301 (2003)

2. E. Lajeunesse, J. Monnier, G. Homsy, Phys. Fluids, 17, 103302 (2005)

3. E. Thompson, H. Huppert, J. Fluid Mech., 575, 177-186 (2007)

4. L. Rondon, O. Pouliquen, P. Aussillous, Phys. Fluids, 23, 073301 (2011)

5. A. Bougouin, L. Lacaze, Phys. Rev. Fluids, 6, 23, (2018).

6. G. Pinzon, M. Cabrera, Phys. Fluids, 31(8), 086603 (2019)

7. V. Topin, Y. Monerie, F. Perales, and F. Radjai, Phys. Rev. Lett. 109, 188001 (2012).

8. Sarlin, W., Morize, C., Sauret, A., & Gondret, P. Phys. Rev. E, 104(6), 064904 (2021).

9. M. Cabrera, N. Estrada, J. Geophys. Res. Solid Earth, 126 (9), e2021JB022589 (2021).

10. O. Polanía, M. Cabrera, M. Renouf, E. Azéma, Phys. Rev. Fluids. (2022).

11. He, K., Shi, H., & Yu, X. Physics of Fluids, 33(10), 103311 (2021).
We observe that the collapse process of polydisperse dry columns is independent of the polydispersity level. For dry monodisperse and polydisperse columns, the collapse initiates at a similar instant and evolves with a similar front velocity. When the collapse comes to an end, the runout slightly varies, being greater for polydisperse columns. On the contrary, for immersed columns, the interaction with a fluid changes the collapse of monodisperse and polydisperse columns. Moreover, the collapse initiation is delayed in polydisperse columns. Then, during the steady propagation stage, the velocity of both monodisperse and polydisperse columns is similar, but polydisperse systems keep this stage for a longer period and their runout considerably increases, compared with monodisperse (see Figure 1).

We show that the differences between monodisperse and polydisperse immersed columns arise from the pore pressure changes $\Delta P$ that columns have beneath them. For the collapse process, we identify the pore pressure evolution in three stages. First, there is a negative change of the pore pressure due to the in-flow of fluid into the column. The fluid comes in into the granular structure when the column dilates to start the motion [4]. Then, during the spreading of the granular flow, the changes of the pore pressure are positive. Finally, pore pressure changes stabilize $\Delta P \approx 0$ (see Figure 2).

Our research provides evidence that the collapse of immersed granular columns strongly depends on the polydispersity level. Our results could be of great use in processes that involve polydisperse granular materials that account for extreme grain size variations, e.g., debris flows, pyroclastic flows, submarine landslides. Our experimental work provides a novel insight of the basal pore pressure variations and the influence that it has on the collapse sequence of immersed flows.

## References

1. S. Courrech du Pont, P. Gondret, B. Perrin, M. Rabaud, Phys. Rev. Lett., 90, 044301 (2003)

2. E. Lajeunesse, J. Monnier, G. Homsy, Phys. Fluids, 17, 103302 (2005)

3. E. Thompson, H. Huppert, J. Fluid Mech., 575, 177-186 (2007)

4. L. Rondon, O. Pouliquen, P. Aussillous, Phys. Fluids, 23, 073301 (2011)

5. A. Bougouin, L. Lacaze, Phys. Rev. Fluids, 6, 23, (2018).

6. G. Pinzon, M. Cabrera, Phys. Fluids, 31(8), 086603 (2019)

7. V. Topin, Y. Monerie, F. Perales, and F. Radjai, Phys. Rev. Lett. 109, 188001 (2012).

8. Sarlin, W., Morize, C., Sauret, A., & Gondret, P. Phys. Rev. E, 104(6), 064904 (2021).

9. M. Cabrera, N. Estrada, J. Geophys. Res. Solid Earth, 126 (9), e2021JB022589 (2021).

10. O. Polanía, M. Cabrera, M. Renouf, E. Azéma, Phys. Rev. Fluids. (2022).

11. He, K., Shi, H., & Yu, X. Physics of Fluids, 33(10), 103311 (2021).
<figure><img src="image_12.png" /><figcaption>Figures</figcaption></figure>

**Figure 1:** Collapse sequence of an immersed granular column with aspect ratio A = H0 / L0 = 2.8 and for two different levels of polydispersity λ = [1.0,20.0] (left and right column). Dots in the base of the experimental setup are pore pressure transducers (PPT). The blue dashed lines indicate the column initial geometry. The horizontal black bar in the right-upper corner of each sequence represents 10 cm.

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

**Figure 2:** Pore pressure changes $\Delta P$ for a column with A = 2.8 and for = [1, 20]. The distance $X_{PPTi}$ corresponds to the position of each PPT. The measure of $\Delta P$ correspond to a range of 0.25 seconds after the arrival of the front to XPPTi. The markers are the median and error bars are the first and third quartile of a range of 0.25 seconds. Note that the first two PPTs are beneath the column before its collapse.
<figure><img src="image_12.png" /><figcaption>Figures</figcaption></figure>

**Figure 1:** Collapse sequence of an immersed granular column with aspect ratio A = H0 / L0 = 2.8 and for two different levels of polydispersity λ = [1.0,20.0] (left and right column). Dots in the base of the experimental setup are pore pressure transducers (PPT). The blue dashed lines indicate the column initial geometry. The horizontal black bar in the right-upper corner of each sequence represents 10 cm.

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

**Figure 2:** Pore pressure changes $\Delta P$ for a column with A = 2.8 and for = [1, 20]. The distance $X_{PPTi}$ corresponds to the position of each PPT. The measure of $\Delta P$ correspond to a range of 0.25 seconds after the arrival of the front to XPPTi. The markers are the median and error bars are the first and third quartile of a range of 0.25 seconds. Note that the first two PPTs are beneath the column before its collapse.
Extending the 3D-H model to snow mechanics: model calibration using oedometer tests

Marie Miot¹, Antoine Wautier², Pit Polfer¹, François Nicot³, Pierre Philippe², Tibor Fulop¹

¹Goodyear Innovation Center Luxembourg, Colmar-Berg, Luxembourg

²INRAE, Aix-Marseille Université, RECOVER

³Université Savoie Mont-Blanc, ISTerre

marie_miot@goodyear.com

**Keywords:** snow mechanics, multi-scale modeling, ice, microstructure

**Abstract**

As a granular material constituted of ice grains, snow has a complex mechanical behaviour which depends on both its microstructure and the behaviour of ice itself. In many science and engineering fields such as avalanche prediction or winter tire design, snow mechanics is a key component to develop predictive and design tools. Several approaches are available to model snow numerically. Discrete element method (DEM), in which the ice grains are modelled explicitly, enable us to capture the detailed influence of the microstructure, but computation time remains too large to address large systems. Continuum approaches are computationally efficient, however, the constitutive laws used for snow modelling do not account well for the various types of snow microstructure. In this work we propose to adopt a multiscale approach by extending a particular micromechanical model (the 3D H-model [Xiong, 2017]) to snow.

The 3D H-model is a multiscale constitutive law that accounts for the microstructure of granular materials (Fig. 1). It relies on the description of a representative elementary volume (REV) as a statistical distribution of mesostructures of a few grains (the H-cells) with different orientations (Fig. 2). In order to provide an analytical description of the mechanical behavior at the mesoscale, the H-cell geometry was chosen as simple as possible in the form of a bi-hexagonal structure of 10 grains. The hypotheses of symmetry in the forces and the geometry of the cell enable us to describe the ten-grain cell with only 6 parameters (4 intergranular distances and 2 opening angles).

The global deformations are obtained at the scale of the REV and the hypothesis of homogeneous strain provides the deformations within the cell. The geometrical description of the cell and the hypothesis of grain force equilibrium make it possible to infer a relationship between local deformations and variations of intergranular distances and opening angles. The contact law at the grain scale gives the intergranular contact forces. Local stresses are deduced from contact forces with the Love-Weber formula, and a statistical homogenization gives the stress at the REV scale.

The snow behaviour is integrated in the 3D-H model by modifying the contact law at microscale to account for the sintering phenomenon. This new ice contact law is adapted from a snow model developed for DEM by [Peters et al., 2021] and [Kabore et al., 20021]. At the beginning of the simulation, two grains in contact are linked by a cylindrical ice bond. The choice of the initial bond radius depends on snow aging (the older the snow, the larger the ice bonds). As long as the bond is active, the interaction between the two grains is modelled using a Maxwell elasto-viscoplastic beam.
Extending the 3D-H model to snow mechanics: model calibration using oedometer tests

Marie Miot¹, Antoine Wautier², Pit Polfer¹, François Nicot³, Pierre Philippe², Tibor Fulop¹

¹Goodyear Innovation Center Luxembourg, Colmar-Berg, Luxembourg

²INRAE, Aix-Marseille Université, RECOVER

³Université Savoie Mont-Blanc, ISTerre

marie_miot@goodyear.com

**Keywords:** snow mechanics, multi-scale modeling, ice, microstructure

**Abstract**

As a granular material constituted of ice grains, snow has a complex mechanical behaviour which depends on both its microstructure and the behaviour of ice itself. In many science and engineering fields such as avalanche prediction or winter tire design, snow mechanics is a key component to develop predictive and design tools. Several approaches are available to model snow numerically. Discrete element method (DEM), in which the ice grains are modelled explicitly, enable us to capture the detailed influence of the microstructure, but computation time remains too large to address large systems. Continuum approaches are computationally efficient, however, the constitutive laws used for snow modelling do not account well for the various types of snow microstructure. In this work we propose to adopt a multiscale approach by extending a particular micromechanical model (the 3D H-model [Xiong, 2017]) to snow.

The 3D H-model is a multiscale constitutive law that accounts for the microstructure of granular materials (Fig. 1). It relies on the description of a representative elementary volume (REV) as a statistical distribution of mesostructures of a few grains (the H-cells) with different orientations (Fig. 2). In order to provide an analytical description of the mechanical behavior at the mesoscale, the H-cell geometry was chosen as simple as possible in the form of a bi-hexagonal structure of 10 grains. The hypotheses of symmetry in the forces and the geometry of the cell enable us to describe the ten-grain cell with only 6 parameters (4 intergranular distances and 2 opening angles).

The global deformations are obtained at the scale of the REV and the hypothesis of homogeneous strain provides the deformations within the cell. The geometrical description of the cell and the hypothesis of grain force equilibrium make it possible to infer a relationship between local deformations and variations of intergranular distances and opening angles. The contact law at the grain scale gives the intergranular contact forces. Local stresses are deduced from contact forces with the Love-Weber formula, and a statistical homogenization gives the stress at the REV scale.

The snow behaviour is integrated in the 3D-H model by modifying the contact law at microscale to account for the sintering phenomenon. This new ice contact law is adapted from a snow model developed for DEM by [Peters et al., 2021] and [Kabore et al., 20021]. At the beginning of the simulation, two grains in contact are linked by a cylindrical ice bond. The choice of the initial bond radius depends on snow aging (the older the snow, the larger the ice bonds). As long as the bond is active, the interaction between the two grains is modelled using a Maxwell elasto-viscoplastic beam.
No other contact law is active at this point. The total bond deformation of the beam $\epsilon_b^{tot}$ can be decomposed into an elastic $\epsilon_b^e$ and a viscous deformation $\epsilon_b^{vp}$:

$$ \epsilon_b^{tot}(t) = \epsilon_b^{vp}(t) + \epsilon_b^e(t) \quad (1) $$

The elastic deformation is directly linked to the stresses in the beam by the oedometer modulus $E_{oedo}$ of ice, as there is no lateral elastic deformation:

$$ \sigma_b(t) = \frac{E(1-\nu)}{(1+\nu)(1-2\nu)} \epsilon_b^e(t) = E_{oedo} \epsilon_b^e(t) \quad (2) $$

The viscous deformation is obtained by an integration over time of the creep rate, which is related to the stress through the Glen's law:

$$ \epsilon_b^{vp}(t) = \epsilon_b^{vp}(t - dt) + \int_{t-dt}^{t} \dot{\epsilon}_b^{vp}(t)dt = \epsilon_b^{vp}(t - dt) + \int_{t-dt}^{t} A(T) \sigma_b^m(t - \delta t)dt \quad (3) $$

This viscous deformation has a lateral component which induces a variation of the bond radius. Assuming bond volume conservation gives a close relationship between the variations of branch vector length and opening angle and the normal and tangential contact forces.

Brittle compression strength, tensile strength and shear strength of ice provide failure threshold values for the normal and tangential bond stress. For practical applications where the time scale is small compared with the sintering characteristic time, bond re-creation does not play an important role and can be omitted.

After bond brittle failure the contact is handled by an elasto-plastic contact between ice grain:

$$ \left\{ \begin{array}{l} \delta N_i = \begin{cases} -k_{n,i} \delta d_i & \text{if } d_i \le 2r_g \\ 0 & \text{else} \end{cases} \\ \\ \delta T_i = \begin{cases} k_{t,i} d_i \delta \alpha_j & \text{if } |T_i + k_{t,i} d_i \delta \alpha_j| < |N_i + \delta N_i| \\ \text{sign}(T_i + k_{t,i} d_i \delta \alpha_j)(N_i + \delta N_i) \tan \varphi_g & \text{else} \end{cases} \end{array} \right. \quad (4) $$

with the normal and tangential contact stiffnesses $k_n$ and $k_t$ depending on the intergranular distance $d_i$:

$$ \begin{cases} k_{n,i} = E_{ice} \frac{\pi}{2} \left( r_g - \frac{d_i^2}{4r_g} \right) \\ k_{t,i} = \nu k_{n,i} \end{cases} \quad (5) $$

Confined compression tests were numerically simulated at REV scale with the developed 3D-H snow model. A parametric study highlights the ability of the model to describe different snow types over a specific range of densities. In particular, the choice of the initial geometry of the 3D-H model cells can be related to the initial density of the numerical sample.

Comparisons with experimental results from [Abele, 1976] show the ability of the model to reproduce the evolution of the axial stress as a function of the axial strain at the REV scale (Fig. 3).

## References

Abele, G., & Gow, A. J. (1975). Compressibility characteristics of compacted snow. *Tech. rep.*, CRREL
No other contact law is active at this point. The total bond deformation of the beam $\epsilon_b^{tot}$ can be decomposed into an elastic $\epsilon_b^e$ and a viscous deformation $\epsilon_b^{vp}$:

$$ \epsilon_b^{tot}(t) = \epsilon_b^{vp}(t) + \epsilon_b^e(t) \quad (1) $$

The elastic deformation is directly linked to the stresses in the beam by the oedometer modulus $E_{oedo}$ of ice, as there is no lateral elastic deformation:

$$ \sigma_b(t) = \frac{E(1-\nu)}{(1+\nu)(1-2\nu)} \epsilon_b^e(t) = E_{oedo} \epsilon_b^e(t) \quad (2) $$

The viscous deformation is obtained by an integration over time of the creep rate, which is related to the stress through the Glen's law:

$$ \epsilon_b^{vp}(t) = \epsilon_b^{vp}(t - dt) + \int_{t-dt}^{t} \dot{\epsilon}_b^{vp}(t) dt = \epsilon_b^{vp}(t - dt) + \int_{t-dt}^{t} A(T) \sigma_b^m(t - \delta t) dt \quad (3) $$

This viscous deformation has a lateral component which induces a variation of the bond radius. Assuming bond volume conservation gives a close relationship between the variations of branch vector length and opening angle and the normal and tangential contact forces.

Brittle compression strength, tensile strength and shear strength of ice provide failure threshold values for the normal and tangential bond stress. For practical applications where the time scale is small compared with the sintering characteristic time, bond re-creation does not play an important role and can be omitted.

After bond brittle failure the contact is handled by an elasto-plastic contact between ice grain:

$$ \left\{ \begin{aligned} \delta N_i &= \begin{cases} -k_{n,i} \delta d_i & \text{if } d_i \le 2r_g \\ 0 & \text{else} \end{cases} \\ \delta T_i &= \begin{cases} k_{t,i} d_i \delta \alpha_j & \text{if } |T_i + k_{t,i} d_i \delta \alpha_j| < |N_i + \delta N_i| \\ \operatorname{sign}(T_i + k_{t,i} d_i \delta \alpha_j)(N_i + \delta N_i) \tan \varphi_g & \text{else} \end{cases} \end{aligned} \right. \quad (4) $$

with the normal and tangential contact stiffnesses $k_n$ and $k_t$ depending on the intergranular distance $d_i$:

$$ \left\{ \begin{aligned} k_{n,i} &= E_{ice} \frac{\pi}{2} \left( r_g - \frac{d_i^2}{4r_g} \right) \\ k_{t,i} &= v k_{n,i} \end{aligned} \right. \quad (5) $$

Confined compression tests were numerically simulated at REV scale with the developed 3D-H snow model. A parametric study highlights the ability of the model to describe different snow types over a specific range of densities. In particular, the choice of the initial geometry of the 3D-H model cells can be related to the initial density of the numerical sample.

Comparisons with experimental results from [Abele, 1976] show the ability of the model to reproduce the evolution of the axial stress as a function of the axial strain at the REV scale (Fig. 3).

## References

* Abele, G., & Gow, A. J. (1975). Compressibility characteristics of compacted snow. *Tech. rep.*, CRREL
Kabore, B. W., Peters, B., Michael, M., Nicot, F. (2021). A discrete element framework for modeling the mechanical behaviour of snow—Part I: Mechanical behaviour and numerical model. *Granular Matter*, 23(2)

Peters, B., Kabore, B. W., Michael, M., Nicot, F. (2021). A discrete element framework for modeling the mechanical behaviour of snow Part II: model validation. *Granular Matter*, 23(2)

Xiong, H., Nicot, F., Yin, Z. Y. (2017). A three-dimensional micromechanically based model. *International Journal for Numerical and Analytical Methods in Geomechanics*, 41(17)

## Figures

<figure><img src="image_14.png" /><figcaption>**Figure 1:** Scheme of the H-model constitutive law extended to snow material</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>**Figure 2:** Statistical distribution of hexagonal cells in 2D (left) and 3D H-model cell (right)</figcaption></figure>
Kabore, B. W., Peters, B., Michael, M., Nicot, F. (2021). A discrete element framework for modeling the mechanical behaviour of snow—Part I: Mechanical behaviour and numerical model. *Granular Matter*, 23(2)

Peters, B., Kabore, B. W., Michael, M., Nicot, F. (2021). A discrete element framework for modeling the mechanical behaviour of snow Part II: model validation. *Granular Matter*, 23(2)

Xiong, H., Nicot, F., Yin, Z. Y. (2017). A three-dimensional micromechanically based model. *International Journal for Numerical and Analytical Methods in Geomechanics*, 41(17)

## Figures

<figure><img src="image_14.png" /><figcaption>**Figure 1:** Scheme of the H-model constitutive law extended to snow material</figcaption></figure>

<figure><img src="image_15.png" /><figcaption>**Figure 2:** Statistical distribution of hexagonal cells in 2D (left) and 3D H-model cell (right)</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>**Figure 3:** Reproduction of experimental oedometer tests from [Abele, 1976] with the 3D-H model at REV sale, for different initial density, temperature, and snow aging.</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>**Figure 3:** Reproduction of experimental oedometer tests from [Abele, 1976] with the 3D-H model at REV sale, for different initial density, temperature, and snow aging.</figcaption></figure>
# Multiscale modeling of oedometer test on crushable tubular particles

Na DENG<fnref n="1" />, Trung-Kien NGUYEN<fnref n="2" />, Vincent RICHEFEU<fnref n="1" />, and Gaël COMBE<fnref n="1" />

**Keywords:** Grain breakage; Bonded spheres; Multiscale modeling; FEM×DEM

## Abstract

This research is dedicated to the numerical modeling of crushable tubular particles under oedometer conditions. First, a numerical tube capable being broken is proposed using a 3D periodic boundary DEM in-house software. The tube is made of bonded spheres Fig.1. It is designed according to the diameter, length and thickness of a real particle, made of baked clay, hereafter called shell. The tube is composed of 8 layers of 22 spheres for its cross section. Two contact models are used: an elasto-plastic frictional contact law and a damageable bonded law corresponding to non-bonded interactions and bonded interactions, respectively. To determine the state of a damaged bonded contact, a criterion based on relative normal ($d_n$) and tangential displacement ($d_t$) and relative rotation ($d_r$) is adopted, following the formulation by (Delenne et al., 2004).

$$G(d_n^b, d_t^b, d_r^b) = \frac{d_n^b}{-d_n^{\text{damage}}} + \left(\frac{|d_t^b|}{d_t^{\text{damage}}}\right)^\alpha + \left(\frac{|d_r^b|}{d_r^{\text{damage}}}\right)^\alpha - 1 = 0 \quad (1)$$

where $d_n^{\text{damage}}$, $d_t^{\text{damage}}$, $d_r^{\text{damage}}$ and $\alpha$ are constants. Compression is counted positive. Parameters of the contact models are calibrated using experimental results of an uniaxial radial compression and an oedometer test on a single shell (Stasiak et al., 2022). The numerical and experimental responses are shown in Fig.3-4. The localisation of primary (vertical) and secondary (horizontal) cracks during the uniaxial compression test are presented in the insets in Fig.3 (a).

Based on the calibrated parameters, we conducted an oedometer loading on a volume element composed of 27 numerical tubes using a coupled FEM×DEM method (Desrues et al., 2019, Nguyen et al., 2022). In the FEM×DEM method, Lagamine FEM code developed in Liège University is used (Charlier, 1989). Fig.5 shows the boundary conditions and one finite element with 8 nodes and 4 integration points. Preliminary results are shown in Fig. 6, compared with the result from a pure DEM simulation and experimental data. Two views of contact forces from one integration point at the initial state and during oedometer test are presented in Fig.7. Qualitatively speaking, the trends of the curves from predictions and the experiment are similar. Quantitatively, the level of predicted stress is lower than the experimental result. This can arise for multiple reasons: (1) the void ratio of the numerical tube is larger than that of the real tube; (2) the orientation of tubes would differ to the reality.

<footnote n="1">¹Université Grenoble Alpes, 3SR, Grenoble, France</footnote>

<footnote n="2">²Hanoi Univ Civil Engn, Fac Bldg & Ind Construction, Hanoi, Vietnam,
na.deng@3sr-grenoble.fr; kiennt3@huce.edu.vn; vincent.richefeu@3sr-grenoble.fr; gael.combe@3sr-grenoble.fr</footnote>
# Multiscale modeling of oedometer test on crushable tubular particles

Na DENG<fnref n="1" />, Trung-Kien NGUYEN<fnref n="2" />, Vincent RICHEFEU<fnref n="1" />, and Gaël COMBE<fnref n="1" />

**Keywords:** Grain breakage; Bonded spheres; Multiscale modeling; FEM×DEM

## Abstract

This research is dedicated to the numerical modeling of crushable tubular particles under oedometer conditions. First, a numerical tube capable being broken is proposed using a 3D periodic boundary DEM in-house software. The tube is made of bonded spheres Fig.1. It is designed according to the diameter, length and thickness of a real particle, made of baked clay, hereafter called shell. The tube is composed of 8 layers of 22 spheres for its cross section. Two contact models are used: an elasto-plastic frictional contact law and a damageable bonded law corresponding to non-bonded interactions and bonded interactions, respectively. To determine the state of a damaged bonded contact, a criterion based on relative normal ($d_n$) and tangential displacement ($d_t$) and relative rotation ($d_r$) is adopted, following the formulation by (Delenne et al., 2004).

$$G(d_n^b, d_t^b, d_r^b) = \frac{d_n^b}{-d_n^{\text{damage}}} + \left(\frac{|d_t^b|}{d_t^{\text{damage}}}\right)^{\alpha} + \left(\frac{|d_r^b|}{d_r^{\text{damage}}}\right)^{\alpha} - 1 = 0 \quad (1)$$

where $d_n^{\text{damage}}$, $d_t^{\text{damage}}$, $d_r^{\text{damage}}$ and $\alpha$ are constants. Compression is counted positive. Parameters of the contact models are calibrated using experimental results of an uniaxial radial compression and an oedometer test on a single shell (Stasiak et al., 2022). The numerical and experimental responses are shown in Fig.3-4. The localisation of primary (vertical) and secondary (horizontal) cracks during the uniaxial compression test are presented in the insets in Fig.3 (a).

Based on the calibrated parameters, we conducted an oedometer loading on a volume element composed of 27 numerical tubes using a coupled FEM×DEM method (Desrues et al., 2019, Nguyen et al., 2022). In the FEM×DEM method, Lagamine FEM code developed in Liège University is used (Charlier, 1989). Fig.5 shows the boundary conditions and one finite element with 8 nodes and 4 integration points. Preliminary results are shown in Fig. 6, compared with the result from a pure DEM simulation and experimental data. Two views of contact forces from one integration point at the initial state and during oedometer test are presented in Fig.7. Qualitatively speaking, the trends of the curves from predictions and the experiment are similar. Quantitatively, the level of predicted stress is lower than the experimental result. This can arise for multiple reasons: (1) the void ratio of the numerical tube is larger than that of the real tube; (2) the orientation of tubes would differ to the reality.

<footnote n="1">¹Université Grenoble Alpes, 3SR, Grenoble, France</footnote>

<footnote n="2">²Hanoi Univ Civil Engn, Fac Bldg & Ind Construction, Hanoi, Vietnam, na.deng@3sr-grenoble.fr; kiennt3@huce.edu.vn; vincent.richefeu@3sr-grenoble.fr; gael.combe@3sr-grenoble.fr</footnote>
<figure><img src="image_17.png" /><figcaption>**Figure 1:** A crushable tubular particle composed of bonded spheres.</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>**Figure 2:** A 3D presentation of the yield surface, $\alpha = 2$.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>**Figure 3:** Deformation and force along uniaxial radial compression: (a) prediction from this research; (b) results of experiments (green) and of a simulation (red) from (Stasiak et al., 2022).</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>Figure 1: A crushable tubular particle composed of bonded spheres.</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>Figure 2: A 3D presentation of the yield surface, $\alpha = 2$.</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>Figure 3: Deformation and force along uniaxial radial compression: (a) prediction from this research; (b) results of experiments (green) and of a simulation (red) from (Stasiak et al., 2022).</figcaption></figure>
<figure><img src="image_20.png" /><figcaption>**Figure 4:** Deformation and force along oedometer loading: (a) prediction from this research (b) results of experiments (green) and of a simulation (red) from (Stasiak, 2019).</figcaption></figure>

<figure><img src="image_21.png" /><figcaption>**Figure 5:** One element with 8 nodes and 4 integration points. At each integration point, the constitutive law is provided by a DEM simulation of tubes.</figcaption></figure>

<figure><img src="image_22.png" /><figcaption>**Figure 6:** Stress-strain curves from simulations using pure DEM and FEM×DEM and from an experiment. $\epsilon_a = \ln \frac{H_0}{H_0 - \Delta H}$, where $\Delta H$ is the deformation along the loading direction; $H_0$ is the initial height. The inset shows the a top view of broken tubes in the oedometer test when $\epsilon_a = 6.8MPa$ (Stasiak, 2019).</figcaption></figure>
<figure><img src="image_20.png" /><figcaption>**Figure 4:** Deformation and force along oedometer loading: (a) prediction from this research (b) results of experiments (green) and of a simulation (red) from (Stasiak, 2019).</figcaption></figure>

<figure><img src="image_21.png" /><figcaption>**Figure 5:** One element with 8 nodes and 4 integration points. At each integration point, the constitutive law is provided by a DEM simulation of tubes.</figcaption></figure>

<figure><img src="image_22.png" /><figcaption>**Figure 6:** Stress-strain curves from simulations using pure DEM and FEM×DEM and from an experiment. $\epsilon_a = \ln \frac{H_0}{H_0 - \Delta H}$, where $\Delta H$ is the deformation along the loading direction; $H_0$ is the initial height. The inset shows the a top view of broken tubes in the oedometer test when $\epsilon_a = 6.8MPa$ (Stasiak, 2019).</figcaption></figure>
<figure><img src="image_23.png" /><figcaption>**Figure 7:** Two front views of contact forces from one integration point at the initial state and during oedometer loading (bonded contact in blue and non-bonded contact in black).</figcaption></figure>

REFERENCES

R. Charlier. Approche unifée de quelques problemes non lineaires de mecanique des milieux continus par la methode des elements finis. 1989.

J.-Y. Delenne, M. S. El Youssoufi, F. Cherblanc, and J.-C. Bénet. Mechanical behaviour and failure of cohesive granular materials. *International Journal for Numerical and Analytical Methods in Geomechanics*, 28(15):1577-1594, 2004.

J. Desrues, A. Argilaga, D. Caillerie, G. Combe, T. K. Nguyen, V. Richefeu, and S. Dal Pont. From discrete to continuum modelling of boundary value problems in geomechanics: An integrated fem-dem approach. *International Journal for Numerical and Analytical Methods in Geomechanics*, 43 (5):919-955, 2019.

T.-K. Nguyen, J. Desrues, T.-T. Vo, and G. Combe. Fem× dem multi-scale model for cemented granular materials: Inter-and intra-granular cracking induced strain localisation. *International Journal for Numerical and Analytical Methods in Geomechanics*, 46(5):1001-1025, 2022.

M. Stasiak. Uniaxial compression of a highly crushable granular material-a 3D DEM study. PhD thesis, Université Grenoble Alpes, 2019.

M. Stasiak, G. Combe, V. Richefeu, G. Armand, and J. Zghondi. High compression of granular assemblies of brittle hollow tubular particles: investigation through a 3d discrete element model. *Computational Particle Mechanics*, 9:1-18, 2022.
<figure><img src="image_23.png" /><figcaption>**Figure 7:** Two front views of contact forces from one integration point at the initial state and during oedometer loading (bonded contact in blue and non-bonded contact in black).</figcaption></figure>

REFERENCES

R. Charlier. Approche unifée de quelques problemes non lineaires de mecanique des milieux continus par la methode des elements finis. 1989.

J.-Y. Delenne, M. S. El Youssoufi, F. Cherblanc, and J.-C. Bénet. Mechanical behaviour and failure of cohesive granular materials. *International Journal for Numerical and Analytical Methods in Geomechanics*, 28(15):1577-1594, 2004.

J. Desrues, A. Argilaga, D. Caillerie, G. Combe, T. K. Nguyen, V. Richefeu, and S. Dal Pont. From discrete to continuum modelling of boundary value problems in geomechanics: An integrated fem-dem approach. *International Journal for Numerical and Analytical Methods in Geomechanics*, 43 (5):919-955, 2019.

T.-K. Nguyen, J. Desrues, T.-T. Vo, and G. Combe. Fem× dem multi-scale model for cemented granular materials: Inter-and intra-granular cracking induced strain localisation. *International Journal for Numerical and Analytical Methods in Geomechanics*, 46(5):1001-1025, 2022.

M. Stasiak. Uniaxial compression of a highly crushable granular material-a 3D DEM study. PhD thesis, Université Grenoble Alpes, 2019.

M. Stasiak, G. Combe, V. Richefeu, G. Armand, and J. Zghondi. High compression of granular assemblies of brittle hollow tubular particles: investigation through a 3d discrete element model. *Computational Particle Mechanics*, 9:1-18, 2022.
# Instability in granular slopes subjected to internal erosion

Pierre-Yves HICHER¹

¹UMR CNRS GeM, Ecole Centrale de Nantes

## Abstract

Fluidised landslides are one of the most dangerous types of mass movements as they can run over long distances at high velocity. A flow-like landslide can occur in both artificially designed and natural slopes, resulting in extensive property damage and significant loss of life. Although various studies have examined the initiation mechanism of this type of landslide, the understanding of the phenomenon of instability which can occur in loose granular slopes deserves further attention. Specific influential factors, such as the changes in the soil microstructure due to water infiltration and seepage forces, still need to be further examined.

In this study, the initiation of a fluidised landslide is investigated through flume tests. The tested material was collected from a coseismic landslide deposit in the 2008 Wenchuan earthquake area. These loose, granular deposits can fail due to intense rainfall and exhibit flow-like movement, potentially evolving into destructive debris flows. Two flume apparatuses have been used (see Fig. 1), which were equipped with multiple sensors to evaluate the pore pressures and the internal displacements. A high-sensitivity seismic accelerometer was installed at the bottom of the flume to record the vibrations produced by the movement of the soil particles during the tests.

The soil mass was progressively wetted either by a uniform inflow from the ground surface at the top of the slope, to simulate the infiltration of run-off water arriving from upslope, or by a gradual rise of the water level at the upper boundary of the slope, to simulate a gradual rise of the groundwater table. Internal erosion of the finer soil fraction driven by the seepage forces is thought to have played a significant role in the slope collapse. It appears that the decrease of the fines content due to internal erosion plays an important role in initiating debris flows by increasing the size of the instability domain of the material constituting the slope as illustrated in Fig. 2 (Hicher, 2013, Yin et al., 2014).

Possible changes in the particle arrangement were recorded by the accelerometers in the form of an increasing high-frequency vibration signal before any significant rise of pore pressure and any significant displacement of the slope could be observed (Hu et al., 2018). Although further parametric investigation is still necessary, and tests at real scale should be devised, the observed behaviour suggests that seismic monitoring of loose, granular deposits could be a promising method for detecting an incipient slope collapse.

## REFERENCES

* P.-Y. Hicher. Modelling the impact of particle removal on granular material behaviour. *Géotechnique*, 63(2):118–128, 2013.

* W. Hu, P.-Y. Hicher, G. Scaringi, Q. Xu, T. Van Asch, and G. Wang. Seismic precursor to instability induced by internal erosion in loose granular slopes. *Géotechnique*, 68(11):989–1001, 2018.

* Z.-Y. Yin, J. Zhao, and P.-Y. Hicher. A micromechanics-based model for sand-silt mixtures. *International journal of solids and structures*, 51(6):1350–1363, 2014.
# Instability in granular slopes subjected to internal erosion

Pierre-Yves HICHER¹

¹UMR CNRS GeM, Ecole Centrale de Nantes

## Abstract

Fluidised landslides are one of the most dangerous types of mass movements as they can run over long distances at high velocity. A flow-like landslide can occur in both artificially designed and natural slopes, resulting in extensive property damage and significant loss of life. Although various studies have examined the initiation mechanism of this type of landslide, the understanding of the phenomenon of instability which can occur in loose granular slopes deserves further attention. Specific influential factors, such as the changes in the soil microstructure due to water infiltration and seepage forces, still need to be further examined.

In this study, the initiation of a fluidised landslide is investigated through flume tests. The tested material was collected from a coseismic landslide deposit in the 2008 Wenchuan earthquake area. These loose, granular deposits can fail due to intense rainfall and exhibit flow-like movement, potentially evolving into destructive debris flows. Two flume apparatuses have been used (see Fig. 1), which were equipped with multiple sensors to evaluate the pore pressures and the internal displacements. A high-sensitivity seismic accelerometer was installed at the bottom of the flume to record the vibrations produced by the movement of the soil particles during the tests.

The soil mass was progressively wetted either by a uniform inflow from the ground surface at the top of the slope, to simulate the infiltration of run-off water arriving from upslope, or by a gradual rise of the water level at the upper boundary of the slope, to simulate a gradual rise of the groundwater table. Internal erosion of the finer soil fraction driven by the seepage forces is thought to have played a significant role in the slope collapse. It appears that the decrease of the fines content due to internal erosion plays an important role in initiating debris flows by increasing the size of the instability domain of the material constituting the slope as illustrated in Fig. 2 (Hicher, 2013, Yin et al., 2014).

Possible changes in the particle arrangement were recorded by the accelerometers in the form of an increasing high-frequency vibration signal before any significant rise of pore pressure and any significant displacement of the slope could be observed (Hu et al., 2018). Although further parametric investigation is still necessary, and tests at real scale should be devised, the observed behaviour suggests that seismic monitoring of loose, granular deposits could be a promising method for detecting an incipient slope collapse.

## REFERENCES

* P.-Y. Hicher. Modelling the impact of particle removal on granular material behaviour. *Géotechnique*, 63(2):118–128, 2013.

* W. Hu, P.-Y. Hicher, G. Scaringi, Q. Xu, T. Van Asch, and G. Wang. Seismic precursor to instability induced by internal erosion in loose granular slopes. *Géotechnique*, 68(11):989–1001, 2018.

* Z.-Y. Yin, J. Zhao, and P.-Y. Hicher. A micromechanics-based model for sand-silt mixtures. *International journal of solids and structures*, 51(6):1350–1363, 2014.
<figure><img src="image_24.png" /><figcaption>Figure 1: Flume test at SKLGP</figcaption></figure>

Influence of fine particle content on soil instability

<figure><img src="image_25.png" /><figcaption>Figure 2: Undrained triaxial tests</figcaption></figure>
<figure><img src="image_24.png" /><figcaption>Figure 1: Flume test at SKLGP</figcaption></figure>

<figure><img src="image_25.png" /><figcaption>Figure 2: Undrained triaxial tests</figcaption></figure>
Field overflow tests on a Rhône levee

S. Nicaise¹, C. Picault², N. Chaouch¹, F. Byron¹, Y. Grémeaux¹, A. Doghmane¹, F. Crampette², F. Golay³, S. Bonelli¹

1. INRAE, Aix Marseille Univ, RECOVER, Aix-en-Provence

2. Compagnie Nationale du Rhône (CNR), Lyon

3. Institut de Mathématique de Toulon (IMATH), Toulon

sylvie.nicaise@inrae.fr, c.picault@cnr.tn.fr, frederic.golay@univ-tln.fr

**Keywords:** earth dikes, levees, overflow, field tests

**Abstract**

The risk of dam or levee failure due to overflowing erosion is still poorly understood: little literature describes the processes, in particular the processes specific to levees. In order to fill some gaps, and to improve their dam safety analysis, Electricité de France (EDF) and the Compagnie Nationale du Rhône (CNR) have launched the OVERCOME project. This project includes multi-scale experiments that will allow to describe and understand the physical processes during overflow.

As part of this work, INRAE performed two in-situ tests on a CNR levee at Montfaucon in May 2022. The levee is 6.2m high, with a slope of 4horizontal/1vertical on the upper part consisting of compacted silt containing some gravels; at the basis of this upper part, a layer of sandy gravel (1m thick) is located. A gravel fill constitutes the lower part, with a slope 3horizontal/1vertical (Figure1). A berm separates the two parts, 6m wide.

<figure><img src="image_26.png" /><figcaption>**Figure 1:** cross-section of the dike and sieve analysis of the different materials</figcaption></figure>

Based on a standardized protocol (D6460, ASTM, 2012), the principle of these tests is to carry out a flow in a channel; several phases are observed, with increase of the flow. Measurements are made during the flow to quantify the hydraulic parameters, and during a pause to measure erosion. Flow rates obtained were from 12 to 140l.s⁻¹ per linear meter, corresponding to a crest water height of 4 to 18cm.

The first test was carried out on intact fill, covered with natural grass. Twelve steps of overflowing, 30min each, were observed on this channel of 1m wide and 25m long.

After 6h of flow over the vegetation, with a maximum velocity of 4m.s⁻¹ corresponding to a hydraulic stress of 20Pa, almost no erosion was observed. However, on the same channel, a scour hole with an approximate volume of 4.2m³ was created on the gravel fill (Figure2).
Field overflow tests on a Rhône levee

S. Nicaise¹, C. Picault², N. Chaouch¹, F. Byron¹, Y. Grémeaux¹, A. Doghmane¹, F. Crampette², F. Golay³, S. Bonelli¹

1. INRAE, Aix Marseille Univ, RECOVER, Aix-en-Provence

2. Compagnie Nationale du Rhône (CNR), Lyon

3. Institut de Mathématique de Toulon (IMATH), Toulon

sylvie.nicaise@inrae.fr, c.picault@cnr.tn.fr, frederic.golay@univ-tln.fr

**Keywords:** earth dikes, levees, overflow, field tests

**Abstract**

The risk of dam or levee failure due to overflowing erosion is still poorly understood: little literature describes the processes, in particular the processes specific to levees. In order to fill some gaps, and to improve their dam safety analysis, Electricité de France (EDF) and the Compagnie Nationale du Rhône (CNR) have launched the OVERCOME project. This project includes multi-scale experiments that will allow to describe and understand the physical processes during overflow.

As part of this work, INRAE performed two in-situ tests on a CNR levee at Montfaucon in May 2022. The levee is 6.2m high, with a slope of 4horizontal/1vertical on the upper part consisting of compacted silt containing some gravels; at the basis of this upper part, a layer of sandy gravel (1m thick) is located. A gravel fill constitutes the lower part, with a slope 3horizontal/1vertical (Figure1). A berm separates the two parts, 6m wide.

<figure><img src="image_26.png" /><figcaption>**Figure 1:** cross-section of the dike and sieve analysis of the different materials</figcaption></figure>

Based on a standardized protocol (D6460, ASTM, 2012), the principle of these tests is to carry out a flow in a channel; several phases are observed, with increase of the flow. Measurements are made during the flow to quantify the hydraulic parameters, and during a pause to measure erosion. Flow rates obtained were from 12 to 140l.s⁻¹ per linear meter, corresponding to a crest water height of 4 to 18cm.

The first test was carried out on intact fill, covered with natural grass. Twelve steps of overflowing, 30min each, were observed on this channel of 1m wide and 25m long.

After 6h of flow over the vegetation, with a maximum velocity of 4m.s⁻¹ corresponding to a hydraulic stress of 20Pa, almost no erosion was observed. However, on the same channel, a scour hole with an approximate volume of 4.2m³ was created on the gravel fill (Figure2).
<figure><img src="image_27.png" /><figcaption>**Figure 2:** channel 1: flow phase (left), surface after the test (middle) and erosion on the gravel fill (right)</figcaption></figure>

The second test concerned soil stripped of vegetation and topsoil, with a channel 60cm wide and 20m long. It lasted less time than the first test: 9 steps of 30min each. We measured a maximum flow velocity of 2.6m.s⁻¹, corresponding to a stress of 110Pa. Erosion occurred from the beginning of the flow; we can consider that during the first phase, erosion corresponds to the leaching of the least resistant materials on the surface of the slope. During the following phases, erosion remained important on the silty slope (Figure3). It led to the creation of steps on the slope. In particular, erosion of the sandy gravel at the base of the upper section formed a step with the silt. This step regressed during the test, the starting mechanism of the silt then seemed to be more mechanical instability than erosion. After 9 steps, the eroded volume was 9m³, we were forced to stop the test for safety reasons, and no measurements were possible.

<figure><img src="image_28.png" /><figcaption>**Figure 3:** channel 2: volume erosion rate during the 8 first phases, for 3 area of interest</figcaption></figure>

One result of these experiments is the remarkable strength of the vegetated slope. A second result is the negative influence of heterogeneities on the silty slope: dramatic erosion was found to come from erosion of the gravel layer, followed by instability of the silt surrounding this layer. These experiments underscore the importance of in-situ testing to better understand overflowing erosion.

## References

* ASTM (2012). Standard test method for determination of Rolled Erosion Control products (RECP) performance in protecting earthen channel from stormwater-induced erosion. Technical Report D 6460.

* Bonelli S. (edt), Erosion of Geomaterials, Wiley/ISTE, 371 p., 2012.
<figure><img src="image_27.png" /><figcaption>**Figure 2:** channel 1: flow phase (left), surface after the test (middle) and erosion on the gravel fill (right)</figcaption></figure>

The second test concerned soil stripped of vegetation and topsoil, with a channel 60cm wide and 20m long. It lasted less time than the first test: 9 steps of 30min each. We measured a maximum flow velocity of 2.6m.s⁻¹, corresponding to a stress of 110Pa. Erosion occurred from the beginning of the flow; we can consider that during the first phase, erosion corresponds to the leaching of the least resistant materials on the surface of the slope. During the following phases, erosion remained important on the silty slope (Figure3). It led to the creation of steps on the slope. In particular, erosion of the sandy gravel at the base of the upper section formed a step with the silt. This step regressed during the test, the starting mechanism of the silt then seemed to be more mechanical instability than erosion. After 9 steps, the eroded volume was 9m³, we were forced to stop the test for safety reasons, and no measurements were possible.

<figure><img src="image_28.png" /><figcaption>**Figure 3:** channel 2: volume erosion rate during the 8 first phases, for 3 area of interest</figcaption></figure>

One result of these experiments is the remarkable strength of the vegetated slope. A second result is the negative influence of heterogeneities on the silty slope: dramatic erosion was found to come from erosion of the gravel layer, followed by instability of the silt surrounding this layer. These experiments underscore the importance of in-situ testing to better understand overflowing erosion.

## References

ASTM (2012). Standard test method for determination of Rolled Erosion Control products (RECP) performance in protecting earthen channel from stormwater-induced erosion. Technical Report D 6460.

Bonelli S. (edt), Erosion of Geomaterials, Wiley/ISTE, 371 p., 2012.
Bonelli S. (edt). Erosion in geomechanics applied to dams and levees, Wiley, 388 p., 2013.

Bonelli S., Nicaise S., Charrier G., Chaouch N., Byron F., Gremeaux Y. (2018), Quantifying the erosion resistance of dikes with the overflowing simulator, 3rd International Conference on Protection against Overtopping, 6-8 June 2018, UK.

Cantré S., Olschewski J., Saathoff F. (2017). Full-Scale Flume Experiments to Analyze the Surface Erosion Re-sistance of Dike Embankments Made of Dredged Mate-rials, J. Waterway Port. Coast. Ocean. Eng. ASCE.

Nerincx N., Bonelli S., Mercier F., Cornacchioli F., Fry J.-J., Herrier G., Richard J.-M., Puiatti D., Tachker P. (2018), DigueElite overflow resistant lime treated soils for dikes and earthdams, 26th International Congress on Large Dam, Vienna, Austria, 4-6 July.
Bonelli S. (edt). Erosion in geomechanics applied to dams and levees, Wiley, 388 p., 2013.

Bonelli S., Nicaise S., Charrier G., Chaouch N., Byron F., Gremeaux Y. (2018), Quantifying the erosion resistance of dikes with the overflowing simulator, 3rd International Conference on Protection against Overtopping, 6-8 June 2018, UK.

Cantré S., Olschewski J., Saathoff F. (2017). Full-Scale Flume Experiments to Analyze the Surface Erosion Re-sistance of Dike Embankments Made of Dredged Mate-rials, J. Waterway Port. Coast. Ocean. Eng. ASCE.

Nerincx N., Bonelli S., Mercier F., Cornacchioli F., Fry J.-J., Herrier G., Richard J.-M., Puiatti D., Tachker P. (2018), DigueElite overflow resistant lime treated soils for dikes and earthdams, 26th International Congress on Large Dam, Vienna, Austria, 4-6 July.
Consequences of scalping and scalping/replacement
procedures on strength properties of coarse-grained
gap-graded soils

March 31, 2023

Ngoc Son NGUYEN¹,* Nadine ALI HASSAN¹ Didier MAROT¹ Fateh BENDAHMANE¹

¹GeM Institute, University of Nantes, 58 rue Michel Ange, BP 420, 44606, Saint-Nazaire Cedex, France

*Email: ngocson.nguyen@univ-nantes.fr

**Key words:** Coarse-grained soils, scalping, scalping/replacement, mechanical behavior, triaxial device, discret element method

Coarse-grained soils are characterized by widely graded particle size distributions. For this particular class of soils, large grains (stone, gravel) coexist with fine particles (sand, silt, clay). Rock fill, alluvial soils, till, sandy and/or silty gravels are examples of coarse-grained soils. They are widely used as construction materials for embankment dams, dikes, subgrade of roads, banks, etc. The measurement of the mechanical properties of these granular soils requires very large triaxial apparatuses or shear boxes whose whose diameter can go up to 1 m. Nevertheless, only a few geotechnical laboratories in the world are equipped with such large devices. It is very tedious to prepare large specimens and perform tests. In addition, it is not easy to guarantee the homogeneity of large samples and the repeatability of their strength properties. Geotechnical laboratories are often equipped with experimental apparatuses such as triaxial cells and shear boxes whose diameter rarely exceeds 300 mm. In this context, with the aim to perform tests on coarse-grained soils with apparatuses available in the laboratory, the scalping and scalping/replacement procedures are often used. The scalping consists in removing all the oversized grains from a given original soil and then reconstituting specimens with the remaining fraction, while the scalping/replacement replaces the oversized grains by an equal mass of smaller grains. While these procedures are usually used to evaluate strength properties of coarse-grained soils, no clear experimental standard has been well established so far. Consequences of these specimen reconstitution procedures on obtained strength properties have not been well understood yet. The scalping procedure has been extensively studied by several authors who showed divergent conclusions (Donaghe and Torrey III, 1985, El Dine, 2007, Reiffsteck et al., 2007). A few studies on the scalping/replacement procedure have been carried out; and, like the scalping one, the conclusions about this procedure are contradictory. It should be noted that the reconstituted soils are different from the original soil so their mechanical properties cannot be expected to be the same as those of the original soil. Nevertheless, when the scalping or scalping/replacement procedure is necessary to evaluate the mechanical properties of a given coarse-grained soil, it is important to minimize this difference and manage consequences of the employed procedure on the obtained strength properties. Therefore, several questions remain still open. (i) Which procedure should be used for a given soil? (ii) Which compactness parameter needs to be controlled for the scalped or scalped/replaced specimen? (iii) If the scalping/replacement procedure is adopted, what is the material that can be used to replace oversized particles, particularly regarding to the
Consequences of scalping and scalping/replacement
procedures on strength properties of coarse-grained
gap-graded soils

March 31, 2023

Ngoc Son NGUYEN¹,* Nadine ALI HASSAN¹ Didier MAROT¹ Fateh BENDAHMANE¹

¹GeM Institute, University of Nantes, 58 rue Michel Ange, BP 420, 44606, Saint-Nazaire Cedex, France

*Email: ngocson.nguyen@univ-nantes.fr

**Key words:** Coarse-grained soils, scalping, scalping/replacement, mechanical behavior, triaxial device, discret element method

Coarse-grained soils are characterized by widely graded particle size distributions. For this particular class of soils, large grains (stone, gravel) coexist with fine particles (sand, silt, clay). Rock fill, alluvial soils, till, sandy and/or silty gravels are examples of coarse-grained soils. They are widely used as construction materials for embankment dams, dikes, subgrade of roads, banks, etc. The measurement of the mechanical properties of these granular soils requires very large triaxial apparatuses or shear boxes whose whose diameter can go up to 1 m. Nevertheless, only a few geotechnical laboratories in the world are equipped with such large devices. It is very tedious to prepare large specimens and perform tests. In addition, it is not easy to guarantee the homogeneity of large samples and the repeatability of their strength properties. Geotechnical laboratories are often equipped with experimental apparatuses such as triaxial cells and shear boxes whose diameter rarely exceeds 300 mm. In this context, with the aim to perform tests on coarse-grained soils with apparatuses available in the laboratory, the scalping and scalping/replacement procedures are often used. The scalping consists in removing all the oversized grains from a given original soil and then reconstituting specimens with the remaining fraction, while the scalping/replacement replaces the oversized grains by an equal mass of smaller grains. While these procedures are usually used to evaluate strength properties of coarse-grained soils, no clear experimental standard has been well established so far. Consequences of these specimen reconstitution procedures on obtained strength properties have not been well understood yet. The scalping procedure has been extensively studied by several authors who showed divergent conclusions (Donaghe and Torrey III, 1985, El Dine, 2007, Reiffsteck et al., 2007). A few studies on the scalping/replacement procedure have been carried out; and, like the scalping one, the conclusions about this procedure are contradictory. It should be noted that the reconstituted soils are different from the original soil so their mechanical properties cannot be expected to be the same as those of the original soil. Nevertheless, when the scalping or scalping/replacement procedure is necessary to evaluate the mechanical properties of a given coarse-grained soil, it is important to minimize this difference and manage consequences of the employed procedure on the obtained strength properties. Therefore, several questions remain still open. (i) Which procedure should be used for a given soil? (ii) Which compactness parameter needs to be controlled for the scalped or scalped/replaced specimen? (iii) If the scalping/replacement procedure is adopted, what is the material that can be used to replace oversized particles, particularly regarding to the
particle size gradation and particle properties? The main objective of this research is to contribute to answering these complex questions for gap-graded cohesionless soils which have not been well studied in the literature. For this study, numerical simulations have been first performed by using the discrete element method (DEM) to better understand factors influencing the behavior of scalped and scalped/replaced soils. Then, experimental tests have been performed on gap-graded soils by using triaxial apparatuses with small, medium and large diameters (50, 100, 200 and 300 mm).

Figure 1 present the grain size distributions (GSD) of a gap-graded soil considered in our study and its corresponding scalped and scalped/replaced ones. It can be seen that the scalping procedure leads to a strong change in the GSD of the tested material; it leads actually to an increase in the fine content of gap-graded soils. For the scalping/replacement procedure, the oversized fraction can be replaced by either the whole remaining coarse fraction (scalped/replaced soil 1) or the class of materials which is closest to the oversized fraction (scalped/replaced soil 2). This procedure allows us to preserve the soil's fine content; however, it reduces the gradation of the coarse fraction.

<figure><img src="image_29.png" /><figcaption>**Figure 1:** GSDs of an original soil with $f_f$ = 30% and its corresponding scalped, scalped/replaced 1 and scalped/replaced 2 soils.</figcaption></figure>

To better understand consequences of the scalping and scalping/replacement procedures, it is vital to understand the effect of the fine content and the gradation of the coarse fraction on the behavior of gap-graded soils. For this purpose, the DEM is used to simulate two series of triaxial tests on bidisperse soils with different fine contents and different gradations of the coarse fraction. The fine particles have diameters $d$ from $d_{min}$ = 1 mm to $d_{max}$ = 2 mm; and the coarse particles have diameters $d$ from $D_{min}$ = 6 mm to $D_{max}$. For the first series, the fine content $f_f$ is varied from 20% to 75% while $D_{max}$ is kept to 12 mm. For the second series, the gradation of the coarse fraction characterized by the ratio $D_{max}/D_{min}$ is varied from 1.5 to 3. Figure 2.a shows the behavior of bidisperse samples with different fine contents $f_f$ when the global void ratio $e$ is controlled to be around 0.32. It can be seen that for $f_f$ < 30% an increase in fine content leads to a decrease in the shear strength, while for $f_f$ > 30% an increase in fine content leads to an increase in the shear strength. This result means that the scalping might underestimate or overestimate the shear strength depending on whether the original fine content $f_f$ is lower or higher than 30%. When the fine content and the global void ratio $e$ are kept constant but the gradation of the coarse
particle size gradation and particle properties? The main objective of this research is to contribute to answering these complex questions for gap-graded cohesionless soils which have not been well studied in the literature. For this study, numerical simulations have been first performed by using the discrete element method (DEM) to better understand factors influencing the behavior of scalped and scalped/replaced soils. Then, experimental tests have been performed on gap-graded soils by using triaxial apparatuses with small, medium and large diameters (50, 100, 200 and 300 mm).

Figure 1 present the grain size distributions (GSD) of a gap-graded soil considered in our study and its corresponding scalped and scalped/replaced ones. It can be seen that the scalping procedure leads to a strong change in the GSD of the tested material; it leads actually to an increase in the fine content of gap-graded soils. For the scalping/replacement procedure, the oversized fraction can be replaced by either the whole remaining coarse fraction (scalped/replaced soil 1) or the class of materials which is closest to the oversized fraction (scalped/replaced soil 2). This procedure allows us to preserve the soil's fine content; however, it reduces the gradation of the coarse fraction.

<figure><img src="image_29.png" /><figcaption>**Figure 1:** GSDs of an original soil with $f_f$ = 30% and its corresponding scalped, scalped/replaced 1 and scalped/replaced 2 soils.</figcaption></figure>

To better understand consequences of the scalping and scalping/replacement procedures, it is vital to understand the effect of the fine content and the gradation of the coarse fraction on the behavior of gap-graded soils. For this purpose, the DEM is used to simulate two series of triaxial tests on bidisperse soils with different fine contents and different gradations of the coarse fraction. The fine particles have diameters $d$ from $d_{min}$ = 1 mm to $d_{max}$ = 2 mm; and the coarse particles have diameters $d$ from $D_{min}$ = 6 mm to $D_{max}$. For the first series, the fine content $f_f$ is varied from 20% to 75% while $D_{max}$ is kept to 12 mm. For the second series, the gradation of the coarse fraction characterized by the ratio $D_{max}/D_{min}$ is varied from 1.5 to 3. Figure 2.a shows the behavior of bidisperse samples with different fine contents $f_f$ when the global void ratio $e$ is controlled to be around 0.32. It can be seen that for $f_f$ < 30% an increase in fine content leads to a decrease in the shear strength, while for $f_f$ > 30% an increase in fine content leads to an increase in the shear strength. This result means that the scalping might underestimate or overestimate the shear strength depending on whether the original fine content $f_f$ is lower or higher than 30%. When the fine content and the global void ratio $e$ are kept constant but the gradation of the coarse
<figure><img src="image_30.png" /><figcaption>**Figure 2:** Stress ratio $q/p$ and volumetric strain $\epsilon_v$ versus axial strain $\epsilon_{11}$ for samples with (a) different fine contents $f_f$ and (b) with different diameter ratios $D_{max}/D_{min}$.</figcaption></figure>

fraction $D_{max}/D_{min}$ is decreased from 3 to 1.5, the shear strength and the dilatancy of bidisperse soils increase as shown in Figure 2.b. This result indicates that the scalping/replacement might lead to an overestimation of the shear strength as the gradation of the coarse fraction is reduced. A micro-mechanical is carried out to explain these macroscopic behaviors.

Experimental triaxial tests are performed on gap-graded soils. Fontainebleau sand is used for the fine fraction with particle diameter $d$ between 0.01 mm and 0.6 mm. To study the effect of the shape and roughness of fine particles, Fontainebleau sand was replaced by glass beads with the same GSD. The coarse fraction is constituted of Palvadeau gravel (subangular shapes) and a natural gravel (angular shapes) and Rhin gravel (rounded shapes). In the following, Fontainebleau sand, glass beads, Palvadeau, natural and Rhin gravels are abbreviated as FS, GB, PG, NG and RG, respectively. For all tested specimens, the ratio $\Phi/d_{max}$ is kept equal to 10; therefore. Each original soil is named by its composition followed by its maximum particle diameter $d_{max}$ in mm. Each scalped or scalped/replaced sample is named by its corresponding original soil's name followed by letter S for scalping, SR1 for scalping/replacement 1 or SR2 for scalping/replacement 2, and followed by its maximum particle diameter ($d_{scalp}$). Several triaxial tests are performed to study effects of different parameters on the behavior of scalped and scalped/replaced soils in comparison with that of the original soil: (i) control compactness parameter, (ii) fine content of the original soil, (iii) particle properties and gradation of the replacement materials and (iv) the gradation of the coarse fraction of the original soil.

The behavior of the scalped soils in comparison with that of the corresponding original soils is shown in Figure 3 for $f_f = 30\%$ (a) and 60% (b). It can be seen that the fine content $f_f$ affects greatly the behavior of scalped and scalped/replaced soils in comparison with that of the original soil. The soil scalped from the original soil with $f_f = 30\%$ (Figure 3.a) is greatly stronger than the original one. On the contrary, for the original soil with $f_f = 60\%$ (Figure 3.b), the scalped soil has
<figure><img src="image_30.png" /><figcaption>**Figure 2:** Stress ratio $q/p$ and volumetric strain $\epsilon_v$ versus axial strain $\epsilon_{11}$ for samples with (a) different fine contents $f_f$ and (b) with different diameter ratios $D_{max}/D_{min}$.</figcaption></figure>

fraction $D_{max}/D_{min}$ is decreased from 3 to 1.5, the shear strength and the dilatancy of bidisperse soils increase as shown in Figure 2.b. This result indicates that the scalping/replacement might lead to an overestimation of the shear strength as the gradation of the coarse fraction is reduced. A micro-mechanical is carried out to explain these macroscopic behaviors.

Experimental triaxial tests are performed on gap-graded soils. Fontainebleau sand is used for the fine fraction with particle diameter $d$ between 0.01 mm and 0.6 mm. To study the effect of the shape and roughness of fine particles, Fontainebleau sand was replaced by glass beads with the same GSD. The coarse fraction is constituted of Palvadeau gravel (subangular shapes) and a natural gravel (angular shapes) and Rhin gravel (rounded shapes). In the following, Fontainebleau sand, glass beads, Palvadeau, natural and Rhin gravels are abbreviated as FS, GB, PG, NG and RG, respectively. For all tested specimens, the ratio $\Phi/d_{max}$ is kept equal to 10; therefore. Each original soil is named by its composition followed by its maximum particle diameter $d_{max}$ in mm. Each scalped or scalped/replaced sample is named by its corresponding original soil's name followed by letter S for scalping, SR1 for scalping/replacement 1 or SR2 for scalping/replacement 2, and followed by its maximum particle diameter ($d_{scalp}$). Several triaxial tests are performed to study effects of different parameters on the behavior of scalped and scalped/replaced soils in comparison with that of the original soil: (i) control compactness parameter, (ii) fine content of the original soil, (iii) particle properties and gradation of the replacement materials and (iv) the gradation of the coarse fraction of the original soil.

The behavior of the scalped soils in comparison with that of the corresponding original soils is shown in Figure 3 for $f_f = 30\%$ (a) and 60% (b). It can be seen that the fine content $f_f$ affects greatly the behavior of scalped and scalped/replaced soils in comparison with that of the original soil. The soil scalped from the original soil with $f_f = 30\%$ (Figure 3.a) is greatly stronger than the original one. On the contrary, for the original soil with $f_f = 60\%$ (Figure 3.b), the scalped soil has
<figure><img src="image_31.png" /><figcaption>**Figure 3:** Stress ratio $q/p$ and volumetric strain $\epsilon_v$ versus axial strain $\epsilon_{11}$ for the original soils with $f_f = (a)$ 30% and (b) 60% and for their scalped (denoted by letter S) and scalped/replaced soils (denoted by letters SR1 and SR2).</figcaption></figure>

almost the same behavior as that of the original soil. The fine content of 30% appears as a turning value under which the scalping induces an underestimation of the shear resistance but over which it leads to an overestimation or good estimation of the shear resistance of coarse-grained soils. The two scalping/replacement methods presented in Figure 1 give different results for the original soil FS30PG40NG30_10 with $f_f = 30\%$ but close results for FS60PG10NG30_10 with $f_f = 60\%$. While the scalping/replacement 1 gives a good estimation of the shear strength of the first soil, the scalping/replacement 2 overestimates its shear strength.

## REFERENCES

- Donaghe, R. T. and Torrey III, V. H. (1985), Strength and deformation properties of earth-rock mixtures., Technical report, Army Engineer Waterway Experiment Station Vicksburg MS Geotechnical Lab.

- El Dine, B. S. (2007), Etude du comportement mécanique de sols grossiers à matrice, PhD thesis, Ecole des Ponts ParisTech.

- Reiffsteck, P., Arbault, J., Sagnard, N., Khay, M., Subrin, D., Chapeau, C. and Levacher, D. (2007), 'Mesures en laboratoire du comportement mécanique des sols hétérogènes', *Bulletin des laboratoires des ponts et chaussées* (268-269).
<figure><img src="image_31.png" /><figcaption>**Figure 3:** Stress ratio $q/p$ and volumetric strain $\epsilon_v$ versus axial strain $\epsilon_{11}$ for the original soils with $f_f = (a)$ 30% and (b) 60% and for their scalped (denoted by letter S) and scalped/replaced soils (denoted by letters SR1 and SR2).</figcaption></figure>

almost the same behavior as that of the original soil. The fine content of 30% appears as a turning value under which the scalping induces an underestimation of the shear resistance but over which it leads to an overestimation or good estimation of the shear resistance of coarse-grained soils. The two scalping/replacement methods presented in Figure 1 give different results for the original soil FS30PG40NG30_10 with $f_f = 30\%$ but close results for FS60PG10NG30_10 with $f_f = 60\%$. While the scalping/replacement 1 gives a good estimation of the shear strength of the first soil, the scalping/replacement 2 overestimates its shear strength.

## REFERENCES

- Donaghe, R. T. and Torrey III, V. H. (1985), Strength and deformation properties of earth-rock mixtures., Technical report, Army Engineer Waterway Experiment Station Vicksburg MS Geotechnical Lab.

- El Dine, B. S. (2007), Etude du comportement mécanique de sols grossiers à matrice, PhD thesis, Ecole des Ponts ParisTech.

- Reiffsteck, P., Arbault, J., Sagnard, N., Khay, M., Subrin, D., Chapeau, C. and Levacher, D. (2007), 'Mesures en laboratoire du comportement mécanique des sols hétérogènes', *Bulletin des laboratoires des ponts et chaussées* (268-269).
# An experimental study on the shear banding phenomenon in granular materials using the orthogonal cutting apparatus

Abhijit HEGDE¹ and Tejas MURTHY¹

¹Indian Institute of Science, Bangalore

**Keywords:** Shear bands, Orthogonal cutting, Granular Materials, Granular chains

## Abstract

When large deformations are imposed on dense granular materials like soil, the resulting strains are confined to narrow regions within the bulk of the material and have been referred to as shear bands or shear zones. The phenomenon of shear banding is commonly encountered during natural disasters such as avalanches, landslides, earthquakes, etc. (Jop et al. (2006), MiDi (2004), Mueth et al. (2000), Nedderman and Laohakul (1980), Rao and Nott (2008))

Though a number of studies do exist which have studied in great detail and have provided great insights into the problem of shear banding (Alshibli and Sture (2000), Arthur et al. (1977), Roscoe (1970), Vardoulakis (1980)), most of the work has been carried out using apparatus such as the bi-axial setup, where the study typically ends with the formation of shear bands in the system. Rather, as a way to study the evolution of these shear bands in time, we have employed the orthogonal cutting configuration, quite popular in the machining literature of metals (Sagapuram et al. (2018)). Cutting configuration allows us to consistently create and sustain shear bands over extended periods of time, hence providing a much larger window for studying the phenomenon of shear banding in granular materials. In our study, we work with three main parameters of the cutting problem, the angle of attack $\alpha$, of the cutting tool, the initial depth of cut D, and the speed of cut, $V_c$ (Fig. 1(b)). Besides varying these parameters corresponding to different boundary conditions, we also perform cutting experiments on a suite of granular materials, varying in size, morphology and density. The experiments are carried out on a CNC machine, inside an aluminum box of dimension 55cm×12cm×2.5cm (Fig. 1(a)). After lowering the cutting tool into the box, granular material is gradually filled via pluviation until the desired depth of cut is reached. Cutting is then carried out a specific rate along with concomitant imaging of the region surrounding the tool, later to be analyzed with a PIV code to extract the deformation parameters such as the velocity of the particles and strain rate. The resolution of captured images is ~ 0.2mm/pixel.

Cutting leads to bifurcation of flow at the tool tip, which results in part of the flow being redirected along the cutting face of the tool while the rest of the material remains unaffected. As a result of this bifurcation, velocity jumps occur within the granular bulk and this leads to the formation of shear bands within the ensemble, as shown in Fig. 1(c), which corresponds to our experiments conducted with angular, quartzitic sand ($d_{50}$ ~ 0.45mm). Formation of localized shear bands is a direct
# An experimental study on the shear banding phenomenon in granular materials using the orthogonal cutting apparatus

Abhijit HEGDE¹ and Tejas MURTHY¹

¹Indian Institute of Science, Bangalore

**Keywords:** Shear bands, Orthogonal cutting, Granular Materials, Granular chains

## Abstract

When large deformations are imposed on dense granular materials like soil, the resulting strains are confined to narrow regions within the bulk of the material and have been referred to as shear bands or shear zones. The phenomenon of shear banding is commonly encountered during natural disasters such as avalanches, landslides, earthquakes, etc. (Jop et al. (2006), MiDi (2004), Mueth et al. (2000), Nedderman and Laohakul (1980), Rao and Nott (2008))

Though a number of studies do exist which have studied in great detail and have provided great insights into the problem of shear banding (Alshibli and Sture (2000), Arthur et al. (1977), Roscoe (1970), Vardoulakis (1980)), most of the work has been carried out using apparatus such as the bi-axial setup, where the study typically ends with the formation of shear bands in the system. Rather, as a way to study the evolution of these shear bands in time, we have employed the orthogonal cutting configuration, quite popular in the machining literature of metals (Sagapuram et al. (2018)). Cutting configuration allows us to consistently create and sustain shear bands over extended periods of time, hence providing a much larger window for studying the phenomenon of shear banding in granular materials. In our study, we work with three main parameters of the cutting problem, the angle of attack $\alpha$, of the cutting tool, the initial depth of cut D, and the speed of cut, $V_c$ (Fig. 1(b)). Besides varying these parameters corresponding to different boundary conditions, we also perform cutting experiments on a suite of granular materials, varying in size, morphology and density. The experiments are carried out on a CNC machine, inside an aluminum box of dimension 55cm×12cm×2.5cm (Fig. 1(a)). After lowering the cutting tool into the box, granular material is gradually filled via pluviation until the desired depth of cut is reached. Cutting is then carried out a specific rate along with concomitant imaging of the region surrounding the tool, later to be analyzed with a PIV code to extract the deformation parameters such as the velocity of the particles and strain rate. The resolution of captured images is ~ 0.2mm/pixel.

Cutting leads to bifurcation of flow at the tool tip, which results in part of the flow being redirected along the cutting face of the tool while the rest of the material remains unaffected. As a result of this bifurcation, velocity jumps occur within the granular bulk and this leads to the formation of shear bands within the ensemble, as shown in Fig. 1(c), which corresponds to our experiments conducted with angular, quartzitic sand ($d_{50}$ ~ 0.45mm). Formation of localized shear bands is a direct
consequence of *mechanical instability* caused in the material as it enters the regime of plastic deformationsDarve et al. (2021), Rice (1976). The corresponding volumetric strain rate maps show strong dilation in the regions corresponding to the shear zones, though the measured dilation angle is found to be much lower than the friction angle of the material. The inclination angle of the shear bands is observed to oscillate in time around a mean value, which roughly corresponds to the Coulomb rupture planes (Roscoe (1970)). Using the inclination angle, we further study the steady state velocity profiles across the shear bands, which are found to be sigmoidal in shape, similar to the wide shear bands of split-couette cells (Fenistein et al. (2004)), and allows us to extract the corresponding length scale associated with the velocity jump/shear band using an error-function fit. To add to this, similar deformation patterns are expected if the experiments are to be carried out for partially or fully saturated soils.

When we perform our cutting experiments on larger glass beads, we notice two major differences with respect to the corresponding deformation fields of sand. Though we observe localized shearing even in the case of glass beads (Fig. 1(d)), we further observe subsequent diffusion of these shear bands, especially as we move farther away from the tool tip. We also observe intense shearing at the free surface of the heap due to avalanching, observed only sporadically in the case of sand. The maximum rate of shearing, measured within the shear bands, also reduces for glass beads. Similar trends are observed for other granular materials, considered in the same size range as the glass beads. The maximum rate of shearing occurring within the shear zones depends primarily on the particle size and follows a power law relationship. Upon re-scaling the velocities with the corresponding shear band widths, we also find that the profiles for different granular materials all collapse onto a single master curve. The width/thickness of shear bands depends on the size of the particle, and less so on other parameters such as the rake angle, depth and the speed of cut (Hegde and Murthy (2022)).

Finally, we perform experiments on granular chains, formed by joining individual metals beads with rigid links made of the same material as the beads(Nickel based alloy). Introducing rigid links between individual particles can be visualized as a means to influence the auto-correlation lengths between the individual grains of the ensemble, defined in the context of force chainsWautier et al. (2017). To this effect, we have employed two kinds of chains, containing 12 and 24 beads each. The effective strain rate field for the case of chains with 12 beads has been presented in Fig. 2(d). We find with chains, that though there is a presence of localized deformation during the initial stages of cutting, with further progress, the resulting deformation is no longer localized but rather diffuse and homogeneous, spreading over the entire ensemble and reaching the boundaries. The effect of chains is even more pronounced when consider the case of 24-beads chains, presented in Fig. 2(b), where no localization is observed even during the initial stages of cutting. There is also significant reduction in the magnitude of the shear rate in the ensemble. One might perhaps employ this strategy of manually altering the auto-correlation length between the particles as means of controlling the width of the shear bands. Thus, one can effectively switch between the localized and diffuse mode of failures in granular materials by providing simple constraints at the particulate level.
consequence of *mechanical instability* caused in the material as it enters the regime of plastic deformationsDarve et al. (2021), Rice (1976). The corresponding volumetric strain rate maps show strong dilation in the regions corresponding to the shear zones, though the measured dilation angle is found to be much lower than the friction angle of the material. The inclination angle of the shear bands is observed to oscillate in time around a mean value, which roughly corresponds to the Coulomb rupture planes (Roscoe (1970)). Using the inclination angle, we further study the steady state velocity profiles across the shear bands, which are found to be sigmoidal in shape, similar to the wide shear bands of split-couette cells (Fenistein et al. (2004)), and allows us to extract the corresponding length scale associated with the velocity jump/shear band using an error-function fit. To add to this, similar deformation patterns are expected if the experiments are to be carried out for partially or fully saturated soils.

When we perform our cutting experiments on larger glass beads, we notice two major differences with respect to the corresponding deformation fields of sand. Though we observe localized shearing even in the case of glass beads (Fig. 1(d)), we further observe subsequent diffusion of these shear bands, especially as we move farther away from the tool tip. We also observe intense shearing at the free surface of the heap due to avalanching, observed only sporadically in the case of sand. The maximum rate of shearing, measured within the shear bands, also reduces for glass beads. Similar trends are observed for other granular materials, considered in the same size range as the glass beads. The maximum rate of shearing occurring within the shear zones depends primarily on the particle size and follows a power law relationship. Upon re-scaling the velocities with the corresponding shear band widths, we also find that the profiles for different granular materials all collapse onto a single master curve. The width/thickness of shear bands depends on the size of the particle, and less so on other parameters such as the rake angle, depth and the speed of cut (Hegde and Murthy (2022)).

Finally, we perform experiments on granular chains, formed by joining individual metals beads with rigid links made of the same material as the beads (Nickel based alloy). Introducing rigid links between individual particles can be visualized as a means to influence the auto-correlation lengths between the individual grains of the ensemble, defined in the context of force chains (Wautier et al. (2017)). To this effect, we have employed two kinds of chains, containing 12 and 24 beads each. The effective strain rate field for the case of chains with 12 beads has been presented in Fig. 2(d). We find with chains, that though there is a presence of localized deformation during the initial stages of cutting, with further progress, the resulting deformation is no longer localized but rather diffuse and homogeneous, spreading over the entire ensemble and reaching the boundaries. The effect of chains is even more pronounced when consider the case of 24-beads chains, presented in Fig. 2(b), where no localization is observed even during the initial stages of cutting. There is also significant reduction in the magnitude of the shear rate in the ensemble. One might perhaps employ this strategy of manually altering the auto-correlation length between the particles as means of controlling the width of the shear bands. Thus, one can effectively switch between the localized and diffuse mode of failures in granular materials by providing simple constraints at the particulate level.
<figure><img src="image_32.png" /><figcaption>**Figure 1:** (a) Experimental Setup. (b) Schematic of the orthogonal cutting problem. (c) Effective strain rate contour field overlapped over the original image showing localized shear deformation in Sand (d) $\dot{\epsilon}_{eff}$ field for glass beads. The shear band tends to diffuse as we move away from the tool tip.</figcaption></figure>

<figure><img src="image_33.png" /><figcaption>**Figure 2:** $\dot{\epsilon}_{eff}$ fields showing diffuse shearing observed in cutting experiments with chains, (a) 12-bead chains (b) 24-bead chains.</figcaption></figure>
<figure><img src="image_32.png" /><figcaption>**Figure 1:** (a) Experimental Setup. (b) Schematic of the orthogonal cutting problem. (c) Effective strain rate contour field overlapped over the original image showing localized shear deformation in Sand (d) $\dot{e}_{eff}$ field for glass beads. The shear band tends to diffuse as we move away from the tool tip.</figcaption></figure>

<figure><img src="image_33.png" /><figcaption>**Figure 2:** $\dot{e}_{eff}$ fields showing diffuse shearing observed in cutting experiments with chains, (a) 12-bead chains (b) 24-bead chains.</figcaption></figure>
REFERENCES

K. A. Alshibli and S. Sture. Shear band formation in plane strain experiments of sand. *Journal of Geotechnical and Geoenvironmental Engineering*, 126(6):495-503, 2000.

J. Arthur, T. Dunstan, Q. Al-Ani, and A. Assadi. Plastic deformation and failure in granular media. *Geotechnique*, 27(1):53-74, 1977.

F. Darve, F. Nicot, A. Wautier, and J. Liu. Slip lines versus shear bands: two competing localization modes. *Mechanics Research Communications*, 114:103603, 2021.

D. Fenistein, J. W. van de Meent, and M. van Hecke. Universal and wide shear zones in granular bulk flow. *Physical Review Letters*, 92(9):094301, 2004.

A. Hegde and T. G. Murthy. Experimental studies on deformation of granular materials during orthogonal cutting. *Granular Matter*, 24(3):1-22, 2022.

P. Jop, Y. Forterre, and O. Pouliquen. A constitutive law for dense granular flows. *Nature*, 441 (7094):727-730, 2006.

G. MiDi. On dense granular flows. *The European Physical Journal E*, 14(4):341-365, 2004.

D. M. Mueth, G. F. Debregeas, G. S. Karczmar, P. J. Eng, S. R. Nagel, and H. M. Jaeger. Signatures of granular microstructure in dense shear flows. *Nature*, 406(6794):385-389, 2000.

R. Nedderman and C. Laohakul. The thickness of the shear zone of flowing granular materials. *Powder Technology*, 25(1):91-100, 1980.

K. K. Rao and P. R. Nott. An introduction to granular flow/k. kesava rao, prabhu r. nott., 2008.

J. R. Rice. Theoretical and applied mechanics. In *Proc. of the 14th IUITAM Congress, North-Holland, Amsterdam, Netherlands*, 1976, pages 207-220, 1976.

K. H. Roscoe. The influence of strains in soil mechanics. *Geotechnique*, 20(2):129-170, 1970.

D. Sagapuram, K. Viswanathan, K. P. Trumble, and S. Chandrasekar. A common mechanism for evolution of single shear bands in large-strain deformation of metals. *Philosophical Magazine*, 98 (36):3267-3299, 2018.

I. Vardoulakis. Shear band inclination and shear modulus of sand in biaxial tests. *International Journal for Numerical and Analytical Methods in Geomechanics*, 4(2):103-119, 1980.

A. Wautier, S. Bonelli, and F. Nicot. Scale separation between grain detachment and grain transport in granular media subjected to an internal flow. *Granular Matter*, 19(2):22, 2017.
REFERENCES

K. A. Alshibli and S. Sture. Shear band formation in plane strain experiments of sand. *Journal of Geotechnical and Geoenvironmental Engineering*, 126(6):495-503, 2000.

J. Arthur, T. Dunstan, Q. Al-Ani, and A. Assadi. Plastic deformation and failure in granular media. *Geotechnique*, 27(1):53-74, 1977.

F. Darve, F. Nicot, A. Wautier, and J. Liu. Slip lines versus shear bands: two competing localization modes. *Mechanics Research Communications*, 114:103603, 2021.

D. Fenistein, J. W. van de Meent, and M. van Hecke. Universal and wide shear zones in granular bulk flow. *Physical Review Letters*, 92(9):094301, 2004.

A. Hegde and T. G. Murthy. Experimental studies on deformation of granular materials during orthogonal cutting. *Granular Matter*, 24(3):1-22, 2022.

P. Jop, Y. Forterre, and O. Pouliquen. A constitutive law for dense granular flows. *Nature*, 441 (7094):727-730, 2006.

G. MiDi. On dense granular flows. *The European Physical Journal E*, 14(4):341-365, 2004.

D. M. Mueth, G. F. Debregeas, G. S. Karczmar, P. J. Eng, S. R. Nagel, and H. M. Jaeger. Signatures of granular microstructure in dense shear flows. *Nature*, 406(6794):385-389, 2000.

R. Nedderman and C. Laohakul. The thickness of the shear zone of flowing granular materials. *Powder Technology*, 25(1):91-100, 1980.

K. K. Rao and P. R. Nott. An introduction to granular flow/k. kesava rao, prabhu r. nott., 2008.

J. R. Rice. Theoretical and applied mechanics. In *Proc. of the 14th IUTAM Congress, North-Holland, Amsterdam, Netherlands*, 1976, pages 207-220, 1976.

K. H. Roscoe. The influence of strains in soil mechanics. *Geotechnique*, 20(2):129-170, 1970.

D. Sagapuram, K. Viswanathan, K. P. Trumble, and S. Chandrasekar. A common mechanism for evolution of single shear bands in large-strain deformation of metals. *Philosophical Magazine*, 98 (36):3267-3299, 2018.

I. Vardoulakis. Shear band inclination and shear modulus of sand in biaxial tests. *International Journal for Numerical and Analytical Methods in Geomechanics*, 4(2):103-119, 1980.

A. Wautier, S. Bonelli, and F. Nicot. Scale separation between grain detachment and grain transport in granular media subjected to an internal flow. *Granular Matter*, 19(2):22, 2017.
# Extending the bifurcation domain concept at the structure scale

Antoine WAUTIER<fnref n="1" />, Anthony MOUYEAU<fnref n="1" />, Richard WAN<fnref n="2" />, François NICOT<fnref n="3" />, and Félix DARVE<fnref n="4" />

## Abstract

The variety of failure modes exhibited by geomaterials are intimately linked to their non-associated plastic behavior. For instance, geomaterials can fail before reaching the plastic limit surface. This is well illustrated by the static liquefaction of loose sand subjected to undrained triaxial tests. The failure of Brumadhino tailing dam in 2019 (Brazil) and Edenville dam in 2020 (USA) are illustrative examples of static liquefaction at the structure scale.

In the wake of the pioneering work of Hill (1958), the second-order work criterion, as proposed by Nicot et al. (2009), Wan et al. (2017), has been shown to be the most general criterion to anticipate failure in geomaterials. This criterion enables us to define a bifurcation domain corresponding to a set of states for which there exists a potential for failure if i) the material is loaded in a critical direction corresponding to the vanishing of the second-order work, ii) the mode of control allows for an inertial transition through a sudden burst of kinetic energy. However so far, this notion of conditional failure at the material point scale cannot be assessed at the structure scale. Only the works of Laouafa et al. (2011), Prunier et al. (2009) provided some clues to answer this question at the scale of an engineering structure, by observing the vanishing of the determinant of the global stiffness matrix that relates to unknown degrees of freedom of a FEM problem. Such an approach is however computationally demanding because the size of the matrix is huge.

In this work, we propose to extend the concept of bifurcation domain from the material point scale to the engineering structure scale. At the structure scale, we define the concept of bifurcation zone that corresponds to the set of material points belonging to their bifurcation domain, i.e. the set of points for which conditional failure may occur. More details can be found in Wautier et al. (2022).

FEM simulations are conducted for a dam made of a non-associated Drucker-Prager elasto-plastic material. For each point of the dam, strain control directional analyses are conducted from an analytical point of view to detect which points of the dam are in the bifurcation domain from a material scale viewpoint. By analyzing the directions of the instability cone thus detected, the existence of potential internal failure surfaces are estimated.

For the points located in the bifurcation zone, there exists a set of incremental directions that will lead to material failure. The cone direction corresponds to an incremental strain tensor $d\boldsymbol{e}_{\text{cone}}$ which is characterized by its principal directions and its principal values. Positive values corresponds to extension along the corresponding

<footnote n="1">¹INRAE, Aix-Marseille Université, UMR RECOVER, Aix-en-Provence, France</footnote>

<footnote n="2">²Civil Engineering Department, University of Calgary, Calgary, AB, Canada</footnote>

<footnote n="3">³Université Savoie Mont Blanc, ISTerre, Chambéry, France</footnote>

<footnote n="4">⁴Université Grenoble Alpes, CNRS, G-INP, Laboratoire 3SR UMR5521, Grenoble, France</footnote>
# Extending the bifurcation domain concept at the structure scale

Antoine WAUTIER<fnref n="1" />, Anthony MOUYEAU<fnref n="1" />, Richard WAN<fnref n="2" />, François NICOT<fnref n="3" />, and Félix DARVE<fnref n="4" />

## Abstract

The variety of failure modes exhibited by geomaterials are intimately linked to their non-associated plastic behavior. For instance, geomaterials can fail before reaching the plastic limit surface. This is well illustrated by the static liquefaction of loose sand subjected to undrained triaxial tests. The failure of Brumadhino tailing dam in 2019 (Brazil) and Edenville dam in 2020 (USA) are illustrative examples of static liquefaction at the structure scale.

In the wake of the pioneering work of Hill (1958), the second-order work criterion, as proposed by Nicot et al. (2009), Wan et al. (2017), has been shown to be the most general criterion to anticipate failure in geomaterials. This criterion enables us to define a bifurcation domain corresponding to a set of states for which there exists a potential for failure if i) the material is loaded in a critical direction corresponding to the vanishing of the second-order work, ii) the mode of control allows for an inertial transition through a sudden burst of kinetic energy. However so far, this notion of conditional failure at the material point scale cannot be assessed at the structure scale. Only the works of Laouafa et al. (2011), Prunier et al. (2009) provided some clues to answer this question at the scale of an engineering structure, by observing the vanishing of the determinant of the global stiffness matrix that relates to unknown degrees of freedom of a FEM problem. Such an approach is however computationally demanding because the size of the matrix is huge.

In this work, we propose to extend the concept of bifurcation domain from the material point scale to the engineering structure scale. At the structure scale, we define the concept of bifurcation zone that corresponds to the set of material points belonging to their bifurcation domain, i.e. the set of points for which conditional failure may occur. More details can be found in Wautier et al. (2022).

FEM simulations are conducted for a dam made of a non-associated Drucker-Prager elasto-plastic material. For each point of the dam, strain control directional analyses are conducted from an analytical point of view to detect which points of the dam are in the bifurcation domain from a material scale viewpoint. By analyzing the directions of the instability cone thus detected, the existence of potential internal failure surfaces are estimated.

For the points located in the bifurcation zone, there exists a set of incremental directions that will lead to material failure. The cone direction corresponds to an incremental strain tensor $d\boldsymbol{e}_{\text{cone}}$ which is characterized by its principal directions and its principal values. Positive values corresponds to extension along the corresponding

<footnote n="1">¹INRAE, Aix-Marseille Université, UMR RECOVER, Aix-en-Provence, France</footnote>

<footnote n="2">²Civil Engineering Department, University of Calgary, Calgary, AB, Canada</footnote>

<footnote n="3">³Université Savoie Mont Blanc, ISTerre, Chambéry, France</footnote>

<footnote n="4">⁴Université Grenoble Alpes, CNRS, G-INP, Laboratoire 3SR UMR5521, Grenoble, France</footnote>
<figure><img src="image_34.png" /><figcaption>**Figure 1:** Bifurcation domain of a dam. For each material point in the bifurcation domain, the cone direction is represented. Such representation corresponds to local deformation mechanisms likely to trigger underlying material instability.</figcaption></figure>

principal direction. A representation of the dam bifurcation domain together with the local deformation mechanisms likely to trigger material instability (i.e. the cone directions) is given in Figure 1.

For the dam considered, the bifurcation domain is quite extensive and looks consistent with the shapes of failure observed in situ. However, one can notice that most of the cone directions have a non negligible off plane component, making them not compatible with plane strain mechanisms. If we restrict the admissible failure mechanisms to plane strain conditions, the bifurcation domain is largely reduced as illustrated in Figure 2. It no longer spans across the dam body.

As a result, we can conclude from these graphs that there exists some plausible mechanisms that can lead to the failure of the dam considered; provided that the plane strain assumption is lifted. These mechanisms cannot develop fully without allowing an extension of the dam in the out-of-plane direction.

The present study provides an illustrative example on how to use the second-order work criterion at the scale of an engineering structure. Such a contribution proposes a modern and innovative view of failure in geomaterials at the engineering scale.

## REFERENCES

* R. Hill. A general theory of uniqueness and stability in elastic-plastic solids. *Journal of the Mechanics and Physics of Solids*, 6(3):236–249, 1958.

* F. Laouafa, F. Prunier, A. Daouadji, H. A. Gali, and F. Darve. Stability in geomechanics, experimental and numerical analyses. *International Journal for Numerical and Analytical Methods in Geomechanics*, 35(2):112–139, 2011.
<figure><img src="image_34.png" /><figcaption>**Figure 1:** Bifurcation domain of a dam. For each material point in the bifurcation domain, the cone direction is represented. Such representation corresponds to local deformation mechanisms likely to trigger underlying material instability.</figcaption></figure>

principal direction. A representation of the dam bifurcation domain together with the local deformation mechanisms likely to trigger material instability (i.e. the cone directions) is given in Figure 1.

For the dam considered, the bifurcation domain is quite extensive and looks consistent with the shapes of failure observed in situ. However, one can notice that most of the cone directions have a non negligible off plane component, making them not compatible with plane strain mechanisms. If we restrict the admissible failure mechanisms to plane strain conditions, the bifurcation domain is largely reduced as illustrated in Figure 2. It no longer spans across the dam body.

As a result, we can conclude from these graphs that there exists some plausible mechanisms that can lead to the failure of the dam considered; provided that the plane strain assumption is lifted. These mechanisms cannot develop fully without allowing an extension of the dam in the out-of-plane direction.

The present study provides an illustrative example on how to use the second-order work criterion at the scale of an engineering structure. Such a contribution proposes a modern and innovative view of failure in geomaterials at the engineering scale.

## REFERENCES

R. Hill. A general theory of uniqueness and stability in elastic-plastic solids. *Journal of the Mechanics and Physics of Solids*, 6(3):236–249, 1958.

F. Laouafa, F. Prunier, A. Daouadji, H. A. Gali, and F. Darve. Stability in geomechanics, experimental and numerical analyses. *International Journal for Numerical and Analytical Methods in Geomechanics*, 35(2):112–139, 2011.
<figure><img src="image_35.png" /><figcaption>**Figure 2:** Bifurcation domain of a dam restricted to plane strain perturbations. For each material point in the bifurcation domain, the plane strain cone direction is represented as in Figure 1.</figcaption></figure>

F. Nicot, L. Sibille, and F. Darve. Bifurcation in granular materials: An attempt for a unified framework. *International Journal of Solids and Structures*, 46(22-23):3938-3947, 2009.

F. Prunier, F. Laouafa, S. Lignon, and F. Darve. Bifurcation modeling in geomaterials: from the second-order work criterion to spectral analyses. *International Journal for Numerical and Analytical Methods in Geomechanics*, 33(9):1169-1202, 2009.

R. Wan, F. Nicot, and F. Darve. *Failure in geomaterials: a contemporary treatise*. Elsevier, 2017.

A. Wautier, A. Mouyeaux, F. Nicot, R. Wan, and F. Darve. Mechanical stability analysis of engineering structures with use of the bifurcation domain concept. In *Multiscale Processes of Instability, Deformation and Fracturing in Geomaterials: Proceedings of 12th International Workshop on Bifurcation and Degradation in Geomechanics*, pages 3-12. Springer, 2022. doi: 10.1007/978-3-031-22213-9_1.
<figure><img src="image_35.png" /><figcaption>**Figure 2:** Bifurcation domain of a dam restricted to plane strain perturbations. For each material point in the bifurcation domain, the plane strain cone direction is represented as in Figure 1.</figcaption></figure>

F. Nicot, L. Sibille, and F. Darve. Bifurcation in granular materials: An attempt for a unified framework. *International Journal of Solids and Structures*, 46(22-23):3938-3947, 2009.

F. Prunier, F. Laouafa, S. Lignon, and F. Darve. Bifurcation modeling in geomaterials: from the second-order work criterion to spectral analyses. *International Journal for Numerical and Analytical Methods in Geomechanics*, 33(9):1169-1202, 2009.

R. Wan, F. Nicot, and F. Darve. *Failure in geomaterials: a contemporary treatise*. Elsevier, 2017.

A. Wautier, A. Mouyeaux, F. Nicot, R. Wan, and F. Darve. Mechanical stability analysis of engineering structures with use of the bifurcation domain concept. In *Multiscale Processes of Instability, Deformation and Fracturing in Geomaterials: Proceedings of 12th International Workshop on Bifurcation and Degradation in Geomechanics*, pages 3-12. Springer, 2022. doi: 10.1007/978-3-031-22213-9_1.
# Capturing onset of bifurcation and post-localisation behaviour in constitutive modelling of partially saturated soils.

Dat G. PHAN<fnref n="1" />, Giang D. NGUYEN<fnref n="1" />*, Ha H. BUI<fnref n="2" />, and Terry BENNETT<fnref n="1" />

**Keywords:** partially saturated soils, bifurcation, localised failure, constitutive modelling

## Abstract

Localised failure in the form of shear bands induce inhomogeneous deformation and distribution of degree of saturation, leading to size effects on macro responses of partially saturated soil samples. The responses inside and outside shear bands are very different: strong inelastic behaviour inside shear bands, while the zone outside shear bands usually undergo relatively small deformation. Classical constitutive models for soils in general, and for partially saturated soils in particular, are formulated based on the assumption of homogenous deformation that is no longer valid once localisation occurs. Therefore, they are not able to describe correctly post-localisation behaviour, due to the lack of a length scale related to the size of the localisation band.

We develop a strategy to embed meso scale behaviour of the localisation band, which is usually missing in classical constitutive models, in a new constitutive structure. This provides a way to connect meso with macro scales, and describe correctly, at the constitutive level, post-localisation behaviour of geomaterials in general and partially saturated soils in particular. The proposed enhancement allows describing correctly both pre- and post-localisation responses of partially saturated soils. This enhancement to the constitutive structure is activated when localisation of deformation is detected, based on the loss of positiveness of the acoustic tensor. Both onset of bifurcation and irreversible hydro-mechanical behaviour inside shear bands are described by a constitutive model possessing a single yield surface in stress-suction space governing the interdependence of two evolution rules for plastic strain and irrecoverable saturation degree. This allows simultaneous activation and evolution of both mechanical and hydraulic yielding responses in different loading and saturation regimes, and facilitates mapping of onset of bifurcation and associated shear band orientations on a single yield surface. The promising features of the proposed approach are demonstrated through experimental data on partially saturated soils under drained triaxial tests.

<footnote n="1">¹School of Architecture and Civil Engineering, University of Adelaide, Australia</footnote>

<footnote n="2">²Department of Civil Engineering, Monash University, Australia</footnote>

<footnote>*Corresponding author E-mail: g.nguyen@adelaide.edu.au</footnote>
# Capturing onset of bifurcation and post-localisation behaviour in constitutive modelling of partially saturated soils.

Dat G. PHAN<fnref n="1" />, Giang D. NGUYEN<fnref n="1" />*, Ha H. BUI<fnref n="2" />, and Terry BENNETT<fnref n="1" />

**Keywords:** partially saturated soils, bifurcation, localised failure, constitutive modelling

## Abstract

Localised failure in the form of shear bands induce inhomogeneous deformation and distribution of degree of saturation, leading to size effects on macro responses of partially saturated soil samples. The responses inside and outside shear bands are very different: strong inelastic behaviour inside shear bands, while the zone outside shear bands usually undergo relatively small deformation. Classical constitutive models for soils in general, and for partially saturated soils in particular, are formulated based on the assumption of homogenous deformation that is no longer valid once localisation occurs. Therefore, they are not able to describe correctly post-localisation behaviour, due to the lack of a length scale related to the size of the localisation band.

We develop a strategy to embed meso scale behaviour of the localisation band, which is usually missing in classical constitutive models, in a new constitutive structure. This provides a way to connect meso with macro scales, and describe correctly, at the constitutive level, post-localisation behaviour of geomaterials in general and partially saturated soils in particular. The proposed enhancement allows describing correctly both pre- and post-localisation responses of partially saturated soils. This enhancement to the constitutive structure is activated when localisation of deformation is detected, based on the loss of positiveness of the acoustic tensor. Both onset of bifurcation and irreversible hydro-mechanical behaviour inside shear bands are described by a constitutive model possessing a single yield surface in stress-suction space governing the interdependence of two evolution rules for plastic strain and irrecoverable saturation degree. This allows simultaneous activation and evolution of both mechanical and hydraulic yielding responses in different loading and saturation regimes, and facilitates mapping of onset of bifurcation and associated shear band orientations on a single yield surface. The promising features of the proposed approach are demonstrated through experimental data on partially saturated soils under drained triaxial tests.

<footnote n="1">¹School of Architecture and Civil Engineering, University of Adelaide, Australia</footnote>

<footnote n="2">²Department of Civil Engineering, Monash University, Australia</footnote>

<footnote>*Corresponding author E-mail: g.nguyen@adelaide.edu.au</footnote>
# Observation of grain-scale erosion/sedimentation process using DEM/CFD coupling simulation

Takashi MATSUSHIMA¹

¹Institute of Systems and Information Engineering, University of Tsukuba

## Abstract

How is geomorphology formed and how does it evolves? This is the central issue of earth & planetary sciences. It is also important in various engineering fields such as civil and environmental engineering, energy resource engineering and agricultural engineering to understand and predict the geomorphological processes quantitatively. On the earth, the surface water flow is a primary driving force for erosion, transportation, and sedimentation of geomaterials. This study attempts to elucidate the grain scale mechanics of sediment transport in riverbed using DEM/CFD simulation.
# Observation of grain-scale erosion/sedimentation process using DEM/CFD coupling simulation

Takashi MATSUSHIMA¹

¹Institute of Systems and Information Engineering, University of Tsukuba

## Abstract

How is geomorphology formed and how does it evolves? This is the central issue of earth & planetary sciences. It is also important in various engineering fields such as civil and environmental engineering, energy resource engineering and agricultural engineering to understand and predict the geomorphological processes quantitatively. On the earth, the surface water flow is a primary driving force for erosion, transportation, and sedimentation of geomaterials. This study attempts to elucidate the grain scale mechanics of sediment transport in riverbed using DEM/CFD simulation.
DEM capabilities with Polyhedral and Level Set shape descriptors

Jérôme Duriez¹, Tarek Mohamed¹, Stéphane Bonelli¹, Frédéric Golay², Cédric Galusinski²

¹INRAE, Aix Marseille Univ, RECOVER, Aix-en-Provence, France

²IMATH, Université de Toulon, CS 60584 83041 Toulon Cedex 9, France

jerome.duriez@inrae.fr

**Keywords:** Granular mechanics; Level Set-Discrete Element Method (LS-DEM); Polyhedra

**Abstract**

Geometrically realistic particles are getting more and more available in discrete (DEM) simulations thanks to a variety of dedicated DEM implementations, enabling one to properly address the known influence of particle shape in geomechanics [1]. In addition to the multi-sphere approach, whereby an appropriate number of overlapping spheres are clumped together and act as one rigid body, one can mention the possibility to define a particle's surface as a polyhedra with its set of vertices, edges and facets or through a Level Set description [2]. The latter relies on a discrete field for the signed distance function to the surface at hand, whose zero level set implicitly describes that surface. For the purpose of contact treatment, a surface discretization in terms of boundary nodes is deduced from the distance field and also enters the method as a second ingredient.

Firstly, the capabilities offered by polyhedral particles lead to propose and validate against experimental data a DEM model for Toyoura sand for a wide range of loading conditions (drained and undrained triaxial compression and extension) and various initial void ratios and confining pressures [3]. Polyhedral Discrete Elements are simply but efficiently defined in 3D from a 2D micrograph of Toyoura sand (Fig. 1). In addition to show satisfactory predictive abilities, the model, through a parametric analysis, illustrates once more the importance for a proper calibration in DEM of initial fabric, in addition to shape or contact parameters.

Second, the polyhedral approach is compared with the Level Set one, in the case of implementations into the YADE open-source code [4,5,6]. While being orders of magnitude heavier than for classical spheres [7], Level Set (LS) computational costs are shown to be possibly lighter than those of polyhedra and may be further optimized be it in terms of memory or time cost. LS-based simulations finally address the mechanical behaviour of superquadric particles (superellipsoids, Fig. 2).

**References**

[1] G.-C. Cho, J. Dodds and J. C. Santamarina, Particle shape effects on packing density, stiffness, and strength: natural and crushed sands, *J. of Geotechnical and Geoenvironmental Engineering*, 132(5), 2006

[2] R. Kawamoto, E. Andò, G. Viggiani, and J. E. Andrade, Level set discrete element method for three-dimensional computations with triaxial case study. *J. of the Mechanics and Physics of Solids*, 91, 2016

[3] T. Mohamed, J. Duriez, G. Veylon and L. Peyras, DEM models using direct and indirect shape descriptions for Toyoura sand along monotonous loading paths, *Computers and Geotechnics*, 142, 2022

[4] V. Šmilauer et al., Yade Documentation 3rd ed. *The Yade Project* (http://yade-dem.org/doc/)

[5] J. Eliáš, Simulation of railway ballast using crushable polyhedral particles, *Powder Technology*, 264, 2014
DEM capabilities with Polyhedral and Level Set shape descriptors

Jérôme Duriez¹, Tarek Mohamed¹, Stéphane Bonelli¹, Frédéric Golay², Cédric Galusinski²

¹INRAE, Aix Marseille Univ, RECOVER, Aix-en-Provence, France

²IMATH, Université de Toulon, CS 60584 83041 Toulon Cedex 9, France

jerome.duriez@inrae.fr

**Keywords:** Granular mechanics; Level Set-Discrete Element Method (LS-DEM); Polyhedra

**Abstract**

Geometrically realistic particles are getting more and more available in discrete (DEM) simulations thanks to a variety of dedicated DEM implementations, enabling one to properly address the known influence of particle shape in geomechanics [1]. In addition to the multi-sphere approach, whereby an appropriate number of overlapping spheres are clumped together and act as one rigid body, one can mention the possibility to define a particle's surface as a polyhedra with its set of vertices, edges and facets or through a Level Set description [2]. The latter relies on a discrete field for the signed distance function to the surface at hand, whose zero level set implicitly describes that surface. For the purpose of contact treatment, a surface discretization in terms of boundary nodes is deduced from the distance field and also enters the method as a second ingredient.

Firstly, the capabilities offered by polyhedral particles lead to propose and validate against experimental data a DEM model for Toyoura sand for a wide range of loading conditions (drained and undrained triaxial compression and extension) and various initial void ratios and confining pressures [3]. Polyhedral Discrete Elements are simply but efficiently defined in 3D from a 2D micrograph of Toyoura sand (Fig. 1). In addition to show satisfactory predictive abilities, the model, through a parametric analysis, illustrates once more the importance for a proper calibration in DEM of initial fabric, in addition to shape or contact parameters.

Second, the polyhedral approach is compared with the Level Set one, in the case of implementations into the YADE open-source code [4,5,6]. While being orders of magnitude heavier than for classical spheres [7], Level Set (LS) computational costs are shown to be possibly lighter than those of polyhedra and may be further optimized be it in terms of memory or time cost. LS-based simulations finally address the mechanical behaviour of superquadric particles (superellipsoids, Fig. 2).

**References**

[1] G.-C. Cho, J. Dodds and J. C. Santamarina, Particle shape effects on packing density, stiffness, and strength: natural and crushed sands, *J. of Geotechnical and Geoenvironmental Engineering*, 132(5), 2006

[2] R. Kawamoto, E. Andò, G. Viggiani, and J. E. Andrade, Level set discrete element method for three-dimensional computations with triaxial case study. *J. of the Mechanics and Physics of Solids*, 91, 2016

[3] T. Mohamed, J. Duriez, G. Veylon and L. Peyras, DEM models using direct and indirect shape descriptions for Toyoura sand along monotonous loading paths, *Computers and Geotechnics*, 142, 2022

[4] V. Šmilauer et al., Yade Documentation 3rd ed. *The Yade Project* (http://yade-dem.org/doc/)

[5] J. Eliáš, Simulation of railway ballast using crushable polyhedral particles, *Powder Technology*, 264, 2014
[6] J. Duriez and C. Galusinski, A Level Set-Discrete Element Method in YADE for numerical, micro-scale, geomechanics with refined grain shapes, *Computers & Geosciences*, 157, 2021

[7] J. Duriez and S. Bonelli, Precision and computational costs of Level Set-Discrete Element Method (LS-DEM) with respect to DEM, *Comp. & Geotechnics*, 134, 2021

[8] B. Li, Effect of fabric anisotropy on the dynamic mechanical behavior of granular materials. Ph.D. thesis, Case Western Reserve University, 2011

## Figures

<figure><img src="image_36.png" /><figcaption>**Figure 1:** Polyhedral model for Toyoura sand [3] Left: micrograph [8] with corresponding Discrete Elements. Right: corresponding numerical packing</figcaption></figure>

<figure><img src="image_37.png" /><figcaption>**Figure 2:** Level Set description of superellipsoids</figcaption></figure>
[6] J. Duriez and C. Galusinski, A Level Set-Discrete Element Method in YADE for numerical, micro-scale, geomechanics with refined grain shapes, *Computers & Geosciences*, 157, 2021

[7] J. Duriez and S. Bonelli, Precision and computational costs of Level Set-Discrete Element Method (LS-DEM) with respect to DEM, *Comp. & Geotechnics*, 134, 2021

[8] B. Li, Effect of fabric anisotropy on the dynamic mechanical behavior of granular materials. Ph.D. thesis, Case Western Reserve University, 2011

## Figures

<figure><img src="image_36.png" /><figcaption>**Figure 1:** Polyhedral model for Toyoura sand [3] Left: micrograph [8] with corresponding Discrete Elements. Right: corresponding numerical packing</figcaption></figure>

<figure><img src="image_37.png" /><figcaption>**Figure 2:** Level Set description of superellipsoids</figcaption></figure>
# Influence of stress state on the initiation and the whole development of suffusion

Bikram OLI¹, Didier MAROT¹, Rachel GELET¹, Fateh BENDAHMANE¹

¹GeM Institut UMR CNRS 6183, Nantes University

58 rue Michel Ange, BP 420, 44606, Saint-Nazaire Cedex, France

Corresponding author: Didier MAROT

didier.marot@univ-nantes.fr

**Key words:** suffusion, stress state, triaxial erosion apparatus, gap-graded soils, hydraulic gradient, expended energy

According to Fell and Fry (2013), four mechanisms of internal erosion can be distinguished: concentrated leak erosion, backward erosion, contact erosion and suffusion. This study focuses on suffusion, which is the selective erosion of fine particles under the effect of seepage flow within the matrix of coarser particles.

For the characterization of suffusion occurrence, Garner and Fannin (2010) as well as Fell and Fry (2013) proposed three criteria to be satisfied. The first criterion concerns the geometry of the porous media and points out that the size of the constrictions between the coarser particles must be larger than the size of detached particles. The second criterion is related to the stress applied on the fine fraction, and the third one considers the hydraulic load. According to Fell and Fry (2013), this third criterion is related to the porewater velocity, whereas Garner and Fannin (2010) considered the variations in porewater pressure. In this context and from results of interface erosion tests, Marot et al. (2011) proposed an analysis based on the total flow power, which is the summation of the power transferred from the fluid to the solid particles, with the power dissipated by the viscous stresses in the bulk. According to Sibille et al. (2015), the first term is negligible in from of the second in the case of suffusion, and the authors suggested characterizing the fluid loading from the total flow power (expressed in W), by the following equation:

$$P_{\text{flow}} = \left( \frac{\gamma_w}{\Delta z + \Delta P} \right) Q \quad (1)$$

where Q is the fluid flow rate (in m³.s⁻¹), γw is the specific weight of water (in N.m⁻³), ΔP = Pₐ - Pₑ is the pressure drop (in Pa) between the upstream section A and the downstream section B, Δz = zₐ - zₑ, and zₐ and zₑ are altitudes (in m) of sections A and B, respectively. Δz > 0 if the flow is downward and Δz < 0 if opposite. The total flow power is equal to QΔP if the flow is horizontal.

The expended energy Eflow (expressed in J) is the time integration of the instantaneous power dissipated by water seepage and is computed until reaching the steady state (i.e., the hydraulic conductivity tends to stabilize and the erosion rate tends to decrease). For the same duration, the cumulative eroded dry mass is determined (expressed in kg) and the erosion resistance index is expressed by:

$$I_a = - \log \left( \frac{\text{Cumulative eroded dry mass}}{E_{\text{flow}}} \right) \quad (2)$$
# Influence of stress state on the initiation and the whole development of suffusion

Bikram OLI¹, Didier MAROT¹, Rachel GELET¹, Fateh BENDAHMANE¹

¹GeM Institut UMR CNRS 6183, Nantes University

58 rue Michel Ange, BP 420, 44606, Saint-Nazaire Cedex, France

Corresponding author: Didier MAROT

didier.marot@univ-nantes.fr

**Key words:** suffusion, stress state, triaxial erosion apparatus, gap-graded soils, hydraulic gradient, expended energy

According to Fell and Fry (2013), four mechanisms of internal erosion can be distinguished: concentrated leak erosion, backward erosion, contact erosion and suffusion. This study focuses on suffusion, which is the selective erosion of fine particles under the effect of seepage flow within the matrix of coarser particles.

For the characterization of suffusion occurrence, Garner and Fannin (2010) as well as Fell and Fry (2013) proposed three criteria to be satisfied. The first criterion concerns the geometry of the porous media and points out that the size of the constrictions between the coarser particles must be larger than the size of detached particles. The second criterion is related to the stress applied on the fine fraction, and the third one considers the hydraulic load. According to Fell and Fry (2013), this third criterion is related to the porewater velocity, whereas Garner and Fannin (2010) considered the variations in porewater pressure. In this context and from results of interface erosion tests, Marot et al. (2011) proposed an analysis based on the total flow power, which is the summation of the power transferred from the fluid to the solid particles, with the power dissipated by the viscous stresses in the bulk. According to Sibille et al. (2015), the first term is negligible in from of the second in the case of suffusion, and the authors suggested characterizing the fluid loading from the total flow power (expressed in W), by the following equation:

$$P_{\text{flow}} = \left( \frac{\gamma_w}{\Delta z + \Delta P} \right) Q \qquad (1)$$

where Q is the fluid flow rate (in m³.s⁻¹), γw is the specific weight of water (in N.m⁻³), ΔP = Pₐ - Pₑ is the pressure drop (in Pa) between the upstream section A and the downstream section B, Δz = zₐ - zₑ, and zₐ and zₑ are altitudes (in m) of sections A and B, respectively. Δz > 0 if the flow is downward and Δz < 0 if opposite. The total flow power is equal to QΔP if the flow is horizontal.

The expended energy Eflow (expressed in J) is the time integration of the instantaneous power dissipated by water seepage and is computed until reaching the steady state (i.e., the hydraulic conductivity tends to stabilize and the erosion rate tends to decrease). For the same duration, the cumulative eroded dry mass is determined (expressed in kg) and the erosion resistance index is expressed by:

$$I_a = - \log \left( \frac{\text{Cumulative eroded dry mass}}{E_{\text{flow}}} \right) \qquad (2)$$
Depending on the values of the erosion resistance index $I_0$, Marot et al. (2016) proposed six categories of suffusion susceptibility from highly erodable to highly resistant.

The two first criteria are mainly related to the grain size distribution and proposals of various geometric assessment methods exist in the literature, mostly based on the particle size distribution. However, the stress state conditions also can influence the suffusion susceptibility. Chang and Zhang (2013) investigated the effect of soil stress state on the values of two hydraulic gradients: an initiation hydraulic gradient associated with a slight increase of the hydraulic conductivity (approach initially proposed by Skempton and Brogan, 1994) and a skeleton-deformation-gradient, referred to sudden increase in eroded soil mass, hydraulic conductivity and deformation. Both hydraulic gradients were slightly lower under triaxial extension in comparison with triaxial compression. However, the hydraulic loading path was not systematically the same for all performed tests, and yet this parameter plays an essential role on the initiation and development of suffusion (Rochim et al., 2017). Since, the small observed differences in hydraulic gradients could be due to the hydraulic loading path so that the influence of the stress state on suffusion needs to be further explored.

In consequence, the objective of the current study is to assess the effect of stress states on the initiation and the development of suffusion by performing suffusion tests on gap-graded soils, under a controlled hydraulic loading path and controlled stress states.

An erosion apparatus was designed to independently control the stress state and the hydraulic loading path. The top cap and base pedestal of the triaxial cell are modified respectively: (i) diffuse the injected flow and (ii) collect the eroded soil particles in a bottom effluent tank. In the effluent tank, a rotational mechanism comprises 16 beakers to collect eroded particles at different times. The discrete eroded mass measurement is favored over a continuous measurement, in order to eliminate the pressure perturbations induced by opening and closing valves. The confining stress is regulated while a deviatoric load is imposed on the specimen via an axial loading jack. The differential head between the upstream and downstream specimen's sections is controlled by an automated upstream reservoir and a stationary downstream reservoir. At the overflow outlet, water is continuously weighted to determine the flow rate.

The four tested soils are all gap-graded and composed by a mixture of sand S1 and gravel G3, both marketed at the Sablière Palvadeau quarry (France). The tested sand percentages are: FC = 15%, 25%, 35% and 40%. Figure 1(a) shows the grain size distribution of the mixture composed by 25% of sand and 75% of gravel (named FC25).

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

<figure><img src="image_39.png" /><figcaption>Figure 1: (a) grain size distribution of FC25, (b) time evolution of the applied hydraulic gradient.</figcaption></figure>

Each experiment was carried out in four phases: (i) saturation (ii) stress state application (iii) erosion and (iv) post-suffusion particle size distribution. First, each specimen was prepared by moist mixing to prevent soil segregation. Thereafter, each tested specimen (100 mm in diameter and 200 mm in height) was compacted in four layers by using the under-compaction moist tamping
Depending on the values of the erosion resistance index $I_α$, Marot et al. (2016) proposed six categories of suffusion susceptibility from highly erodable to highly resistant.

The two first criteria are mainly related to the grain size distribution and proposals of various geometric assessment methods exist in the literature, mostly based on the particle size distribution. However, the stress state conditions also can influence the suffusion susceptibility. Chang and Zhang (2013) investigated the effect of soil stress state on the values of two hydraulic gradients: an initiation hydraulic gradient associated with a slight increase of the hydraulic conductivity (approach initially proposed by Skempton and Brogan, 1994) and a skeleton-deformation-gradient, referred to sudden increase in eroded soil mass, hydraulic conductivity and deformation. Both hydraulic gradients were slightly lower under triaxial extension in comparison with triaxial compression. However, the hydraulic loading path was not systematically the same for all performed tests, and yet this parameter plays an essential role on the initiation and development of suffusion (Rochim et al., 2017). Since, the small observed differences in hydraulic gradients could be due to the hydraulic loading path so that the influence of the stress state on suffusion needs to be further explored.

In consequence, the objective of the current study is to assess the effect of stress states on the initiation and the development of suffusion by performing suffusion tests on gap-graded soils, under a controlled hydraulic loading path and controlled stress states.

An erosion apparatus was designed to independently control the stress state and the hydraulic loading path. The top cap and base pedestal of the triaxial cell are modified respectively: (i) diffuse the injected flow and (ii) collect the eroded soil particles in a bottom effluent tank. In the effluent tank, a rotational mechanism comprises 16 beakers to collect eroded particles at different times. The discrete eroded mass measurement is favored over a continuous measurement, in order to eliminate the pressure perturbations induced by opening and closing valves. The confining stress is regulated while a deviatoric load is imposed on the specimen via an axial loading jack. The differential head between the upstream and downstream specimen's sections is controlled by an automated upstream reservoir and a stationary downstream reservoir. At the overflow outlet, water is continuously weighted to determine the flow rate.

The four tested soils are all gap-graded and composed by a mixture of sand S1 and gravel G3, both marketed at the Sablière Palvadeau quarry (France). The tested sand percentages are: FC = 15%, 25%, 35% and 40%. Figure 1(a) shows the grain size distribution of the mixture composed by 25% of sand and 75% of gravel (named FC25).

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

<figure><img src="image_39.png" /><figcaption>Figure 1: (a) grain size distribution of FC25, (b) time evolution of the applied hydraulic gradient.</figcaption></figure>

Each experiment was carried out in four phases: (i) saturation (ii) stress state application (iii) erosion and (iv) post-suffusion particle size distribution. First, each specimen was prepared by moist mixing to prevent soil segregation. Thereafter, each tested specimen (100 mm in diameter and 200 mm in height) was compacted in four layers by using the under-compaction moist tamping
technique. Afterward, the sample was confined with a confining pressure of 20 kPa to prevent preferential flow between the specimen and the membrane during the saturation phase. Subsequently, CO₂ was injected at a slow rate for 20 minutes and tap water was upwardly injected at the rate of 0.25 mm/min to limit erosion prior to erosion phase. Upon completion of the saturation, the sample was consolidated to the target mean effective stress p' of 70 kPa and, afterwards, vertically loaded to a target value of deviatoric stress q. The corresponding values of shear stress ratio η = q/p' vary between -0.63 and 0.79. Thereafter, suffusion was triggered by applying a controlled multi-stage hydraulic gradient. For all tests, the first stage of applied hydraulic gradient was 0.1 and it was automatically increased by stages every 20 min (see Fig. 1(b)). In the final phase of the test, the specimen was divided into four layers to assess the spatial variability of its post-suffusion gradation. To enhance the confidence in our results, several tests were repeated and additionally, a specimen was examined under oedometric conditions, by keeping the mold around the membrane after the compaction, for comparison with data from literature.

Figures 2(a) and 2(b) show the time evolution of hydraulic conductivity and erosion rate respectively, for the tests performed with 25% of sand (FC25). Four phases can be identified: during the first phase, the hydraulic conductivity increases slightly and remains stable. The predominant process during the second phase is a self-filtration, as evidenced by the simultaneous decrease of the hydraulic conductivity and the erosion rate. The third phase is characterized by a strong increase in the erosion rate followed by an increase in hydraulic conductivity. This phase is referred to as the blow-out event and is predominantly caused by the detachment and the transport of solid particles. In the fourth phase, the hydraulic conductivity tends to stabilize while the erosion rate decreases. This can be explained by the presence of one or more preferential flow paths created by the erosion process, leading to a final steady state.

<figure><img src="image_40.png" /><figcaption>Figure 2: Time evolution of (a) hydraulic conductivity, (b) erosion rate, for FC25 specimens.</figcaption></figure>

The results for all FC25 specimens are analyzed thanks to several interpretative methods. First, by using the approach outlined by Skempton and Brogan (1994), we observe that the initiation hydraulic gradient ($i_{SB}$) for all tests is 0.2. This indicates that there is no influence of deviatoric stress variations on the critical hydraulic gradient when the mean effective stress is constant. The initiation of the blow-out can be described by the increase in hydraulic conductivity or in erosion rate. The comparison of the corresponding hydraulic gradient values shows that the hydraulic gradient associated with the turning point of erosion rate ($i_{RE}$) is consistently lower than that associated with the strong increase in hydraulic conductivity ($i_{HC}$) across all test configurations. This suggests that the use of the erosion rate time evolution permits to detect earlier the initiation of the blowout, than by using the hydraulic conductivity time evolution. From the energy-based approach, the energy corresponding to the
technique. Afterward, the sample was confined with a confining pressure of 20 kPa to prevent preferential flow between the specimen and the membrane during the saturation phase. Subsequently, CO₂ was injected at a slow rate for 20 minutes and tap water was upwardly injected at the rate of 0.25 mm/min to limit erosion prior to erosion phase. Upon completion of the saturation, the sample was consolidated to the target mean effective stress p' of 70 kPa and, afterwards, vertically loaded to a target value of deviatoric stress q. The corresponding values of shear stress ratio η = q/p' vary between -0.63 and 0.79. Thereafter, suffusion was triggered by applying a controlled multi-stage hydraulic gradient. For all tests, the first stage of applied hydraulic gradient was 0.1 and it was automatically increased by stages every 20 min (see Fig. 1(b)). In the final phase of the test, the specimen was divided into four layers to assess the spatial variability of its post-suffusion gradation. To enhance the confidence in our results, several tests were repeated and additionally, a specimen was examined under oedometric conditions, by keeping the mold around the membrane after the compaction, for comparison with data from literature.

Figures 2(a) and 2(b) show the time evolution of hydraulic conductivity and erosion rate respectively, for the tests performed with 25% of sand (FC25). Four phases can be identified: during the first phase, the hydraulic conductivity increases slightly and remains stable. The predominant process during the second phase is a self-filtration, as evidenced by the simultaneous decrease of the hydraulic conductivity and the erosion rate. The third phase is characterized by a strong increase in the erosion rate followed by an increase in hydraulic conductivity. This phase is referred to as the blow-out event and is predominantly caused by the detachment and the transport of solid particles. In the fourth phase, the hydraulic conductivity tends to stabilize while the erosion rate decreases. This can be explained by the presence of one or more preferential flow paths created by the erosion process, leading to a final steady state.

<figure><img src="image_40.png" /><figcaption>Figure 2: Time evolution of (a) hydraulic conductivity, (b) erosion rate, for FC25 specimens.</figcaption></figure>

The results for all FC25 specimens are analyzed thanks to several interpretative methods. First, by using the approach outlined by Skempton and Brogan (1994), we observe that the initiation hydraulic gradient ($i_{SB}$) for all tests is 0.2. This indicates that there is no influence of deviatoric stress variations on the critical hydraulic gradient when the mean effective stress is constant. The initiation of the blow-out can be described by the increase in hydraulic conductivity or in erosion rate. The comparison of the corresponding hydraulic gradient values shows that the hydraulic gradient associated with the turning point of erosion rate ($i_{RE}$) is consistently lower than that associated with the strong increase in hydraulic conductivity ($i_{HC}$) across all test configurations. This suggests that the use of the erosion rate time evolution permits to detect earlier the initiation of the blowout, than by using the hydraulic conductivity time evolution. From the energy-based approach, the energy corresponding to the
turning intersection of cumulative lost mass vs cumulative energy ($E_{MVE}$) is determined. The quite close values of $E_{MVE}$ and of erosion resistance index $I_α$ for all tests performed under different deviatoric stresses highlight the little influence of the deviatoric stress on these parameters.

The comparison between oedometric and triaxial stress states shows that the occurrence of a blow-out event is systematically earlier in the oedometric configuration. It was also noted that the corresponding hydraulic gradients ($i_{RE}$, $i_{HC}$) and the cumulative expended energy ($E_{MVE}$) are lower in oedometric stress conditions than in triaxial stress ones, for all tested sand %. Consequently, the oedometric configuration leads to a lower resistance. In detail, post-suffusion grain size distributions exhibit that oedometric stress conditions seem to favor a circumferential preferential flow that shortens the time required to reach the final steady state.

Now considering the results with other sand percentages, it is worth noting that the initiation hydraulic gradient, from the Skempton and Brogan (1994) method cannot always be determined. In contrast, the initiation of the blow-out based on the M versus E relationship could be measured for all tests. This parameter could serve as threshold since it can always be estimated both in the laboratory and on-site by multiplying the flow rate and difference of head.

Finally, if we consider the parameters related to the initiation of the blow-out (values of hydraulic gradient or expended energy: $i_{RE}$, $i_{HC}$, or $E_{MVE}$), the most resistant mixture seems to be the soil with 35% of fines. Whereas by considering the whole development and the following parameters: axial strain, percentage of eroded mass or $I_α$, the 15FC mixture appears as the most resistant among the tested mixtures.

## References

* Chang, D.S., and L.M. Zhang (2013). Critical hydraulic gradients of internal erosion under complex stress states. J. Geotech. Geoenviron. Eng., Vol. 139, No. 9, pp. 1454–1467

* Fell, R., and J. J. Fry. (2013). Erosion in geomechanics applied to dams and levees. Bonelli S. Editor., ISTE-Wiley., pp. 1-99.

* Garner, S.J., and R.J. Fannin (2010). Understanding internal erosion: a decade of research following a sinkhole event. The International Journal on Hydropower and Dams, 17, pp. 93–98

* Marot, D., P. L. Regazzoni, and T. Wahl (2011). Energy based method for providing soil surface erodibility rankings. J. Geotech. Geoenviron. Eng., Vol. 137, No. 12, pp. 1290-1293.

* Marot, D., A. Rochim, H. H. Nguyen, F. Bendahmane, and L. Sibille (2016). Assessing the susceptibility of gap graded soils to internal erosion characterization: proposition of a new experimental methodology. Nat. Hazards, Vol. 83, No. 1, pp. 365-388.

* Rochim, A., D. Marot, L. Sibille, and V. T. Le (2017). Effect of hydraulic loading history on suffusion susceptibility of cohesionless soils. J. Geotech. Geoenviron. Eng., Vol. 143, No. 7,

* Sibille, L., F. Lominé, P. Poullain, Y. Sail, and D. Marot (2015). Internal erosion in granular media: direct numerical simulations and energy interpretation. Hydrological Processes, Vol. 29, No. 9, pp. 2149-2163.

* Skempton, A. W., and J. M. Brogan (1994). Experiments on piping in sandy gravels. Géotechnique, Vol. 44, No. 3, pp. 440-460.
turning intersection of cumulative lost mass vs cumulative energy ($E_{MVE}$) is determined. The quite close values of $E_{MVE}$ and of erosion resistance index $I_α$ for all tests performed under different deviatoric stresses highlight the little influence of the deviatoric stress on these parameters.

The comparison between oedometric and triaxial stress states shows that the occurrence of a blow-out event is systematically earlier in the oedometric configuration. It was also noted that the corresponding hydraulic gradients ($i_{RE}$, $i_{HC}$) and the cumulative expended energy ($E_{MVE}$) are lower in oedometric stress conditions than in triaxial stress ones, for all tested sand %. Consequently, the oedometric configuration leads to a lower resistance. In detail, post-suffusion grain size distributions exhibit that oedometric stress conditions seem to favor a circumferential preferential flow that shortens the time required to reach the final steady state.

Now considering the results with other sand percentages, it is worth noting that the initiation hydraulic gradient, from the Skempton and Brogan (1994) method cannot always be determined. In contrast, the initiation of the blow-out based on the M versus E relationship could be measured for all tests. This parameter could serve as threshold since it can always be estimated both in the laboratory and on-site by multiplying the flow rate and difference of head.

Finally, if we consider the parameters related to the initiation of the blow-out (values of hydraulic gradient or expended energy: $i_{RE}$, $i_{HC}$, or $E_{MVE}$), the most resistant mixture seems to be the soil with 35% of fines. Whereas by considering the whole development and the following parameters: axial strain, percentage of eroded mass or $I_α$, the 15FC mixture appears as the most resistant among the tested mixtures.

## References

* Chang, D.S., and L.M. Zhang (2013). Critical hydraulic gradients of internal erosion under complex stress states. J. Geotech. Geoenviron. Eng., Vol. 139, No. 9, pp. 1454–1467

* Fell, R., and J. J. Fry. (2013). Erosion in geomechanics applied to dams and levees. Bonelli S. Editor., ISTE-Wiley., pp. 1-99.

* Garner, S.J., and R.J. Fannin (2010). Understanding internal erosion: a decade of research following a sinkhole event. The International Journal on Hydropower and Dams, 17, pp. 93–98

* Marot, D., P. L. Regazzoni, and T. Wahl (2011). Energy based method for providing soil surface erodibility rankings. J. Geotech. Geoenviron. Eng., Vol. 137, No. 12, pp. 1290-1293.

* Marot, D., A. Rochim, H. H. Nguyen, F. Bendahmane, and L. Sibille (2016). Assessing the susceptibility of gap graded soils to internal erosion characterization: proposition of a new experimental methodology. Nat. Hazards, Vol. 83, No. 1, pp. 365-388.

* Rochim, A., D. Marot, L. Sibille, and V. T. Le (2017). Effect of hydraulic loading history on suffusion susceptibility of cohesionless soils. J. Geotech. Geoenviron. Eng., Vol. 143, No. 7,

* Sibille, L., F. Lominé, P. Poullain, Y. Sail, and D. Marot (2015). Internal erosion in granular media: direct numerical simulations and energy interpretation. Hydrological Processes, Vol. 29, No. 9, pp. 2149-2163.

* Skempton, A. W., and J. M. Brogan (1994). Experiments on piping in sandy gravels. Géotechnique, Vol. 44, No. 3, pp. 440-460.
Role of soil preparation in controlling the erodibility behavior of cohesive soils when subjected to surface erosion

Shadi Youssef*, Nadia Benahmed*, Pierre Philippe*, Sylvie Nicaise*, Abdelkrim Bennabi**, Adrien Poupardin**

* INRAE, RECOVER Research Unit, Aix Marseille University, 3275 route de Cézanne, CS 40061, 13182 Aix-en-Provence Cedex 5, France

** Université Paris-Est, Institut de Recherche en Constructibilité, ESTP, 28 avenue du Président Wilson, 94234 Cachan, France

Shadi.youssef@inrae.fr

Nadia.benahmed@inrae.fr

Pierre.phillipe@inrae.fr

Sylvie.nicaise@inrae.fr

abennabi@estp-paris.eu

apoupardin@estp-paris.eu

**Keywords:** Internal erosion, soil preparation method, hole erosion test, erosion function apparatus, cohesive soils, erosion properties, x-ray tomography.

Abstract

Hydraulic earth structures are built with soils that are susceptible to erosion, which can have diverse impacts on the structures. Erosion is the process whereby water flow removes and carries away specific particles of soil [1]. External erosion takes place on the surface of the structure in contact with water, while internal erosion occurs within the structure. Fry et al. (2012) [2] and Bonelli (2013) [3] have identified four primary mechanisms that initiate internal erosion: suffusion, leakage concentration, regressive erosion, and contact erosion. This section focuses primarily on the characterization of concentrated leak erosion. When a leakage concentration zone exists, it promotes the detachment of particles from the soil surface in that zone. Several researchers have made efforts to devise techniques and assessments to evaluate the erodibility of soil. Some of these methods ascertain soil erosion qualitatively, while others aim to establish a correlation between the application of hydraulic shear stresses and the resistance of soil to erosion. When water flows over soil, hydraulic shear stresses occur at the interface between them, which can cause the detachment and transportation of soil particles beyond a certain threshold. The extent of erosion depends on various parameters that can be categorized into two groups: (i) hydraulic parameters linked to the potency of the erosive agents, and (ii) parameters associated with a specific soil's resistance to applied erosive forces. The correlation between the rate of soil erosion and hydraulic shear stress is deemed to be indicative of the soil's erodibility. The erosion law can therefore be expressed by the relationship given by the equation below:
Role of soil preparation in controlling the erodibility behavior of cohesive soils when subjected to surface erosion

Shadi Youssef*, Nadia Benahmed*, Pierre Philippe*, Sylvie Nicaise*, Abdelkrim Bennabi**, Adrien Poupardin**

* INRAE, RECOVER Research Unit, Aix Marseille University, 3275 route de Cézanne, CS 40061, 13182 Aix-en-Provence Cedex 5, France

** Université Paris-Est, Institut de Recherche en Constructibilité, ESTP, 28 avenue du Président Wilson, 94234 Cachan, France

Shadi.youssef@inrae.fr

Nadia.benahmed@inrae.fr

Pierre.phillipe@inrae.fr

Sylvie.nicaise@inrae.fr

abennabi@estp-paris.eu

apoupardin@estp-paris.eu

**Keywords:** Internal erosion, soil preparation method, hole erosion test, erosion function apparatus, cohesive soils, erosion properties, x-ray tomography.

Abstract

Hydraulic earth structures are built with soils that are susceptible to erosion, which can have diverse impacts on the structures. Erosion is the process whereby water flow removes and carries away specific particles of soil [1]. External erosion takes place on the surface of the structure in contact with water, while internal erosion occurs within the structure. Fry et al. (2012) [2] and Bonelli (2013) [3] have identified four primary mechanisms that initiate internal erosion: suffusion, leakage concentration, regressive erosion, and contact erosion. This section focuses primarily on the characterization of concentrated leak erosion. When a leakage concentration zone exists, it promotes the detachment of particles from the soil surface in that zone. Several researchers have made efforts to devise techniques and assessments to evaluate the erodibility of soil. Some of these methods ascertain soil erosion qualitatively, while others aim to establish a correlation between the application of hydraulic shear stresses and the resistance of soil to erosion. When water flows over soil, hydraulic shear stresses occur at the interface between them, which can cause the detachment and transportation of soil particles beyond a certain threshold. The extent of erosion depends on various parameters that can be categorized into two groups: (i) hydraulic parameters linked to the potency of the erosive agents, and (ii) parameters associated with a specific soil's resistance to applied erosive forces. The correlation between the rate of soil erosion and hydraulic shear stress is deemed to be indicative of the soil's erodibility. The erosion law can therefore be expressed by the relationship given by the equation below:
$$ \dot{\epsilon}_t = \begin{cases} C_e(\tau_t - \tau_c), & \tau_t > \tau_c \\ 0, & \tau_t < \tau_c \end{cases} $$

Where, $\dot{\epsilon}_t$ is the rate of erosion per unit surface area of the hole at time t [Kg/(s/m²)], $C_e$ is a constant named by the authors as the Coefficient of Soil Erosion [s/m] (it is then noted $k_d$ and defined such that $k_{er} = \rho_s k_d, \rho_s$ being the dry density of the soil [4]), $\tau_t$ is the hydraulic shear stress along the hole at time t [N/m²], $\tau_c$ is the Critical Shear Stress by a number of authors [N/m²].

Yin et al. (2019) have suggested that the method used for preparing soil mixtures can impact their inherent properties and internal soil structure [5]. The authors have proposed three protocols for mixing soil and have found that the first protocol results in the most uniform structure. Therefore, the current study aims to examine how soil preparation affects erosion parameters, including critical shear stress ($\tau_c$) and erosion coefficient ($K_{er}$). To achieve this objective, HET and EFA tests were conducted on compacted soil samples (five samples per mixture). All samples were comprised of 50% fine cohesive material (KA) and 50% coarse material (HN), with a water content of 15.6% at optimum level, based on the following methodology:

* Method P1: water and sand are firstly mixed, and then clay is progressively added.

* Method P2: sand and clay are firstly mixed, and then water is sprayed all over the mixture.

* Method P3: clay and water are firstly mixed, and then sand is progressively added.

Five different degrees of compaction, ranging from approximately 88% to 100% of Standard Proctor maximum dry density, were applied to the soil of each mixture, resulting in dry densities ranging from 1.50 $t/m^3$ to 1.71 $t/m^3$. The results obtained were used to complementarily discuss the effect of compaction on erosion resistance. Additionally, X-ray tomography was utilized to interpret the 3D mesoscopic structure of KA-HN mixtures. For this purpose, three soil samples, measuring 20 mm in height and 50 mm in diameter, were dynamically compacted at the same density ($\rho_d = 1.623 t/m^3$), which is 95% of the maximum dry density, with a water content corresponding to the optimum ($W - opt = 15.6\%$). These samples were composed of 50% Armoricaine kaolinite (KA) and 50% Hostun sand (HN0.4/0.8), and were prepared using three different methods, namely Method P1, Method P2, and Method P3. The corresponding samples were labeled as S1, S2, and S3, respectively.

The erosion results (Figure 1 and Figure 2) show that, the samples produced by preparation method P3 were the most resistant to critical shear stress, but paradoxically, they also eroded the most rapidly with the highest $K_{er}$ values, especially during EFA tests. On the other hand, the samples prepared using method P2 had the lowest critical shear stress values and systematically smaller erosion coefficient values than those produced by method P3. Samples produced by method P1 appeared to be the least erodible, with $K_{er}$ values often similar to those obtained with method P2, and their resistance to erosion fell between the two previous cases. In terms of microstructural properties (Figure 3), the sample prepared by method P1 had the most heterogeneous particle distribution with aggregates of sand grains separated by clods of pure clay and very little interstitial void. The sample prepared by method P2 had a relatively homogeneous distribution of sand grains with some small empty spots, and the sample prepared by method P3 also had strong heterogeneities, perhaps slightly less than for method P1, with the presence of a little more interstitial void.

## References

[1] S. Bonelli, *Erosion of Geomaterials*. Wiley, 2012.

[2] J. Fry, "Introduction to the Process of Internal Erosion in Hydraulic Structures: Embankment Dams and Dikes," in *Erosion of Geomaterials*, Wiley, 2012, pp. 1–37.
$$ \dot{\epsilon}_t = \begin{cases} C_e(\tau_t - \tau_c), & \tau_t > \tau_c \\ 0, & \tau_t < \tau_c \end{cases} $$

Where, $\dot{\epsilon}_t$ is the rate of erosion per unit surface area of the hole at time t [Kg/(s/m²)], $C_e$ is a constant named by the authors as the Coefficient of Soil Erosion [s/m] (it is then noted $k_d$ and defined such that $k_{er} = \rho_s k_d$, $\rho_s$ being the dry density of the soil [4]), $\tau_t$ is the hydraulic shear stress along the hole at time t [N/m²], $\tau_c$ is the Critical Shear Stress by a number of authors [N/m²].

Yin et al. (2019) have suggested that the method used for preparing soil mixtures can impact their inherent properties and internal soil structure [5]. The authors have proposed three protocols for mixing soil and have found that the first protocol results in the most uniform structure. Therefore, the current study aims to examine how soil preparation affects erosion parameters, including critical shear stress ($\tau_c$) and erosion coefficient ($K_{er}$). To achieve this objective, HET and EFA tests were conducted on compacted soil samples (five samples per mixture). All samples were comprised of 50% fine cohesive material (KA) and 50% coarse material (HN), with a water content of 15.6% at optimum level, based on the following methodology:

* Method P1: water and sand are firstly mixed, and then clay is progressively added.

* Method P2: sand and clay are firstly mixed, and then water is sprayed all over the mixture.

* Method P3: clay and water are firstly mixed, and then sand is progressively added.

Five different degrees of compaction, ranging from approximately 88% to 100% of Standard Proctor maximum dry density, were applied to the soil of each mixture, resulting in dry densities ranging from 1.50 $t/m^3$ to 1.71 $t/m^3$. The results obtained were used to complementarily discuss the effect of compaction on erosion resistance. Additionally, X-ray tomography was utilized to interpret the 3D mesoscopic structure of KA-HN mixtures. For this purpose, three soil samples, measuring 20 mm in height and 50 mm in diameter, were dynamically compacted at the same density ($\rho_d = 1.623 t/m^3$), which is 95% of the maximum dry density, with a water content corresponding to the optimum ($W - opt = 15.6\%$). These samples were composed of 50% Armoricaine kaolinite (KA) and 50% Hostun sand (HN0.4/0.8), and were prepared using three different methods, namely Method P1, Method P2, and Method P3. The corresponding samples were labeled as S1, S2, and S3, respectively.

The erosion results (Figure 1 and Figure 2) show that, the samples produced by preparation method P3 were the most resistant to critical shear stress, but paradoxically, they also eroded the most rapidly with the highest $K_{er}$ values, especially during EFA tests. On the other hand, the samples prepared using method P2 had the lowest critical shear stress values and systematically smaller erosion coefficient values than those produced by method P3. Samples produced by method P1 appeared to be the least erodible, with $K_{er}$ values often similar to those obtained with method P2, and their resistance to erosion fell between the two previous cases. In terms of microstructural properties (Figure 3), the sample prepared by method P1 had the most heterogeneous particle distribution with aggregates of sand grains separated by clods of pure clay and very little interstitial void. The sample prepared by method P2 had a relatively homogeneous distribution of sand grains with some small empty spots, and the sample prepared by method P3 also had strong heterogeneities, perhaps slightly less than for method P1, with the presence of a little more interstitial void.

## References

[1] S. Bonelli, Erosion of Geomaterials. Wiley, 2012.

[2] J. Fry, "Introduction to the Process of Internal Erosion in Hydraulic Structures: Embankment Dams and Dikes," in Erosion of Geomaterials, Wiley, 2012, pp. 1–37.
[3] S. Bonelli, *Erosion in Geomechanics Applied to Dams and Levees*. Hoboken, NJ, USA: John Wiley & Sons, Inc., 2013.

[4] N. Benahmed and S. Bonelli, "Investigating concentrated leak erosion behaviour of cohesive soils by performing hole erosion tests," *European Journal of Environmental and Civil Engineering*, vol. 16, no. 1, pp. 43–58, 2012.

[5] K. Yin et al., "Influence of sample preparation on the multi scale structure of sand-clay mixtures," *E3S Web of Conferences*, vol. 92, p. 01007, 2019.

## Figures

<figure><img src="image_41.png" /><figcaption>**Figure 1:** Influence of the soil preparation method on the erosion parameters, i.e. critical shear stress $\tau_c$ and erosion coefficient $K_{er}$, for KA50HN50 mixtures with different compaction degrees, using HET.</figcaption></figure>

<figure><img src="image_42.png" /><figcaption>**Figure 2:** Influence of soil preparation method on the erosion parameters, i.e. critical shear stress $\tau_c$ and erosion coefficient $K_{er}$, for KA50HN50 mixtures with different compaction degrees, using EFA.</figcaption></figure>
[3] S. Bonelli, *Erosion in Geomechanics Applied to Dams and Levees*. Hoboken, NJ, USA: John Wiley & Sons, Inc., 2013.

[4] N. Benahmed and S. Bonelli, "Investigating concentrated leak erosion behaviour of cohesive soils by performing hole erosion tests," *European Journal of Environmental and Civil Engineering*, vol. 16, no. 1, pp. 43–58, 2012.

[5] K. Yin et al., "Influence of sample preparation on the multi scale structure of sand-clay mixtures," *E3S Web of Conferences*, vol. 92, p. 01007, 2019.

## Figures

<figure><img src="image_41.png" /><figcaption>**Figure 1:** Influence of the soil preparation method on the erosion parameters, i.e. critical shear stress $\tau_c$ and erosion coefficient $K_{er}$, for KA50HN50 mixtures with different compaction degrees, using HET.</figcaption></figure>

<figure><img src="image_42.png" /><figcaption>**Figure 2:** Influence of soil preparation method on the erosion parameters, i.e. critical shear stress $\tau_c$ and erosion coefficient $K_{er}$, for KA50HN50 mixtures with different compaction degrees, using EFA.</figcaption></figure>
Method P1

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

Method P2

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

Method P3

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

**Figure 3:** X-ray tomography scans from the KA50HN50 compacted soil samples prepared following the three different methods (P1, P2 and P3). Left: raw tomography scans; Middle: images after threshold 1 to isolate the voids; Right: images after threshold 2 to isolate the grains.
Method P1

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

Method P2

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

Method P3

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

**Figure 3:** X-ray tomography scans from the KA50HN50 compacted soil samples prepared following the three different methods (P1, P2 and P3). Left: raw tomography scans; Middle: images after threshold 1 to isolate the voids; Right: images after threshold 2 to isolate the grains.
# Experimental study of the hydromechanical instabilities of a cemented granular material

Abbas FARHAT<fnref n="1" />,<fnref n="2" />, Pierre PHILIPPE<fnref n="1" />, Li Hua LUU<fnref n="1" />, Pablo CUELLAR<fnref n="3" />, Nadia BENHAMED<fnref n="1" />, and Torsten WICHTMANN<fnref n="2" />

**Keywords:** Cemented granular material, mechanical characterisation, fluidisation and hydromechanical instabilities.

## Abstract

Civil engineering hydraulic earthen structures, such as dams, levees and embankments, are built for various functions such as retaining water, possibly to produce electricity, protecting against flooding, or canalising a river. Such hydraulic structures are susceptible to erosion of soil particles, which can be caused by either seepage flow or overtopping flow, potentially resulting in hydromechanical instabilities. Erosion is by far the most common cause for failure of hydraulic earthen structures. Two different situations can be distinguished: External erosion, the degradation of the outer surface of the hydraulic structure in the event of overflow or over-topping, and internal erosion, which takes place inside the structure or the foundation by seepage flow. When infiltrating through the soil or rock foundation of the structure, water can create a pressure gradient capable to erode soil particles and, eventually, lead the structure to fail. This phenomenon is exacerbated by the presence of fractures or discontinuities in the foundation, which can increase the rate of seepage and result in further erosion.

The hydromechanical instability of hydraulic structures due to seepage flow is a complex phenomenon that can be influenced by a range of factors, including the hydraulic gradient, soil properties, geometry of the structure, and loading conditions. The understanding of the intrinsic physics of the mechanism underlying the specific processes of soil erosion has been extensively studied at the sample scale for the specific case of granular materials Brunier-Coulin et al. (2017), Mena et al. (2017, 2018), Philippe and Badiane (2013), but very few for coherent soilsCui (2013), Cui et al. (2012, 2014), for which many questions are still open.

This work deals with water flow-induced erosion of geomaterials more complex than purely frictional soil. More specifically, we focus here on cemented granular materials made of grains connected by solid inter-particle bonds. We aim to study the mechanical characterization of an artificial cemented granular system at different scales and the development of hydro-mechanical instabilities within this model soil by a localized hydraulic loading exerted by an upward flow. This multidisciplinary research work involves different scientific communities, notably those of soil mechanics, geotechnics, fluid mechanics, and geophysics.

Our homemade artificial cemented granular material is prepared from a mixture of paraffin and spherical glass beads Farhat et al.. The material's cementation strength is quantified utilizing laboratory tests designed and carried out for investigation

<footnote n="1">¹RECOVER, INRAE, Aix Marseille Univ, Aix-en-Provence, France</footnote>

<footnote n="2">²Ruhr University Bochum Univ, Bochum, Germany</footnote>

<footnote n="3">³BAM Federal Institute for Materials Research and Testing, Berlin, Germany</footnote>
Experimental study of the hydromechanical instabilities of
a cemented granular material

Abbas FARHAT<fnref n="1" />,<fnref n="2" />, Pierre PHILIPPE<fnref n="1" />, Li Hua LUU<fnref n="1" />, Pablo CUELLAR<fnref n="3" />, Nadia BENHAMED<fnref n="1" />, and Torsten WICHTMANN<fnref n="2" />

***Keywords:*** Cemented granular material, mechanical characterisation, fluidisation and hy-
dromechanical instabilities.

Abstract

Civil engineering hydraulic earthen structures, such as dams, levees and embank-
ments, are built for various functions such as retaining water, possibly to produce
electricity, protecting against flooding, or canalising a river. Such hydraulic structures
are susceptible to erosion of soil particles, which can be caused by either seepage flow
or overtopping flow, potentially resulting in hydromechanical instabilities. Erosion
is by far the most common cause for failure of hydraulic earthen structures. Two
different situations can be distinguished: External erosion, the degradation of the outer
surface of the hydraulic structure in the event of overflow or over-topping, and internal
erosion, which takes place inside the structure or the foundation by seepage flow.
When infiltrating through the soil or rock foundation of the structure, water can create
a pressure gradient capable to erode soil particles and, eventually, lead the structure to
fail. This phenomenon is exacerbated by the presence of fractures or discontinuities in
the foundation, which can increase the rate of seepage and result in further erosion.

The hydromechanical instability of hydraulic structures due to seepage flow is a complex phenomenon that can be influenced by a range of factors, including the hydraulic gradient, soil properties, geometry of the structure, and loading conditions. The understanding of the intrinsic physics of the mechanism underlying the specific processes of soil erosion has been extensively studied at the sample scale for the specific case of granular materials Brunier-Coulin et al. (2017), Mena et al. (2017, 2018), Philippe and Badiane (2013), but very few for coherent soilsCui (2013), Cui et al. (2012, 2014), for which many questions are still open.

This work deals with water flow-induced erosion of geomaterials more complex
than purely frictional soil. More specifically, we focus here on cemented granular
materials made of grains connected by solid inter-particle bonds. We aim to study
the mechanical characterization of an artificial cemented granular system at different
scales and the development of hydro-mechanical instabilities within this model soil
by a localized hydraulic loading exerted by an upward flow. This multidisciplinary
research work involves different scientific communities, notably those of soil mechanics,
geotechnics, fluid mechanics, and geophysics.

Our homemade artificial cemented granular material is prepared from a mixture of
paraffin and spherical glass beads Farhat et al.. The material's cementation strength
is quantified utilizing laboratory tests designed and carried out for investigation

<footnote n="1">¹RECOVER, INRAE, Aix Marseille Univ, Aix-en-Provence, France</footnote>

<footnote n="2">²Ruhr University Bochum Univ, Bochum, Germany</footnote>

<footnote n="3">³BAM Federal Institute for Materials Research and Testing, Berlin, Germany</footnote>
at different scales. We determined the tensile yield strength at the macroscale of a cemented sample as well as at the interparticle microscale, between two bonded particles extracted from the preceding macroscale sample. We also investigated the yield micro-strength for solid bonds ruptured by shear, bending, and torsion. Despite the huge dispersion of the results, we discovered a satisfactory one-to-one agreement between the tensile force and the other loadings in the experiments. Furthermore, we provide a theoretical model for the micro tensile force required to rupture a bond, based on the assumption of uniform adhesive strength, where a good agreement with the experimental results is reached.

We carried out several experiments on submerged cemented granular materials subjected to a localized upward flow from bottom inlet. After a static regime where no motion is visible in the sample, three distinct modes of hydromechanical destabilization of cemented granular material were observed (see Figure 1): (i) a block rupture, characterized by the appearance of a median crack above the flow inlet; (ii) a fluidized path rupture, by progressive burrowing along the walls towards the inlet; (iii) a block uplift rupture, when the sample slides upward at the lateral walls after prior bonds breakage at its boundaries. Globally, all of these experiments showed that cemented beds made of larger beads and higher paraffin contents required a greater flow rate to be destabilized. Finally, in consistency with previous data from experiments on purely granular localized fluidization Mena et al. (2017), it was possible to satisfactorily gather all the data, approximately on Ergun's law, when considering the relevant dimensionless numbers, namely the inlet particle Reynolds number and the Archimedes number, the latter being adapted to the cemented situation by considering the micro-tensile force instead of the usual grain buoyant weight.

REFERENCES

Brunier-Coulin, Cuéllar, and Philippe. Erosion onset of a cohesionless granular medium by an immersed impinging round jet. *Phys. Rev. Fluids*, 2:034302, 2017.

Cui. *Numerical simulation of internal fluidisation and cavity evolution due to a leaking pipe using the coupled DEM-LBM technique*. PhD thesis, University of Birmingham, 2013.

Cui, Li, Chan, and Chapman. A 2d dem-lbm study on soil behaviour due to locally injected fluid. *Particuology*, 10(2):242-252, 2012.

Cui, Li, Chan, and Chapman. Coupled dem-lbm simulation of internal fluidisation induced by a leaking pipe. *Powder Technology*, 254:299-306, 2014.

Farhat, Luu, Philippe, and Cuéllar. Multi-scale cohesion force measurements for cemented granular materials.

Mena, Luu, Cuéllar, Philippe, and Curtis. Parameters affecting the localized fluidization in a particle medium. *AIChE Journal*, 63(5):1529-1542, 2017.

Mena, Brunier, Curtis, and Philippe. Experimental observation of the initial stages of localized fluidization. In *2018 AIChE Annual Meeting*. AIChE, 2018.

Philippe and Badiane. Localized fluidization in a granular medium. *Physical Review E*, 87(4): 042206, 2013.
at different scales. We determined the tensile yield strength at the macroscale of a cemented sample as well as at the interparticle microscale, between two bonded particles extracted from the preceding macroscale sample. We also investigated the yield micro-strength for solid bonds ruptured by shear, bending, and torsion. Despite the huge dispersion of the results, we discovered a satisfactory one-to-one agreement between the tensile force and the other loadings in the experiments. Furthermore, we provide a theoretical model for the micro tensile force required to rupture a bond, based on the assumption of uniform adhesive strength, where a good agreement with the experimental results is reached.

We carried out several experiments on submerged cemented granular materials subjected to a localized upward flow from bottom inlet. After a static regime where no motion is visible in the sample, three distinct modes of hydromechanical destabilization of cemented granular material were observed (see Figure 1): (i) a block rupture, characterized by the appearance of a median crack above the flow inlet; (ii) a fluidized path rupture, by progressive burrowing along the walls towards the inlet; (iii) a block uplift rupture, when the sample slides upward at the lateral walls after prior bonds breakage at its boundaries. Globally, all of these experiments showed that cemented beds made of larger beads and higher paraffin contents required a greater flow rate to be destabilized. Finally, in consistency with previous data from experiments on purely granular localized fluidization Mena et al. (2017), it was possible to satisfactorily gather all the data, approximately on Ergun's law, when considering the relevant dimensionless numbers, namely the inlet particle Reynolds number and the Archimedes number, the latter being adapted to the cemented situation by considering the micro-tensile force instead of the usual grain buoyant weight.

## REFERENCES

* Brunier-Coulin, Cuéllar, and Philippe. Erosion onset of a cohesionless granular medium by an immersed impinging round jet. *Phys. Rev. Fluids*, 2:034302, 2017.

* Cui. *Numerical simulation of internal fluidisation and cavity evolution due to a leaking pipe using the coupled DEM-LBM technique*. PhD thesis, University of Birmingham, 2013.

* Cui, Li, Chan, and Chapman. A 2d dem-lbm study on soil behaviour due to locally injected fluid. *Particuology*, 10(2):242-252, 2012.

* Cui, Li, Chan, and Chapman. Coupled dem-lbm simulation of internal fluidisation induced by a leaking pipe. *Powder Technology*, 254:299-306, 2014.

* Farhat, Luu, Philippe, and Cuéllar. Multi-scale cohesion force measurements for cemented granular materials.

* Mena, Luu, Cuéllar, Philippe, and Curtis. Parameters affecting the localized fluidization in a particle medium. *AIChE Journal*, 63(5):1529-1542, 2017.

* Mena, Brunier, Curtis, and Philippe. Experimental observation of the initial stages of localized fluidization. In *2018 AIChE Annual Meeting*. AIChE, 2018.

* Philippe and Badiane. Localized fluidization in a granular medium. *Physical Review E*, 87(4): 042206, 2013.
<figure><img src="image_46.png" /><figcaption>**Figure 1:** Four distinct regimes observed for a cemented granular material under localized hydromechanical load: (a) Static regime, (b) Fluidized path regime, (c) Block rupture regime and (d) Block uplift regime.</figcaption></figure>
<figure><img src="image_46.png" /><figcaption>**Figure 1:** Four distinct regimes observed for a cemented granular material under localized hydromechanical load: (a) Static regime, (b) Fluidized path regime, (c) Block rupture regime and (d) Block uplift regime.</figcaption></figure>