# Discontinuous shear thickening in dry granular materials induced by non-Coulombian friction

Denis Dumont, Francisco Rocha, Maxime Nicolas, Olivier Pouliquen

► To cite this version:

Denis Dumont, Francisco Rocha, Maxime Nicolas, Olivier Pouliquen. Discontinuous shear thickening in dry granular materials induced by non-Coulombian friction. Physical Review Fluids, 2025, 10 (6), pp.064302. 10.1103/PhysRevFluids.10.064302. hal-05099976

HAL Id: hal-05099976

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

Submitted on 5 Jun 2025

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
# Discontinuous Shear Thickening in Dry Granular Materials induced by non-Coulombian Friction

Denis Dumont, Francisco M. Rocha, Maxime Nicolas, and Olivier Pouliquen*

Aix Marseille Univ, CNRS, IUSTI, Marseille, France.

(Dated: May 5, 2025)

Over the last decade, the physics of discontinuous shear thickening (DST) in dense suspensions has undergone a paradigm shift. Growing evidence indicates that it stems from a frictional transition at the particle scale induced by repulsive forces between the suspended grains. However, in dry granular materials, interparticle friction can also exhibit a transition when the grains are coated, typically with polymers or lubricants widely used in industry. In this work, we use discrete-element-method simulations to demonstrate that DST occurs in dry granular materials when considering non-Coulombian behavior, where the interparticle friction varies with the contact normal force to mimic the influence of a coating. By performing pressure- and volume-imposed measurements, we reveal all regimes of dry shear thickening and present a mean-field model that predicts the observed scenarios by relating the bulk rheology to the microscopic friction law.

## I. INTRODUCTION

Dense suspensions, such as the mixture of cornstarch particles in water, can exhibit the spectacular phenomenon of discontinuous shear thickening (DST). When prepared at high enough packing fractions, these suspensions display an abrupt increase in effective viscosity when sheared, which may completely solidify the material [1–7]. The phenomenon of a violent increase in viscosity was first observed by Freundlich and Röder [8] in the 30s when studying the drag felt by an object within a dense suspension of quartz powders. However, its physical origin has remained an open question until recent years.

Over the last decade, many experimental and numerical studies have shown that DST stems from a frictional transition at the particle scale [9–13]. The key physical ingredient is the presence of a short-range repulsive force between the grains, which prevents solid contact at low stresses. The suspension then behaves as an ensemble of frictionless grains. However, when the material is subjected to high stresses and the particle pressure overcomes the repulsive stress, grains come into solid contact, activating the interparticle friction [9, 10, 14]. These observations have been used to rationalize DST theoretically in a mean-field model that considers a transition from a frictionless to a frictional state as the shear stress increases [14]. Although the studies mentioned above considered the presence of an interstitial fluid, they also revealed that in such dense states, hydrodynamic interactions play a minor role, and the change between different types of contact is responsible for the change in effective viscosity. This scenario suggests that a similar phenomenon can emerge even in the complete absence of a suspending fluid.

Dong and Trulsson [7], using discrete-element-method (DEM) simulations, have shown that the DST transition is indeed a very general phenomenon. It can be observed in suspensions across a wide range of Stokes numbers, from the viscous regime controlled by dissipation in the suspending fluid to inertial suspensions, where the fluid plays no role and the rheology is controlled by collisions between particles. The only necessary ingredient for DST is a frictionless-frictional transition at the interparticle scale.

In this paper, we follow Dong and Trulsson's work by studying the role of a variable friction coefficient and the occurrence of DST in a completely dry granular materials, without any suspending fluid. In various industrial applications, particles are often coated with liquids, lubricants, or polymers [15–17]. This coating alters the interparticle interactions, which are no longer governed solely by Coulombian friction. Tribology research has established that the presence of polymers or lubricants at particle contact points modifies the friction coefficient, making it dependent on both the normal force and the sliding velocity in intricate ways [18–23]. While the impact of non-Coulombian friction on rheology is well-documented for suspensions [24, 25], it is much less understood for dry granular materials [26, 27].

The rheology of dry, Coulombian granular materials, interacting only by collision has been extensively studied in the past two decades [28]. Under constant confining pressure, the leading-order dynamics is well described by the so-called $\mu(I)$ rheology, in which the bulk friction coefficient $\mu(I, \mu_p)$ and the solids volume fraction $\phi(I, \mu_p)$ are functions of the dimensionless inertial number $I$ and of the interparticle friction coefficient $\mu_p$ [29, 30]. In configurations where the volume of the sample is fixed, the constitutive laws are Bagnoldians, i.e. shear stress and pressure scale with the shear

<footnote>* olivier.pouliquen@univ-amu.fr</footnote>
rate squared [29, 30]. Both approaches are entirely equivalent and the constitutive relations in the volume-imposed configuration can be computed from the ones at imposed pressure and vice versa [29]. This description of the rheology of dry granular materials only holds for constant friction coefficient and in the limit of rigid particles.

In this study, we investigate how the dependence of the friction coefficient on the normal force affects the rheology of granular media. Through DEM simulations, we demonstrate that a discontinuous shear thickening behavior can be observed in non-Coulombian dry granular materials, analogous to that seen in dense suspensions, but without the presence of a suspending fluid. We conduct numerical experiments in both pressure- and volume-imposed configurations to show that all regimes of shear thickening, observed both numerically [9, 10] and experimentally [11–13] in dense suspensions, are also present in dry systems due to the non-constant friction coefficient. We present a mean-field model that incorporates the normal-force dependence of interparticle friction into the bulk friction coefficient, successfully capturing all the regimes observed in the simulations.

