# A Numerical Model for the Elasto-Plastic Behavior of Cellular and Granular Materials

Yohann Trivino, Vincent Richefeu, Farhang Radjai, Komlanvi Lampoh, Jean-Yves Delenne

► To cite this version:

Yohann Trivino, Vincent Richefeu, Farhang Radjai, Komlanvi Lampoh, Jean-Yves Delenne. A Numerical Model for the Elasto-Plastic Behavior of Cellular and Granular Materials. IOP Conference Series: Earth and Environmental Science, 2025, 1480, pp.012037. 10.1088/1755-1315/1480/1/012037. hal-05493359

HAL Id: hal-05493359

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

Submitted on 4 Feb 2026

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

PAPER • OPEN ACCESS

# A Numerical Model for the Elasto-Plastic Behavior of Cellular and Granular Materials

To cite this article: Y. Trivino et al 2025 IOP Conf. Ser.: Earth Environ. Sci. 1480 012037

View the [article online](#) for updates and enhancements.

You may also like

- Atom Fountains at NTSC
Jun Ruan, Xinliang Wang, Hui Zhang et al.

- Effect of thermal deformation on the microstructure and tensile property of the TiZrAl medium entropy alloy
Rong Guo, Jiaxiang Cui, Bang Dou et al.

- Development and performance research of flash welding process for DH32 marine bulb flat steel
Su Li, Shiyi Gao, Zhenlong Li et al.
# A NUMERICAL MODEL FOR THE ELASTO-PLASTIC BEHAVIOR OF CELLULAR AND GRANULAR MATERIALS

Y. Trivino<fnref n="1" />,<fnref n="2" />, V. Richefeu<fnref n="3" /> F. Radjai<fnref n="1" />, K. Lampoh<fnref n="2" /> and JY. Delenne<fnref n="2" />✉

**Abstract:** We present a 2D numerical model for simulating soft deformable cells and particles. To allow large deformations, the particle boundary is represented by an array of elements interacting at their nodes through elasto-plastic force laws governing their linear and angular relative displacements. Thus, the particle shape change is controlled by a combination of boundary stiffness and an internal particle pressure, which is coupled with the particle volume through a bulk modulus. The model is verified by comparing the results of our simulations of single particle deformation under load with theoretical predictions. We also showcase the potential power of our model in application to packings of particles involving both large particle shape change and dynamic rearrangements such as the compaction of a granular assembly of ultrasoft particles or the fracture of plant cell tissues.

## 1. Introduction

The field of soft particle dynamics represents a fascinating and broad area of materials science. Soft particles encompass a wide variety of materials, from cells and polymers to colloids, emulsions and liquid crystals with structures intermediate between crystalline solids and conventional fluids [1, 2]. This class of soft particle materials includes clays, which can be modeled as assemblies composed of flexible platelets at the nano scale [3, 4].

What is common to all these materials, despite their diversity, is their shared property of undergoing significant deformation under the effect of external load [5]. Exploring the behaviors and properties inherent in soft-particle materials has profound implications for many disciplines. In materials research, manipulating these flexible particles enables the creation of novel materials with tailored properties. For example, since the surface area between particles can change under applied pressure, such materials can be tailored with modulable adhesion in response to environmental actions [6]. Modeling the effects of particle deformability and shape change is essential for a better understanding of the fracture of soft materials. For example, for the comminution of cellular materials, the extensive cell deformations before breakage are expected to affect the process [7]. In the same way, deformation mechanisms on the scales of the aggregates in clays and the effect of composition on these mechanisms remain poorly understood [8].

In this context, we need new simulation methods going beyond the Discrete Element Method (DEM), which is a common numerical approach for the simulation of granular materials. In this paper, we introduce a numerical method based on frictional contact interactions between individual particles, as in DEM, but with deformable particles. In this method, the particle boundary is modeled by an array of elements interacting at their nodes through elasto-plastic force laws. The volume change behavior is also controlled through an internal pressure that is coupled with particle volume. We present an application of this method to the uniaxial compaction of a granular packing.

## 2. The core of the model

The model is based on the discrete element method (DEM) where the rigid objects are mass points or equivalently vertices. The novelty in the proposed approach lies in various interactions that

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

<footnote n="1">¹ Univ. Montpellier, CNRS, LMGC, Montpellier, France</footnote>

<footnote n="2">² Institut Agro Montpellier, INRAE, IATE, Montpellier, France</footnote>

<footnote n="3">³ Univ. Grenoble Alpes, CNRS, 3SR Lab., Grenoble, France</footnote>

<footnote>✉ jean-yves.delenne@inrae.fr</footnote>
contribute to the resultant forces at the vertices. The mass points are interconnected to form open or closed cells or particles, with a given wall in-plane thickness $\phi$ or radius $r = \phi/2$.

## 2.1. Indexing cells and their vertices

The numbering system used is of great concern for the comprehensiveness of the model, but also for the data structure in the computer code. Figure 1.a) illustrates how cells are organised by holding interconnected vertices. Thanks to a numbering of the vertices within each cell, called local numbering, a vertex is designated by a pair of index numbers $C_i|V_k$, where $i$ is the cell number and $v$ is the local vertex number.

