# Analyse numérique de la fragmentation par la méthode des éléments discrets dans un système multi-contacts

Mohamed Guessasma, Patrice Coorevits, Jérôme Fortin, Emmanuel Bellenger

► To cite this version:

Mohamed Guessasma, Patrice Coorevits, Jérôme Fortin, Emmanuel Bellenger. Analyse numérique de la fragmentation par la méthode des éléments discrets dans un système multi-contacts. 9e Colloque national en calcul des structures, CSMA, May 2009, Giens, France. hal-01413840

HAL Id: hal-01413840

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

Submitted on 11 Dec 2016

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

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

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

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Analyse numérique de la fragmentation par la MED dans un système multi-contacts

M. Guessasma¹, P. Coorevits¹, J. Fortin², E. Bellenger¹

¹ Laboratoire des Technologies Innovantes - EA 3899 (IUT de l'Aisne)
48 rue d'Ostende 02100 Saint-Quentin
mohamed.guessasma@u-picardie.fr, patrice.coorevits@u-picardie.fr,
emmanuel.bellenger@u-picardie.fr

² Laboratoire des Technologies Innovantes - EA 3899 (INSSET)
48 rue Raspail 02100 Saint-Quentin
jerome.fortin@u-picardie.fr

**Résumé** — On propose dans ce papier une modélisation discrète de la fragmentation dans un système multi-contacts parfaitement rigide. La formulation du contact unilatéral à frottement sec est basée sur la méthode du bipotentiel. Par ailleurs, on a également introduit un modèle de cohésion pour l'interface de contact dans le cas d'un milieu cohésif, en utilisant un critère de résistance fonction du type de sollicitation. Le modèle de cohésion proposé est implémenté dans le code de calcul par éléments discrets MULTICOR.

**Mots clés** — DEM, Mécanique du contact, Mécanique de la rupture, Fragmentation, Cohésion.

## 1 Introduction

Face aux besoins d'améliorations de la sécurité et de la fiabilité des équipements industriels, la maîtrise du frottement est devenue aujourd'hui un enjeu incontournable dans l'industrie et les transports. Plusieurs études ont révélé la complexité des phénomènes de frottement et la nécessité de raisonner en termes de couplage multi-physiques et multi-échelles. Dans ce contexte, on propose dans ce papier de développer une modélisation originale de l'interface de contact intégrant un modèle de fragmentation en utilisant la Méthode des Eléments Discrets (MED). Cette dernière est du type dynamique du contact [1] et elle est basée sur la méthode du bipotentiel pour la formulation du contact unilatéral à frottement sec [2]. Dans un premier temps, nous nous sommes intéressés à la modélisation de la fragmentation (sous sollicitation de compression) d'un milieu discret parfaitement rigide et non cohésif (Figure 1) en s'appuyant sur des travaux antérieurs [3, 4]. Les résultats numériques issus de cette modélisation ont permis de mieux comprendre les phénomènes mis en jeu lors de la fragmentation et présentent un intérêt pour certaines applications industrielles (ex. le broyage). Dans un second temps, nous avons considéré un milieu discret parfaitement rigide en tenant compte cette fois-ci du caractère cohésif. Par ce biais, il sera possible de mieux appréhender les mécanismes de fragmentation par décohésion dans un système multi-contacts. En effet, la prise en compte de l'aspect cohésif, pourrait apporter des réponses pour les systèmes avec contact frottant dont les surfaces de contact subissent une dégradation par décohésion de fines particules.
## 2 Approche discrète

La modélisation bidimensionnelle (2D) du système multi-contacts est basée sur la MED. Les particules sont modélisées par des disques parfaitement rigides et non-pénétrables. En choisissant comme paramètre de configuration q pour chaque disque (Figure 1), l'équation à modéliser s'écrit sous la forme suivante (1) :

$$ M\ddot{q} = F_{ext}(q, \dot{q}, t) + R^{\alpha}, \qquad (1) $$

où $F_{ext}$ représentent les efforts extérieurs explicitement connus et $R^{\alpha}$ les efforts intérieurs associés aux réactions de contact à priori inconnues.

<figure><img src="image_2.png" /><figcaption>Figure 1 – Système mutli-contacts</figcaption></figure>

A chaque paire de grains $\alpha$, constituée des candidats au contact $\Omega_i$ et $\Omega_j$, on lui associe un repère local dont les axes sont définis par la normale $n$ et la tangente $t$ (Figure 1). La normale $n$ est orientée de $\Omega_i$ vers $\Omega_j$. Les variables mises en dualité sont la vitesse relative $\dot{u}$ et la réaction de contact $r$ (de $\Omega_i$ par rapport à $\Omega_j$). Ces variables se décomposent dans un repère local de la façon suivante (2) :

$$ \dot{u} = \dot{u}_t + \dot{u}_n \cdot n, \quad r = r_t + r_n \cdot n \qquad (2) $$

où $\dot{u}_t$ représente la vitesse de glissement, $\dot{u}_n$ la vitesse de séparation normale, $r_t$ la force d'adhérence et $r_n$ la réaction normale de contact. L'introduction d'un frottement sec de Coulomb conduit à un problème de complémentarité non linéaire dont la résolution est basée sur la méthode du bipotentiel [2]. L'algorithme utilisé est celui d'Usawa avec un schéma de prédiction-correction donné par l'équation (3) :

$$ \text{prediction : } \tau^{ij} = r^{ij} - \gamma[\dot{u}_t^{ij} + (\dot{u}_n^{ij} + \mu\| - \dot{u}_t^{ij}\| \cdot n)], \quad \text{correction : } r^{ij} = \text{proj}(K_{\mu}, \tau^{ij}) \qquad (3) $$

où $\gamma$ est un paramètre numérique et $K_{\mu}$ est le cône de Coulomb. La projection de la prédiction $\tau^{ij}$ sur le cône de Coulomb conduit à l'un des états suivants : relâchement, contact avec adhérence ou glissement.

## 3 Champ de contraintes et mécanisme de rutture

Pour comprendre les mécanismes qui gouvernent la fragmentation dans un système multi-contacts sous sollicitation de compression, il convient de déterminer dans un premier temps le tenseur des contraintes dans chaque particule. La connaissance du champ de contraintes permettra ainsi de prédire la fragmentation de la particule par rapport à un critère de rupture donné. Considérons une particule $p$ sous l'action de plusieurs forces de contact réparties aléatoirement (Figure 2(a)), le tenseur des contraintes qui en découle est calculé de la manière suivante (4) :

$$ \sigma_p = \frac{a_p}{V_p} \sum_{\alpha} n_{c\alpha} \otimes r_{\alpha} \qquad (4) $$
avec $a_p$ et $V_p$ respectivement le rayon et le volume de la particule $p$, $r_\alpha$ les forces de contact appliquées à la particule, $n_{c\alpha}$ le vecteur unitaire normal et $\alpha$ le nombre de contacts (Figure 2(a)).

<figure><img src="image_3.png" /><figcaption>Figure 2 – Répartition aléatoire des forces de contact (a) - Compression biaxiale (b)</figcaption></figure>

## 3.1 Mécanisme de rupture sous sollicitation de compression biaxiale

Connaissant le champ de contrainte $\sigma_p$, on peut facilement calculer les contraintes majeure et mineure respectivement $\sigma_{max}$ et $\sigma_{min}$ selon les directions principales. Cela revient à considérer la particule $p$ sous un chargement de compression biaxiale (Figure 2(b)). Les efforts de compression $r_{max}$ et $r_{min}$ sont calculés à partir des expressions suivantes (5) :

$$r_{max,min} = \pi \frac{a_p}{2} \sigma_{max,min} \qquad (5)$$

En s'appuyant sur les techniques de calcul lors d'un essai brésilien, les contraintes de traction $\sigma_{traction}$ et de compression $\sigma_{compression}$, agissant au voisinage du centre de la particule, sont calculées en supposant les effets des contraintes majeure $\sigma_{max}$ et mineure $\sigma_{min}$. Ces contraintes se calculent comme suit (6) :

$$\sigma_{traction} = \frac{\sigma_{max} - 3\sigma_{min}}{2} , \quad \sigma_{compression} = \frac{\sigma_{min} - 3\sigma_{max}}{2} \qquad (6)$$

A partir des constatations observées lors des essais de compression sur les matériaux fragiles, il s'en suit que $\sigma_{traction}$ est à l'origine de la fragmentation. Cette contrainte donne lieu dans un premier temps à une amorce de fissure qui se propage par la suite le long de la direction principale majeure conduisant ainsi à la fragmentation. Différentes études ont été menées afin de définir la contrainte critique à partir de laquelle la rupture se produit. Ces études ont conclu que cette valeur seuil dépend de plusieurs facteur, notamment, la géométrie, la dimension du problème et la contrainte [3, 4].

## 3.2 Contrainte critique en rupture et modèle numérique de fragmentation

Il existe plusieurs formulations pour la détermination de la contrainte critique à la rupture en fonction de la prise en compte ou pas de l'existence de fissures initiales au sein de la particule. Nous avons utilisé pour le calcul de cette contrainte seuil une formulation basée sur la théorie de Weibull. L'expression de la contrainte critique $\sigma_{crit}$ est donnée par l'équation suivante (7) :

$$\sigma_{crit} = \beta K_0 a_p^{\frac{d}{m}} , \quad \beta = \frac{2^\alpha}{\pi} , \quad \alpha = d - 1 - \frac{d}{m} \qquad (7)$$

avec $m$ le module de Weibull, $d$ la dimension du problème et $K_0$ le préfacteur. Quand $m \to \infty$, $K_0$ est équivalent à une contrainte. Dans le cas particulier où $m = 2$, $K_0$ est équivalent à un facteur d'intensité de contrainte $K_{IC}$ en mode I (mode en traction). Quand la condition $\sigma_{traction} > \sigma_{crit}$ est satisfaite, il y aura donc fragmentation de la particule par traction.
Plusieurs modèles numériques sont proposés dans la littérature afin de décrire le mode de fragmentation [3, 5, 6]. Dans le cadre de ce travail, on a utilisé le modèle proposé par Tsoungui et al [3], qui considère la particule initiale comme un assemblage de 12 particules (Figure 3). De part et d'autre des directions majeure et mineure on a respectivement 2 particules de rayon $a_p/2$ et 10 particules de rayons $a_p/3$, $a_p/6$ et $a_p/11$.

<figure><img src="image_4.png" /><figcaption>Figure 3 – Modèle numérique de fragmenation proposé par Tsoungui et al [3]</figcaption></figure>

Le processus de fragmentation est arrêté pour une taille de particule critique, $d_{crit}$, donnée. L'expression qui permet de calculer $d_{crit}$ est issue des travaux de Kendall [7] et s'écrit comme suit : $d_{crit} = \frac{32}{2} (\frac{K_{IC}}{\sigma_e})^2$, avec $K_{IC}$ le facteur d'intensité de contrainte et $\sigma_e$ la contrainte élastique limite.

# 4 Simulation numérique d'un essai œdométrique

Le modèle numérique de fragmentation est mis en œuvre à travers un essai œdométrique réalisé avec le code de calcul MULTICOR [8]. On a considéré un empilement de 211 particules de rayons 1.25, 0.75 et 0.5mm distribués uniformément. Le matériau choisi est le plâtre avec une masse volumique de 1200kg/m³, une résistance en traction de 10MPa, un coefficient de friction de 0.2 et un module de Weibull de 11 (Figure 4).

<figure><img src="image_5.png" /><figcaption>Figure 4 – Essai œdométrique</figcaption></figure>

L'empilement a une taille de $2.5 \times 2.7 \text{cm}^2$. Des conditions périodiques sont imposées sur les faces latérales. La paroi du haut a une vitesse constante de $-0.1 \text{m/s}$.
La Figure 5 est une représentation des contraintes de traction et de compression moyennées par couche en fonction de la hauteur de l'empilement. On remarque sur les deux courbes le niveau élevé des contraintes dans la zone proche de la paroi fixe. Cette différence pourrait être expliquée par le fait que les particules se trouvant dans la partie inférieure supporteraient un chargement hydrostatique. La Figure 6 permet de suivre l'évolution du nombre de grains au cours de l'essai œdométrique.
<figure><img src="image_6.png" /><figcaption>Figure 5 – Variations des contraintes de traction et de compression en fonction de l'épaisseur</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Figure 6 – Evolution du nombre de grains</figcaption></figure>

## 5 Introduction d'un modèle de cohésion

Le modèle de cohésion utilisé est basé sur une loi de contact contraignant le mouvement relatif des particules. Ce modèle doit vérifier à la fois les critères de résistance en traction, au cisaillement et au roulement [9]. La réaction de contact entre deux particules (2), donnée précédemment, s'écrit donc comme suit (8) :

$$r = r_t + (r_n + r_{cohe}).n \quad (8)$$

où $r_{cohe}$ est la force de cohésion calculée par rapport à une contrainte limite en traction $\sigma_{cohe}$ (carréristique intrinsèque du matériau) : $r_{cohe} = \sigma_{cohe} \cdot a_{moy}$ ; $a_{moy}$ étant le rayon moyen des particules en contact. La décohésion survient au point de contact si l'une des deux conditions suivantes est satisfaite : $r_n + r_{cohe} < 0$ ou $|r_t| > \mu \cdot (r_n + r_{cohe})$ ; $r_t$ étant la composante tangentielle générée par cisaillement ou roulement.

Pour tester le modèle de cohésion, on a réalisé un essai de compression uniaxiale sur une particule élémentaire représentée par un agglomérat de petites particules faiblement cohésives (Figure 7(a)). Les Figures 7(b) et 7(c) montrent la dégradation de l'agglomérat sous l'effet de l'effort de compression par la formation de plusieurs lignes de glissement, caractérisant ainsi le mode de rupture de matériaux faiblement rigides (ex. comportement lâche de certains géomatériaux).
<figure><img src="image_8.png" /><figcaption>Figure 7 – Agglomérat de particules cohésives en compression uniaxiale</figcaption></figure>

# 6 Conclusion

A travers cette étude, nous avons implémenté dans le code de calcul MULTICOR un modèle de fragmentation dans le cas d'un système multi-contacts non cohésif et parfaitement rigide. Sous un chargement œdométrique, le mécanisme de fragmentation, à l'échelle de la particule, est provoqué par la contrainte de traction conséquence de la superposition des efforts de contact. La fragmentation a lieu pour une valeur de contrainte seuil calculée d'après une formulation statistique. La prise en compte de la cohésion vient compléter notre travail sur la modélisation de la rupture pour un système multi-contacts aggloméré. L'exemple du test de compression unixaxial réalisé sur un agglomérat de particules cohésives permet de valider les conditions de contact en présence de forces adhésives. D'autres test mécaniques avec plusieurs agglomérats devront nous aider à mieux comprendre la complexité du comportement d'un milieu discret cohésif.

# Références

[1] M. Jean. The Non Smooth Contact Dynamics Method. *Compt. Methods Appl. Math. Engrg.*, 177, 235-257, 1999.

[2] G. De Saxcé, Z.Q. Feng. The bipotential method : a constructive approach to design the complete contact law with friction and improved numerical algorithms. *Int. J. Math. Comput. Model.*, 28, 225-245, 1998.

[3] O. Tsoungui, D. Vallet, J.C. Charmet. Numerical model of crushing of grains inside two-dimensional branular materials. *Powder Technology*, 105, 190-198, 1999.

[4] A.W. Larson. *The effect of size and contact geometry in single particle crushing*, Master's thesis, MIT, Boston, 1986.

[5] J.A. Aström, H.J. Herrmann. Fragmentation of grains in a two-dimensional packing. *Eur. Phys. J.*, 5, 551-554, 1998.

[6] S. Lobo-Guerrero, L.E. Vallejo. Influence of pile shape and pile interaction on the crushable behaviour of granular material around driven piles : DEM analyses. *Granular Matter*, 9, 241-250, 2007.

[7] K. Kendall. The impossibility of comminuting small particles by compression. *Nature*, 272, 710-711, 1978.

[8] J. Fortin, O. Millet, G. de Saxce, Numerical simulation of granular materials by an improved discrete element method. *Int. J. For Num. Meth. in Ing.*, 62, 639-663, 2004.

[9] A. Taboada, N. Estrada, F. Radjaï. Additive decomposition of shear strength in cohesive granular media from grain-scale intraction. *Phys. Rev. Lett.*, 97, 2006.