It is worth noting that the DST transition reported in this paper and in Dong and Trulsson [7], which is controlled by a frictional transition at the interparticle scale, is distinct from the ‘elastic’ DST reported in previous numerical studies of dry granular materials [31–34]. Those studies demonstrated that, with a constant friction coefficient but finite stiffness, a shear thickening between a flowing branch and a jammed state is observed in a narrow range of volume fractions near the jamming transition, linked to the elasticity of the grains.

## II. SIMULATION METHODOLOGY

The DEM simulations are performed using the open source software LAMMPS [35]. The granular medium consists of spherical beads of diameter $d = 10^{-3}$ m and density $\rho_p = 2500$ kg/m$^3$. The interaction forces between grains are described using the Hertz-Mindlin model [36], with a Young’s modulus $E \in [10^5, 10^7]$ Pa, and a Poisson ratio $\nu = 0.3$. The restitution coefficient is equal to $e = 0.2$.

The novelty of this work is to consider a non-Coulombian frictional behavior. Typically, in presence of a coating, $\mu_p$ depends on both $F_n$ and the relative velocity. In this study we focus solely on the role of $F_n$ and consider that the interparticle friction $\mu_p$ is an explicit function of the contact normal force $F_n$. We assume that a transition from a low friction limit $\mu_1$ to a high friction limit $\mu_2$ occurs as $F_n$ increases (see Fig. 1a). To model such phenomenology, we choose a friction coefficient:

$$ \mu_p(F_n) = \frac{\mu_2 + \mu_1}{2} + \frac{\mu_2 - \mu_1}{2} \tanh \left( \frac{F_n - F^*}{S^*} \right). \quad (1) $$

As shown in Fig. 1(a), the frictional transition between the two limit cases occurs at a critical force $F^*$ and spreads over a force $S^*$. Hereafter in the paper, we choose $S^*/F^* \ll 1$, meaning that the transition is sufficiently sharp such that results are independent of $S^*$. In addition, the lower and upper friction limits are set to $\mu_1 = 0.2$ and $\mu_2 = 0.7$, respectively.

To investigate how the rheology of the granular medium is affected by the normal-force-dependent friction coefficient, we consider the plane shear configuration of Fig. 1(b). The granular material is confined between two rough and rigid plates made of bonded grains. The bottom plate is held fixed while the top plate moves in the horizontal direction with a fixed velocity $V_{top}$, prescribing the shear rate $\dot{\gamma}$. The surface area of the plates is $(L_x, L_y) = (20d, 20d)$ and the thickness of the layer is typically $20d$ at rest. Periodic boundary conditions are imposed in the flow ($x$) and lateral ($y$) directions.

To perform pressure-imposed measurements, a constant pressure $P$ is applied to the top plate, whilst for volume-imposed experiments the vertical position of the top plate is prescribed at all times. Simulations are carried for long enough times such that the system achieves a steady state, and we have checked that in both configurations the flow remains uniformly sheared, with a uniform shear rate $\dot{\gamma}$ and volume fraction $\phi$ across the layer [37, 38]. To analyze the constitutive laws, we extract from each simulation the shear stress $\tau$, the normal stress $P$, the volume fraction $\phi$, and the shear rate $\dot{\gamma}$. In the following, we analyze the influence of the non-Coulombian particle friction coefficient in both pressure- and volume-imposed configurations.

## III. PRESSURE-IMPOSED RHEOLOGY

We first consider the granular medium sheared under constant pressure $P$. In the Coulombian case of constant $\mu_p$ and in the limit of rigid particles, the granular rheology reduces to two constitutive laws: the bulk friction, i.e. the ratio of shear to normal stresses $\tau/P = \mu_{Coul}(I, \mu_p)$, and the volume fraction $\phi_{Coul}(I, \mu_p)$, both depending on
the inertial number $I = \dot{\gamma}d/\sqrt{P/\rho_p}$ and the interparticle friction $\mu_p$ [39, 40]. The red and blue squares in Fig. 2 correspond to the constitutive laws for the two limit values of $\mu_p$ prescribed by Eq. (1), $\mu_p = \mu_1$ in blue and $\mu_p = \mu_2$ in red.

When considering a normal-force-dependent friction coefficient, the new force scale $F^*$ controlling the transition between low and high interparticle friction introduces a new dimensionless number $F^*/Pd^2$, which compares $F^*$ to the typical normal force between the particles $F_n \approx Pd^2$. We thus expect the rheology to be controlled by a friction law $\mu(I, F^*/Pd^2)$ and a volume fraction law $\phi(I, F^*/Pd^2)$. The results of the simulations carried out at different $F^*$, $P$, and $\dot{\gamma}$ are depicted by the circles in Fig. 2, where the different colours correspond to runs at same values of $F^*/Pd^2$.