<figure><img src="image_4.png" /><figcaption>Figure 1. a) Identification of vertices and bars in each cell. A vertex has a local number within the cell; to refer to it globally, the cell number followed by the local vertex number is required. For example, the bottom right vertex is identified as $C_j|V_4$, and any vertex in the cell $i$ is referred to as $C_i|V_k$. b) Representation of a cell $C_i|V_k$. Sliders are not displayed on this sketch</figcaption></figure>

## 2.2. Cell-internal forces and moments

We introduce longitudinal linear elastic springs between adjacent nodes and angular springs at each node between adjacent bars. These springs are plastified when they reach thresholds $F_b^Y$ and $M_v^Y$, respectively. The notation $\Delta\mathcal{X} \pm \mathcal{P}$ means that $\Delta\mathcal{X} = 0$ if $|\mathcal{X} + \Delta\mathcal{X}| \le \mathcal{P}$. In the longitudinal direction of the bar, we have

$$ \Delta F_b = \left[ -k_b \dot{l} dt \right] \pm F_b^Y \quad (1a) $$

$$ \Delta M_v = \left[ -k_v \dot{\theta} dt \right] \pm M_v^Y \quad (1b) $$

where $\dot{l} = (\ell_{t+\delta t} - \ell_t)/\delta t$, $k_b$ is the longitudinal stiffness, and $k_v$ is the angular stiffness.

## 2.3. Cell internal pressure

In order to represent the true physical behavior of systems such as cellular tissues, another important feature is the presence of a fluid in the cell. Let $V$ be the volume of the cell, $V_0$ the initial volume and $p$ the pressure. We set

$$ p = -K_{\Omega} \frac{V - V_0}{V_0}, \quad (2) $$

where $K_{\Omega}$ is the bulk modulus of the cell due to the presence of the fluid. The total bulk modulus $K$ is a combination of $K_{\Omega}$ and a contribution from the cell contour $K_{\partial\Omega}$, so that $K = K_{\partial\Omega} + K_{\Omega}$. It can
<figure><img src="image_5.png" /><figcaption>Figure 2. Snapshots of a simulated assembly of soft particles under uniaxial conditions inside a box.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Figure 3. Compression/consolidation curves for uniaxial compression-decompression simulations with a plastic yield moment $M^Y/(k_v) = 6 \times 10^{-4}$.</figcaption></figure>

be shown that in the case of a regular polygon composed of $N$ elements we have

$$K_{\partial\Omega} = \frac{k_b}{s} \sin \frac{\pi}{N}, \qquad (3)$$

where $s$ is the depth of the domain. For $N \to \infty$, $K_{\partial\Omega}$ tends to $K_{\partial\Omega}^{\infty} = ES/(2Rs)$, where $R$ is the radius of the particle and $S$ is the sectional area of the bars.

### 3. Results & Discussion

To showcase the performance and applicability of this model, we simulated the uniaxial compression of an assembly of soft particles as shown in figure 2. The sample is loaded vertically at a constant low speed. To create the initial sample, we used gravity to fill a box with regular elastic hexagons of high volume stiffness to avoid shape change during sample preparation. The particles are allowed to relax to static equilibrium. We used 500 particles whose internal pressure, the yield moment and the fluid bulk modulus were varied. The friction coefficient is zero between particles.

Figure 3 shows the void ratio $e$ of the sample as a function of the vertical stress $\sigma$ for $M^Y/k_v = 6 \times 10^{-4}$ and four values of the ratio $K_{\partial\Omega}/K_{\partial\Omega} \in \{0.001, 0.01, 0.1, 1\}$. We see that, for all $K$-ratios and starting from the same value of void ratio $e_i \simeq 0.2$, $e$ declines during the compression phase as a
function of the vertical stress $\sigma$. The value $e_{min}$ at $\sigma = \sigma_{max}$ slightly decreases for increasing $K$-ratio. During decompression, $e$ increases again from $e_{min}$ towards a constant final value $e_f$ as $\sigma$ vanishes. This relaxed void ratio $e_f$ increases with $K$-ratio. The lower level of the decompression plateau $e_f$ compared to the initial plateau $e_i$ indicates the plastification of the sample due to two mechanisms: 1) particle rearrangements and 2) permanent particle shape change due to node plasticity. Only in the case $K_{\Omega} = K_{\partial\Omega}$, the initial level of void ratio is recovered during decompression. This is due to the high level of the internal fluid bulk modulus $K_{\Omega}$, which implies that the particles tend to deform by shape change, thereby filling well the space (low values of $e_{min}$). However, this shape change is nearly uniform on the particles with small deviations of the curvature from a circular shape so that the angular plastic threshold $M^Y$ is not reached and the deformation is almost fully reversed upon unloading. In all other cases, the plastic threshold of the particles is reached and the particles partially keep their deformed shape upon unloading, depending on the value of the $K$-ratio.

Another related feature of the consolidation curves, which is worth to be briefly mentioned here, is that, at low $K$-ratios the particles tend to deform isotropically at low stress so that the particle diameters decrease. This effect is enhanced by zero friction between particles, allowing efficient rearrangement of the particles and their filling of the available space. Therefore, only at high values of $\sigma$ particle shapes begin to change. In contrast to the consolidation curve of rigid undeformable particles, in the case of soft particles the void ratio and solid fraction are independent parameters and a full description of the volume-change behavior requires also the variation of the solid fraction with stress.

## 4. Conclusion

In this paper, we introduced a model based on the discrete element method for the mechanical simulation of soft particles. The model is based on a set of equations describing the motion of nodes in a two-dimensional space. Initially designed for plant cell tissues, this model opens up new perspectives for the simulation of complex structures at the microscopic scale. For example, once calibrated, it can also be employed for the simulation of fibrous materials, soft aggregates in clays, and metallic powders.

## References

[1] Vlassopoulos, D. and Cloitre, M. (2014). Tunable rheology of dense soft deformable colloids. *Current opinion in colloid & interface science*, **19**, 561–574.

[2] Höhler, R. and Weaire, D. (2019). Can liquid foams and emulsions be modeled as packings of soft elastic particles? *Advances in colloid and interface science*, **263**, 19–37.

[3] Mitchell, J. K., Soga, K., et al. (2005). *Fundamentals of soil behavior*, vol. 3. John Wiley & Sons New York.

[4] Gea Pagano, A., Alonso-Marroquin, F., Ioannidou, K., Radjaï, F., and O'Sullivan, C. (2023). Clay micromechanics: mapping the future of particle-scale modelling of clay. *International society for soil mechanics and geotechnical engineering*, Porto, Portugal, Sep.

[5] Nezamabadi, S., Nguyen, T. H., Delenne, J.-Y., and Radjai, F. (2017). Modeling soft granular materials. *Granular Matter*, **19**, 1–12.

[6] Liu, K.-K. (2006). Deformation behaviour of soft particles: a review. *Journal of Physics D: Applied Physics*, **39**, R189.

[7] Guida, G., Casini, F., Viggiani, G. M., Ando, E., and Viggiani, G. (2018). Breakage mechanisms of highly porous particles in 1d compression revealed by x-ray tomography. *Géotechnique Letters*, **8**, 155–160.

[8] Carrier, B., Vandamme, M., Pellenq, R. J.-M., and Van Damme, H. (2014). Elastic properties of swelling clay particles at finite temperature upon hydration. *The Journal of Physical Chemistry C*, **118**, 8933–8943.