At this point, three observations can be made. First, the different points obtained for different values of pressure $P$ and different critical forces $F^*$ but keeping the same ratio $F^*/Pd^2$ (same color) collapse onto single curves, showing that the rheology is indeed controlled by this ratio as predicted by dimensional analysis. Secondly, changing $F^*/Pd^2$ from 0.01 to 5, the $\mu(I, F^*/Pd^2)$ and $\phi(I, F^*/Pd^2)$ laws evolve from the high-friction curve, corresponding to the constant $\mu_p = \mu_2$, to the low-friction curve associated with $\mu_p = \mu_1$. This is qualitatively expected, as when $Pd^2$ is much larger than $F^*$, most contacts involve normal forces $F_n > F^*$, which leads to a high friction coefficient in the majority of the contacts. Conversely, when $Pd^2$ is much smaller than $F^*$, most contacts lie in the lower frictional state. Thirdly, for intermediate values of $F^*/Pd^2$, the curves $\mu(I, F^*/Pd^2)$ and $\phi(I, F^*/Pd^2)$ tend to start close to the low-friction limit and approach the high-friction limit as $I$ increases. This occurs because as $I$ increases, the system dilates, decreasing the average number of contacts per particle. Hence, to balance the constant confining pressure, the average normal force between particles increases, leading the system to the high-friction state. The transition observed between the low-friction and the high-friction curves is reminiscent of what is observed with shear thickening suspensions studied under pressure-imposed condition [11, 12].

## IV. VOLUME-IMPOSED RHEOLOGY

Now, we turn to the volume-imposed configuration, where the positions of the top and bottom plates are kept fixed at all time, fixing the average volume fraction $\phi$. In contrast with the pressure-imposed case, the constitutive laws are expressed in terms of how the pressure $P$ and the shear stress $\tau$ depend on the shear rate $\dot{\gamma}$ and volume fraction $\phi$. For Coulombian materials with constant interparticle friction coefficient $\mu_p$ and in the limit of rigid particles, dimensional analysis imposes that $P$ and $\tau$ vary with the square of the shear rate, i.e. $P = \rho_p d^2 \dot{\gamma}^2 f_1(\phi)$ and $\tau = \rho_p d^2 \dot{\gamma}^2 f_2(\phi)$. The functions $f_1(\phi)$ and $f_2(\phi)$ are equivalent to the shear and normal dimensionless viscosities in dense suspensions and they diverge when $\phi$ approaches the maximum volume fraction $\phi_c$. This scaling breaks down at high shear rates

<figure><img src="image_2.png" /><figcaption>FIG. 1. (a) The interparticle friction coefficient $\mu_p$ as a function of the normal force $F_n$ with $\mu_1 = 0.2$, $\mu_2 = 0.7$, $F^* = 30N$ and $S^* = 3N$ given by Eq. (1). (b) Snapshot of the simulated shear cell, the walls are in blue, the flowing granular layer is in red.</figcaption></figure>
when stresses are sufficiently high to deform the particles. The system then enters into an elastic regime, where the pressure varies with the square root of the shear rate [41, 42].

In the following, we present the influence of the non-Coulombian interparticle friction on the rescaled pressure. Equivalent behaviour is observed for the shear stress $\tau$, so it is omitted here. In Fig. 3(a), the normalized pressure $P/\rho_p d^2 \dot{\gamma}^2$ is plotted as a function of the normalized shear rate $\dot{\gamma} d(\rho_p/E)^{1/2}$ for $\phi = 0.576$ (main plot) and $0.587$ (inset). The blue and red square symbols correspond to simulations with a constant interparticle friction coefficients $\mu_p = \mu_1$ and $\mu_p = \mu_2$, respectively. For $\phi = 0.576$, both curves are constant as predicted by the dimensional analysis but they are shifted in pressure, showing the influence of the interparticle friction $\mu_p$. For $\phi = 0.587$ (inset), both curves present a plateau at low shear rates corresponding to the rigid limit. For $\mu_p = 0.7$ (red squares), a shear-thinning branch is observed with slope $-3/2$, which corresponds to the elastic regime.

If we now consider the case of a normal-force-dependent interparticle friction, the system transitions from the low-to the high-friction curve when increasing the shear rate, as shown by the circles in Fig. 3(a), where the colours correspond to different values of the critical force $F^*$. The higher the critical force $F^*$, the higher the associated critical shear rate at which shear thickening occurs. From dimensional analysis, this critical shear rate should scale with $\sqrt{F^*/\rho_p/d^2}$.

These different curves from Fig. 3a, corresponding to different $F^*$ indeed collapse onto a master curve, when the pressure is plotted as a function of the dimensionless shear rate $\dot{\gamma}d^2/\sqrt{F^*/\rho_p}$, as shown in Fig 3(b). A deviation is observed for the two smallest values of $F^*$ (shown in light blue), which can be attributed to the non negligible width of the frictional transition $S^*$ (Fig. 1a) at such small $F^*$ values. It is also important to note that this collapse holds only in the rigid limit. As soon as the elasticity of the grains starts to play a role, another force scale intervenes in the rescaling of the shear rate. This is observed in the limit of high shear rates for $\phi = 0.587$ [inset of Fig. 3(b)].

We proceed to investigate the influence of the volume fraction $\phi$ in Fig. 3(c), where $P/\rho_p d^2 \dot{\gamma}^2$ is plotted against $\dot{\gamma}d^2/\sqrt{F^*/\rho_p}$ for different values of $\phi$, which is equivalent to the behavior of a suspension normal viscosity as a function of the shear rate [9, 10]. At low volume fractions, the system transitions continuously from the lower to the higher frictional branch, and $P/\rho_p d^2 \dot{\gamma}^2$ is a monotonically increasing function of the shear rate, leading to a continuous shear-thickening regime. When $\phi = 0.583$, the slope of the curve becomes vertical, which leads to a drastic jump in viscosity around $\dot{\gamma}_c d^2/\sqrt{F^*/\rho_p} \approx 0.06$. This threshold corresponds to the onset of discontinuous shear thickening in dense suspensions. For all higher volume fractions, the flow curve becomes S-shaped with a shear-weakening branch connecting the lower and higher frictional branches. This branch is equivalent to the shear-weakening regime in colloidal suspensions, which is intrinsically unstable, and the reason for the discontinuous jump in shear and normal stresses at a critical shear rate [14]. Note that a second shear-weakening regime is observed in the limit of high shear rates, after the material has already transitioned to the higher frictional state. This shear-weakening

<figure><img src="image_3.png" /><figcaption>FIG. 2. (a) Bulk friction coefficient $\mu(I)$ and (b) volume fraction $\phi(I)$ laws for different values of $P$, $F^*$ and $\dot{\gamma}$. Square symbols: data obtained for constant friction coefficient $\mu_p = \mu_1$ (blue) and $\mu_p = \mu_2$ (red). Circles are data obtained for normal-force-dependent friction coefficient (Eq. (1)) colored as function of $F^*/Pd^2$. Dashed lines are fit given by Eqs. (4) and (5). Continuous lines are predictions of the mean-field model (see text) colored as function of $F^*/Pd^2$.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>FIG. 3. (a) Rescaled pressure, $P/(\rho_p \dot{\gamma}^2 d^2)$ as a function of the shear rate rescaled by the Young's modulus, $E$ for $\phi = 0.576$ and different critical forces $F^*/Ed^2 = 5 \times 10^{-7}, 10^{-6}, 2.5 \times 10^{-6}, 5 \times 10^{-6}, 10^{-5}, 2.5 \times 10^{-5}$ (from light to dark purple circles). Inset: same graph for $\phi = 0.587$ and $F^*/Ed^2 = 2.5 \times 10^{-6}, 5 \times 10^{-6}, 10^{-5}$ (from light to dark orange circles). Blue and red squares correspond to the case of constant $\mu_p = 0.2$ and $\mu_p = 0.7$. (b) Rescaled pressure, $P/(\rho_p \dot{\gamma}^2 d^2)$ as a function of the dimensionless shear rate $\dot{\gamma}d^2/\sqrt{F^*/\rho_p}$ for the same data than in (a). (c) Rescaled pressure as a function of rescaled shear rate for a fixed $F^*/Ed^2 = 10^{-5}$ and different $\phi$. Data are represented by disks, the predictions of the mean-field model for rigid particles by dotted lines, predictions of the model for soft particles by the solid line curves.</figcaption></figure>

regime is, however, induced by the elasticity of the grains, where $P/\rho_p d^2 \dot{\gamma}^2 \propto (\dot{\gamma}d\sqrt{\rho_p/E})^{-3/2}$ (black line). This plot is a direct equivalent of the shear-thickening curves observed experimentally and numerically in dense colloidal suspensions [9, 10, 13, 14], and all the shear-thickening regimes are recovered with the non-Coulombian nature of the contact.

## V. MEAN-FIELD MODEL

To explain the observations, we propose hereafter a mean field approach to predict the macroscopic shear-thickening rheology from the microscopic interparticle friction law $\mu_p(F_n)$. The key assumption of our model is that the rheology of a non-Coulombian granular medium can be described by the Coulombian granular rheological laws, where the interparticle friction is replaced by an effective value $\tilde{\mu}_p$, which depends on $I$ and $F^*/Pd^2$. This hypothesis is close to the one used in [25] and similar to the main assumption of the Wyart & Cates model [14] for describing DST in suspensions. In their case, the rheology is assumed to be given by the newtonian rheology with an effective maximum volume fraction $\phi_c$ that depends on the ratio between the repulsive force and the average contact force.

To do so, we first characterize the Coulombian rheology $\mu_{Coul}(I, \mu_p)$ and $\phi_{Coul}(I, \mu_p)$ for all value of constant $\mu_p$. Secondly, we infer what would be the distribution of friction coefficient in the case of a normal-force-dependent friction coefficient $\mu_p(F_n)$ assuming that the distribution of $F_n$ is the one from the Coulombian system. Then, an effective average friction coefficient $\tilde{\mu}_p$ is defined, which depends on $I$ and $F^*/Pd^2$. Injecting $\tilde{\mu}_p(I, F^*/Pd^2)$ into the Coulombian rheology yields the non-Coulombian constitutive laws:

$$ \mu(I, F^*/Pd^2) = \mu_{Coul}(I, \tilde{\mu}_p(I, F^*/Pd^2)), \qquad (2) $$

$$ \phi(I, F^*/Pd^2) = \phi_{Coul}(I, \tilde{\mu}_p(I, F^*/Pd^2)). \qquad (3) $$

We thus start by systematically computing the friction $\mu_{Coul}(I, \mu_p)$ and volume fraction $\phi_{Coul}(I, \mu_p)$ for a constant $\mu_p \in [0.05, 0.9]$ (the detailed results are presented in Appendix A). The two constitutive relations are fitted according to the following functional forms proposed in the literature [30]:

$$ \mu_{Coul}(I, \mu_p) = \mu_c(\mu_p) + \frac{\mu_\infty(\mu_p) - \mu_c(\mu_p)}{I_0/I + 1}, \qquad (4) $$

$$ \phi_{Coul}(I, \mu_p) = \phi_c(\mu_p) - (\phi_c(\mu_p) - \phi_m(\mu_p))I, \qquad (5) $$

where $I_0 = 0.5$ and the functions $\mu_c(\mu_p)$, $\mu_\infty(\mu_p)$, $\phi_c(\mu_p)$, and $\phi_m(\mu_p)$ are obtained from a best fit procedure and are well approximated by an exponential dependence with $\mu_p$ (see Appendix A for details).
Then, we define the effective average friction coefficient $\tilde{\mu}_p$ for the non-Coulombian system as

$$ \tilde{\mu}_p = \int_0^\infty \mu_p(F_n) \mathcal{P}(F_n) dF_n, \qquad (6) $$

where $\mu_p(F_n)$ is the interparticle friction law given by Eq. (1) and $\mathcal{P}(F_n)$ is the probability distribution of normal forces within the material. In this mean field approach, we assume that the probability distribution $\mathcal{P}(F_n)$ can be taken from the Coulombian system. In Appendix A, we show that a good approximation of the normal force distribution for constant $\mu_p$ is an exponential distribution:

$$ \mathcal{P}(F_n) = \frac{1}{F_0} \exp(-F_n/F_0) \qquad (7) $$

characterized by a single parameter $F_0$, which is the mean normal force. It should be noted that, as mentioned in previous studies [43–45], the normal force distribution deviates from a purely exponential distribution when normal forces are low (particularly when $I$ is low), a property neglected in what follows. The parameter $F_0$ a priori is proportional to $Pd^2$, with a coefficient depending on both $I$ and $\mu_p$. In the Appendix A, we compute the average normal force $\langle F_n \rangle/Pd^2$ as a function of $I$ and $\mu_p$ to show that the role of $\mu_p$ is weak. We then propose the following expression for $F_0/Pd^2$ as a function of $I$ only:

$$ \frac{F_0}{Pd^2} = 1.2 + 23I^2. \qquad (8) $$

Substituting the normal force distribution, Eq. (7), with $F_0$ given by Eq. (8) in the definition of $\tilde{\mu}_p$, Eq. (6), gives the expression for the effective friction $\tilde{\mu}_p(I, F^*/Pd^2)$. To compute the non-Coulombian rheology, we inject $\tilde{\mu}_p$ in Eqs (2) and (3) where the expressions for $\mu_{Coul}$ and $\phi_{Coul}$ are given by Eqs. (4) and (5).

The predictions of this mean field approach for $\mu$ and $\phi$ are shown as solid lines in Fig. 2(a,b) for the pressure imposed configuration. The model captures well the transition from the low- to the high-friction branch when $F^*/Pd^2$ decreases and $I$ increases, and matches quantitatively the DEM numerical data. A discrepancy between the theoretical model and the data is observed at low inertial number $I$, probably due to the assumption of the exponential distribution for the normal force, which is less accurate at low $I$.

To predict the constitutive equations in the volume-imposed case, we just need to invert (numerically) the relation $\phi(I, F^*/Pd^2)$ knowing the definition of $I$, to get $P/\rho d^2 \dot{\gamma}^2$ as a function of $\phi$ and $\dot{\gamma}d^2/\sqrt{F^*/\rho_p}$. The predictions are shown as dotted lines in Fig. 3(b,c). The discontinuous shear thickening and the transition between low to high friction regime are well captured by the model, with typical S-shaped curves. Since the DEM simulations are conducted at fixed shear rate, the unstable branches of the S-shape curves are not accessible and one observes sudden jumps in pressure when increasing $\dot{\gamma}$. To capture the shear thinning observed at high shear rates corresponding to the elastic regime [see inset of Fig. 3(a)] the model needs to be improved. We thus use the following empirical relation for the pressure capturing the crossover between the rigid and soft regime:

$$ \frac{P}{E} = a \left( \left( \hat{\gamma}^2 \right)^{-1} + \left( b \hat{\gamma}^{1/2} \right)^{-1} \right)^{-1}, \qquad (9) $$

where $\hat{\gamma} = \rho_p^{1/2} d\dot{\gamma}/E^{1/2}(\phi_c(\mu_p) - \phi)$. The numerical parameters $a = 0.025$ and $b = 0.029$ are obtained via a fitting procedure (see Appendix A for details). Substituting $\mu_p$ by $\tilde{\mu}_p$ in Eq. (9) we obtain numerically $P$ as a function of $\phi$ and $\dot{\gamma}d^2/\sqrt{F^*/\rho_p}$. The results are plotted in Fig. 3(b,c) as solid lines and capture the behaviour over the entire range of shear rates, from the shear thickening induced by the frictional transition, to the shear thinning due to the particle softness.

## VI. CONCLUSION

To conclude, we have shown that discontinuous shear thickening can indeed occur in dry granular media if the interparticle friction transitions from a low to a high frictional state when subjected to higher contact normal loads. Performing discrete-element-method simulations in both pressure- and volume-imposed flow configurations, we demonstrate the effect of the non-Coulombian interparticle friction in the macroscopic bulk rheology. In the pressure-imposed configuration, the normal-force-dependent interparticle friction induces a transition from a lower to a higher frictional state as the material is subjected to higher normal stresses, similar to what is observed in dense colloidal suspensions. When describing the system in a volume-imposed framework, the system transitions smoothly between frictional
states when the volume fraction is low, analogously to a continuous shear-thickening regime in suspensions. On the other hand, when the non-Coulombian material is prepared in a dense packing, the flow curves become S-shaped, leading to a discontinuous shear-thickening regime due to the presence of a shear-weakening, unstable region across the transition. Our simulations also capture the elastic regime in the limit of high shear rates when the level of stresses is sufficient to deform the particles. This leads to a second velocity-weakening regime, even after the system has already become fully frictional. Here, we propose a phenomenological mean-field model inspired by the Wyart & Cates model [14] for colloidal suspensions. In our model, we depart from the rheology of dry, Coulombian granular materials and introduce an explicit dependence on the interparticle friction coefficient, which depends on the stress imposed on the material. We also present an extension of the model to account for the transition in scaling when the material enters the elastic regime. Our model can capture quantitatively all the regimes observed in the simulations, from continuous to discontinuous shear-thickening, as well as the elastic regime.

Importantly, our work demonstrates that discontinuous shear thickening does not require a suspending fluid to occur. Thus, shear thickening is not a unique feature of granular suspensions and can be observed in a variety of systems [46, 47]. Although the current numerical investigation uses a simplified model of granular material, we believe these findings can be of interest in potential applications. To enhance the flowability of powders, coatings are often used to lubricate and reduce friction at low stresses [48]. A parallel can be drawn with suspensions, where the use of plasticizers in concrete or lecithin in chocolate allows for better flowability and prevents jamming at low stresses while maintaining a high concentration of particles and the desired resistance at high stresses [49]. Controlling particle interactions through surface modifications to enhance flow properties in dry granular media is thus of interest for industrial applications. Further investigation into how surface coatings influence the rheology of granular materials would be a valuable avenue for future research.

## ACKNOWLEDGMENT

This project has received funding from the European Research Council (ERC) under the European Union's Horizon (Grant 101097842). The authors also acknowledge financial support from International Fine Particle Research Institute (IFPRI) and the Agence Nationale de la Recherche under RheoCom (ANR-22-CE06-0020) grant. This work has benefited from fruitful discussions with Vincent Bertin.

## Appendix A: Flow properties for constant friction coefficient

To obtain the two constitutive laws $\mu(I, \mu_p)$ and $\phi(I, \mu_p)$ necessary to derive the mean field approach, we systematically varied $\mu_p$ in the range 0.05 to 0.9. Results are presented in Fig. 4(a,b). By fitting $\mu(I)$ and $\phi(I)$ using Eqs. (2) and (3) for the different $\mu_p$, we obtained the parameters $\mu_c(\mu_p)$, $\mu_\infty(\mu_p)$, $\phi_c(\mu_p)$, $\phi_m(\mu_p)$ as shown in Fig. 4(c-f), where we also plot results from previous studies [39, 40]. The parameter $\phi_m$ appears to be independent of $\mu_p$, whereas the other parameters are well fitted by an exponential function of $\mu_p$ shown as solid lines (see the caption of the figure for the quantitative expression of $\mu_c(\mu_p)$, $\mu_\infty(\mu_p)$, $\phi_c(\mu_p)$).

To derive the mean field approach and compute the average friction coefficient between the particles, we need to know the distribution of normal forces in a granular assembly. Fig. 5(a) shows how the distribution $P(F_n)$ varies when increasing the inertial number $I$. The distribution is close to an exponential distribution (Eq. (7)) characterised by its mean value $F_0$. From dimensional analysis, one expect $F_0$ to scale with $Pd^2$, the only force scale in the system. In Fig. 5(b), the normal force averaged in all the contacts $(F_n)/Pd^2$ is plotted as a function of the inertial number $I$ for various values of $\mu_p$. To first order, one can consider that $F_0$ is independent of $\mu_p$ and well fitted by Eq. (8), as shown by the magenta solid line.

At high stresses as reached after the shear thickening transition, the elasticity of the grains can no longer be neglected. An empirical law (Eq. (9)) is proposed for the evolution of the pressure with the shear rate and the volume fraction. In Fig. 5(c), the evolution of the normalised pressure $P/E$ with the normalised shear rate $\dot{\gamma}$ is plotted for different interparticle friction coefficients $\mu_p = 0.2, 0.5, 0.7$ and Young's modulus $E$. The dotted line represents the inertial regime ($P \sim \dot{\gamma}^2$) and the dashed line the elastic regime ($P \sim \dot{\gamma}^{1/2}$). The red line is Eq. (9).

<footnote>[1] E. Brown and H. M. Jaeger, *Shear thickening in concentrated suspensions: phenomenology, mechanisms and relations to jamming*, Reports on Progress in Physics, **77**, no. 4, 046602 (2014).</footnote>
<figure><img src="image_5.png" /><figcaption>FIG. 4. (a) $\mu(I)$ and (b) $\phi(I)$ for constant interparticle friction $\mu_p$ ranging between 0.05 and 0.9. Solid lines correspond to best fits obtained from Eqs. (2) and (3). Evolution of (c) $\mu_c$, (d) $\mu_\infty$ (e) $\phi_c$ and (f) $\phi_m$ as a function of $\mu_p$, with $I_0 = 0.5$. The solid lines represent the exponential fits $\mu_c = 0.355 - 0.207\exp(-6.956\mu_p)$, $\mu_\infty = 0.742 - 0.300\exp(-3.869\mu_p)$, $\phi_c = 0.584 + 0.055\exp(-3.868\mu_p)$. The dashed line corresponds to $\phi_m = 0.454$. Data from previous works of Estrada et al. [39] and Man et al. [40] are shown for comparison.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>FIG. 5. (a) Rescaled distribution of the normal forces $F_0P(F_n)$ as a function of $F_n/F_0$, where $F_0$ is the mean force and the solid line represents the exponential distribution obtained for $\mu_p = 0.5$ and different $I$; Inset: Same data without the rescaling. (b) Evolution of the normal force averaged over all the contacts, $\langle F_n \rangle$, as a function of the inertial number $I$ for different values of $P$ and $\mu_p$. The magenta curve corresponds to the best fit $F_0 = Pd^2(1.2 + 23I^2)$ used in the mean field model. (c) Evolution of the normalized pressure $P/E$ with the normalised shear rate $\hat{\gamma} = \rho_p^{1/2}d\dot{\gamma}/E^{1/2}(\phi_c(\mu_p) - \phi)$ for $\mu_p = 0.2, 0.5, 0.7$ (upside-down triangle, circle, triangle) and Young's modulus $E$ (colour). The dotted and dashed lines represents the inertial regime ($P \sim \dot{\gamma}^2$) and the elastic regime ($P \sim \dot{\gamma}^{1/2}$) respectively. The red line is given by the Eq. (9) and the value of $a = 0.025$ and $b = 0.029$ is obtained from the best fit.</figcaption></figure>

[2] J. F. Morris, *Shear thickening of concentrated suspensions: Recent developments and relation to other phenomena*, Annual Review of Fluid Mechanics, **52**, no. 1, 121-144 (2020).

[3] C. Ness, R. Seto, and R. Mari, *The physics of dense suspensions*, Annual Review of Condensed Matter Physics, **13**, no. 1, 97-117 (2022).

[4] Holmes, C. B., M.E. Cates, M. Fuchs, and P. Sollich, *Glass transitions and shear thickening suspension rheology*, Journal of Rheology, **49**, no. 1, 237-269 (2005).
[5] W.T. Kranz, F. Frahsa, A. Zippelius, M. Fuchs, and M. Sperl, Rheology of inelastic hard spheres at finite density and shear rate, Physical Review Letters, 121, no. 14, 148002 (2018).

[6] C. Ness, and J. Sun, Shear thickening regimes of dense non-Brownian suspensions, Soft matter, 12, no. 3, 914-924 (2016).

[7] J. Dong, and M. Trulsson, Unifying viscous and inertial regimes of discontinuous shear thickening suspensions, Journal of Rheology, 64, no. 2, 255-266 (2020).

[8] H. Freundlich and H. L. Röder, Dilatancy and its relation to thixotropy, Transactions of the Faraday Society, 34, 308-316 (1938)

[9] R. Mari, R. Seto, J. F. Morris, and M. M. Denn, Shear thickening, frictionless and frictional rheologies in non-Brownian suspensions, Journal of Rheology, 58, no. 6, 1693-1724 (2014).

[10] R. Seto, R. Mari, J. F. Morris, and M. M. Denn, Discontinuous shear thickening of frictional hard-sphere suspensions, Physical Review Letters, 111, no. 21, 218301 (2013).

[11] J. Dong and M. Trulsson, Analog of discontinuous shear thickening flows under confining pressure, Physical Review Fluids, 2, no. 8, 081301 (2017).

[12] B. Etcheverry, Y. Forterre, and B. Metzger, Capillary-stress controlled rheometer reveals the dual rheology of shear-thickening suspensions, Physical Review X, 13, no. 1, 011024 (2023).

[13] C. Clavaud, A. Bérut, B. Metzger, and Y. Forterre, Revealing the frictional transition in shear-thickening suspensions, Proceedings of the National Academy of Sciences, 114, no. 20, 5147-5152 (2017).

[14] M. Wyart and M. E. Cates, Discontinuous shear thickening without inertia in dense non-Brownian suspensions, Physical Review Letters, 112, no. 9, 098302 (2014).

[15] A. Gans, O. Pouliquen, and M. Nicolas, Cohesion-controlled granular materials, Physical Review E, 101, 032904 (2020).

[16] A. Gans, A. Abramian, P. Y. Lagrée, M. Gong, A. Sauret, O. Pouliquen, and M. Nicolas, Collapse of a cohesive granular column, Journal of Fluid Mechanics, 959, A41 (2023).

[17] R. S. Sharma, W. Sarlin, L. Xing, C. Morize, P. Gondret, and A. Sauret, Effects of interparticle cohesion on the collapse of granular columns, Physical Review Fluids, 9, no. 7, 074301 (2024).

[18] J. Klein, E. Kumacheva, D. Mahalu, D. Perahia, and L. J. Fetters, Reduction of frictional forces between solid surfaces bearing polymer brushes, Nature, 370, no. 6491, 634-636 (1994).

[19] L. Bureau and L. Léger, Sliding friction at a rubber/brush interface, Langmuir, 20, no. 11, 4523-4529 (2004).

[20] J. H. H. Bongaerts, K. Fourtouni, and J. R. Stokes, Soft-tribology: Lubrication in a compliant PDMS-PDMS contact, Tribology International, 40, no. 10-12, 1531-1542 (2007).

[21] Z. Ye, P. Egberts, G. H. Han, A. T. Johnson, C. R. Carpick, and A. Martini, Load-dependent friction hysteresis on graphene, ACS Nano, 10, no. 5, 5161-5168 (2016).

[22] N. Kumar, A. T. Kozakov, T. R. Ravindran, S. Dash, and A. K. Tyagi, Load dependent friction coefficient of crystalline graphite and anomalous behavior of wear dimension, Tribology International, 88, 280-289 (2015).

[23] C. E. Tsai and J. C. Tsai, Dynamical force measurements for contacting soft surfaces upon steady sliding: Fixed-depth tribology, Physical Review E, 109, no. 6, 064802 (2024).

[24] N. Fernandez, R. Mani, D. Rinaldi, D. Kadau, M. Mosquet, H. Lombois-Burger, J. Cayer-Barrioz, H. J. Herrmann, N. D. Spencer, and L. Isa, Microscopic mechanism for shear thickening of non-Brownian suspensions, Physical Review Letters, 111, no. 10, 108301 (2013).

[25] L. Lobry, E. Lemaire, F. Blanc, S. Gallier, and F. Peters, Shear thinning in non-Brownian suspensions explained by variable friction between particles, Journal of Fluid Mechanics, 860, 682-710 (2019).

[26] J. C. Tsai, G. H. Huang, and C. E. Tsai, Signature of transition between granular solid and fluid: Rate-dependent stick slips in steady shearing, Physical Review Letters, 126, no. 12, 128001 (2021).

[27] C. E. Tsai, W. C. Li, H. C. Fan-Chiang, P. Y. Hsiao, and J. C. Tsai, Two types of quaking and shear unjamming: State diagram for soft granular particles in quasistatic shear, Physical Review Research, 6, no. 2, 023065 (2024).

[28] K. Kamrin, K. M. Hill, D. I. Goldman, and J. E. Andrade, Advances in modeling dense granular media, Annual Review of Fluid Mechanics, 56, no. 1, 215-240 (2024).

[29] GDR MiDi, On dense granular flows, The European Physical Journal E, 14, 341 (2004).

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

[31] C. Ness, and S.M. Fielding, Nonmonotonic Constitutive Curves and Shear Banding in Dry and Wet Granular Flows, Physical Review Letters, 134, no. 6, 038201 (2025).

[32] M. Otsuki, and H. Hayakawa, Avalanche contribution to shear modulus of granular materials, Physical Review E, 90, no. 4, 042202 (2014).

[33] M. Grob, C. Heussinger, and A. Zippelius, Jamming of frictional particles: A nonequilibrium first-order phase transition, Physical Review E, 89, no. 5, 050201 (2014).

[34] M. Grob, A. Zippelius, and C. Heussinger, Rheological chaos of frictional grains, Physical Review E, 93, no. 3, 030901 (2016).

[35] A. P. Thompson, H. M. Aktulga, R. Berger, D. S. Bolintineanu, W. M. Brown, P. S. Crozier, and S. J. Plimpton, LAMMPS-a flexible simulation tool for particle-based materials modeling at the atomic, meso, and continuum scales, Computer Physics Communications, 271, 108171 (2022).

[36] R. D. Mindlin and H. Deresiewicz, Elastic spheres in contact under varying oblique forces, J. Appl. Mech., 20, 327-344 (1953).

[37] F. Da Cruz, S. Emam, M. Prochnow, J. N. Roux, and F. Chevoir, Rheophysics of dense granular materials: Discrete simulation of plane shear flows, Physical Review E, 72, no. 2, 021309 (2005).

[38] G. Lois, A. Lemaître, and J. M. Carlson, Numerical tests of constitutive laws for dense granular flows, Physical Review E,
72, 051303 (2005).

[39] N. Estrada, A. Taboada, and F. Radjai, *Shear strength and force transmission in granular media with rolling resistance*, Physical Review E, **78**, no. 2, 021301 (2008).

[40] T. Man, P. Zhang, Z. Ge, S. A. Galindo-Torres, and K. M. Hill, *Friction-dependent rheology of dry granular systems*, Acta Mechanica Sinica, **39**, no. 1, 722191 (2023).

[41] T. Hatano, *Critical scaling of granular rheology*, Progress of Theoretical Physics Supplement, **184**, 143-152 (2010).

[42] S. Chialvo, J. Sun, and S. Sundaresan, *Bridging the rheology of granular flows in three regimes*, Physical Review E, **85**, no. 2, 021305 (2012).

[43] F. Radjai, M. Jean, J. J. Moreau, and S. Roux, *Force distributions in dense two-dimensional granular systems*, Physical Review Letters, **77**, no. 2, 274 (1996).

[44] D. M. Mueth, H. M. Jaeger, and S. R. Nagel, *Force distribution in a granular medium*, Physical Review E, **57**, no. 3, 3164 (1998).

[45] F. Radjai, S. Roux, and J. J. Moreau, *Contact forces in a granular packing*, Chaos, **9**, no. 3, 544-550 (1999).

[46] M. J. Hertaeg, S. M. Fielding, and D. Bi, *Discontinuous shear thickening in biological tissue rheology*, Physical Review X, **14**, no. 1, 011027 (2024).

[47] D. Hernandez-Delfín, A. García, and M. Ellero, *Rheology of a crowd: from faster-is-slower to shear thickening*, Communications Physics, **7**, no. 1, 152 (2024).

[48] A. A. Tracton, *Coatings technology handbook*, CRC press (2005)

[49] F. Toussaint, C. Roy and P-H. Jézéquel Reducing shear thickening of cement-based suspensions, Rheologica acta, **48**, no 8, 883-895 (2009)