# Micromechanism Underlying Nonlinear Stress-DependentK

Authors: Xiang-yu Shang; Chen Yang; Guo-qing Zhou; Xiu-zhong Zheng
Journal: Advances in Materials Science and Engineering
 (2015)

---

## Abstract

In order to investigate the mechanism underlying the reported nonlinear at-rest coefficient of earth pressure,K
0 of clays at high pressure, a particle-scale model which can be used to calculate vertical and horizontal repulsion between clay particles has been proposed. This model has two initial states which represent the clays at low pressure and high pressure, and the particles in this model can undergo rotation and vertical translation. The computation shows that the majority of particles in a clay sample at high pressure state would experience rotation during one-dimensional compression. In addition, rotation of particles which tends to form a parallel structure causes an increase of the horizontal interparticle force, while vertical translation leads to a decrease in it. Finally, the link between interparticle force, microstructure, and macroscopic K
0 is analyzed and it can be used to interpret well the nonlinear changes in K
0 with both vertical consolidation stress and height-diameter ratio.


---

## Body

## 1. Introduction

The coefficient of earth pressure at rest,K
0, defined as the ratio of horizontal effective stress to vertical effective stress under the condition of zero horizontal deformation representing the in situ stress state of the ground, is a fundamental parameter in the analysis and design of geotechnical structures.Numerous studies in the past have indicated thatK
0 of a given soil is a constant depending on its strength parameter [1]. However, accumulated evidence over recent two decades demonstrated that K
0 is not necessarily a constant but generally a function of void ratio, stress level, and critical state friction angle even for given normally consolidated clay [2–5]. In particular, K
0 of clay increases nonlinearly with consolidation stresses over a wide range of pressures [6–9].In addition, previous investigation has indicated that, during one-dimensional compression, the distance between clay particles decreases continuously and the orientations of clay particles tend to be parallel to each other with their normal line pointing to the vertical direction [3, 10, 11].It is well known that there exist noncontact forces such as repulsion between clay particles due to the electric charge on the surface of clay particle, and these interparticle forces usually dominate the mechanical behaviors of clayey material. In addition, interparticle forces, which are balanced with macroscopic stress in the clays, determine the arrangements and orientations (i.e., microscopic structure) of clay particles. The macroscopic stress and deformation can be readily measured during mechanical tests on the clay specimens, but the information related to microscopic structure such as pore size, arrangement, and orientation of clay particles only can be analyzed after stopping the mechanical tests. Since the evolution of microscopic structure of clay during tests is usually unknown, it is difficult to establish the relation between macroscopic mechanical behavior and microscopic structure. Moreover, it is almost impossible to measure the interaction forces between clay particles during macroscopic tests. Therefore, the studies on the links between microscopic structure, interparticle forces, and macroscopic behavior, which is important for the thorough understanding of the intrinsic mechanism of macroscopic mechanical properties of clayey material, are seldom reported. It is not surprising that the micromechanism relating above nonlinearK
0 of clay to its microscopic structure during high pressure one-dimensional compression has not been exploited.This study aims to reveal the links between microscopic structure, interparticle forces, and macroscopicK
0 using numerical method and get insight into the micromechanism underlying the reported nonlinear K
0 of clay at high pressure.


## 2. Numerical Study on the Links between Interparticle Forces and Microscopic Structure

### 2.1. Repulsive Forces between Clay Particles

Since it has been revealed that the calculated relations between void ratio and vertical pressure, based on the double layer repulsive forces between clay particles, agree well with the measured results of both low and high pressure one-dimensional compression tests [12, 13], the noncontact forces between clay particles in this study are limited to the repulsive forces. The repulsion consists of osmotic pressure Π and electrical stress T
e. Π and T at arbitrary positon x in a system of clay-electrolyte are given as follows [14]:(1)
Π
x
=
R
T
∑
i
=
1
N
c
i
x
(2)
T
e
x
=
1
2
ε
E
x
·
E
x
I
-
ε
E
x
E
x
,where x is the position vector, c
i is the concentration of ith ion in the electrolyte, R
=
8.3114 J/(mol·K) is the universal gas constant, T is the absolute temperature, ε is the permittivity of the electrolyte, and E is electric field intensity. The repulsive pressure between two clay particles is(3)
F
=
Π
x
1
-
Π
x
2
-
T
e
x
1
-
T
e
x
2
,where x
2 is placed at any position between two particles and x
1 is an arbitrary position outside the region between these two particles.Bothc
i in (1) and E in (2) are related to the electric potential:(4)
c
i
=
c
0
cosh
⁡
v
e
ψ
k
T
(5)
E
=
-
∇
ψin which c
0 is the concentration of background electrolyte, e is the electronic unit charge, v is the valence of ion in the electrolyte, k
=
1.3806505
×
10
-
23 J/K is Boltzmann’s constant, ∇ denotes the gradient operator, and ψ is the electric potential. It should be noted that (4) only applies to binary monovalent electrolytes which is the common case in the study of saturated clay.In order to calculate the repulsion between clay particles, it is necessary to obtain the potential distribution in the clay-electrolyte system. Dimensionless potentialϕ can be determined by solving the following well-known Poisson-Boltzmann equation:(6)
∂
2
ϕ
∂
ζ
2
+
∂
2
ϕ
∂
η
2
=
sinh
⁡
ϕ
(7)
ϕ
=
v
e
ψ
k
T
(8)
ζ
=
K
x
(9)
η
=
K
y
,where(10)
K
=
2
n
v
e
2
ϵ
k
Tis the reciprocal of double-layer thickness in which n
=
c
0
×
A
v
o
g
a
d
r
o’s number.Combined with the following boundary conditions, the particular solution of (7) in a given domain Ω can be determined. Consider(11)
a
ϕ
+
h
∂
ϕ
∂
n
=
b
on
Γ
1
ϕ
=
g
on
Γ
2in which Γ
=
Γ
1
+
Γ
2 is the boundary enclosing the domain Ω and a, b, g, and h are constant parameters used to determine the boundary conditions.There are few analytical solutions of (6) due to its strong nonlinearity, and it is necessary to solve it using numerical method. Using a partial differential equation solver FreeFem++ which is open finite element software [14], (6) can be solved numerically. The first step to solve (6) in this software is to derive its weak form as follows: (12)
J
ϕ
,
w
=
∬
∂
w
∂
ζ
∂
ϕ
∂
ζ
+
∂
w
∂
η
∂
ϕ
∂
η
d
Ω
+
∫
w
a
ϕ
-
w
b
d
Γ
1
+
∫
w
sinh
⁡
ϕ
d
Γ
=
0
,in which w is an arbitrary weight function which is forced to be zero on the boundary Γ
2. Equation (12) is nonlinear and can be solved using iterative method. Assuming that there is only Dirichlet boundary condition, the following linear iterative equation based on Newton method can be obtained:(13)
J
ϕ
i
,
w
+
∬
∂
w
∂
ζ
∂
δ
ϕ
i
∂
ζ
+
∂
w
∂
η
∂
δ
ϕ
i
∂
η
d
Ω
+
∫
w
δ
ϕ
i
cosh
⁡
ϕ
i
d
Γ
=
0
,
ϕ
i
+
1
=
ϕ
i
+
δ
ϕ
i
,where the superscript i denotes the number of iteration steps.Figure1 presents the finite element mesh used for computation of potential around two inclined charged plates with finite length and the corresponding potential nephogram. Substituting the calculated potential into (4), the repulsion between clay particles can be solved.Figure 1
Finite element mesh and corresponding potential nephogram.


### 2.2. Numerical Model for Analysis

Actual microstructure of clay is complex and the amount of particles in a clay sample is huge as shown schematically in Figure2(a). Because the very fine mesh near the particle surface as shown in Figure 1 is required to accurately calculate the electrical repulsion, it is impractical to simulate actual saturated clay samples. The numerical model shown in Figure 2(b) can be seen as a local point in the clay sample. This particle-scale model consists of four clay particles. Each particle can feel both the vertical and horizontal forces from the neighbor particles. At the same time, the effects of macroscopic vertical compression and lateral constraint on the motion of particles can be considered in this model.Figure 2
Diagram of numerical model.
(a)
(b)The particle is 100 nm long and 1 nm thick and the potential on the surface of particle is 248.2 mV, which are the properties of typical montmorillonite clays [15]. The pore space is occupied by binary monovalent electrolyte solution, and the background electrolyte concentration is 1 mol/m3. According to (7)–(10), the corresponding dimensionless length of clay particle and surface potential are 10 and 9.94, respectively.Two kinds of initial states of the numerical model will be used in this study. The first one represents the clays consolidated at low pressure, and the second is corresponding to the compacted clays at high pressure. According to the reported microscopic structure of clays during one-dimensional compression [3, 11], edge-to-edge distance dH1 and face-to-face distance dV1 between clay particles shown in Figure 2(b) at high pressure are smaller than those at low pressure and the orientation of clay particle θ shown in Figure 2(b) at high pressure is also smaller than that at low pressure. θ corresponding to high and low pressure is assumed to be 10° and 30°, respectively. Based on the following relation between pore size and void ratio [16], dH1 and dV1 can be determined. Consider(14)
t
=
2
e
S
d
s
,where S is the specific surface area, d
s is the relative density for deep clay, t is the average distance between two neighbor parallel clay plates which represents the pore size, and e is void ratio. It is assumed that S
=
1
2
0 m2/g and d
s
=
2.7 g/m3 [13]. Initial void ratios related to high and low pressure are 0.49 and 0.97, respectively, and the calculated initial distances between clay particles using (14) are 3 nm and 6 nm.The basic models of motion of rigid clay particles are translation and rotation. Due to the effects of macroscopic vertical compression and lateral constraint on the motion of particles, the initial configuration of numerical model shown in Figure2(b) can change into the two new configurations. The first one is shown by the dotted lines in Figure 3(a), which results from the rotation of particles. The second as shown by the dotted lines in Figure 3(b) results from the vertical translation of particles.Figure 3
Possible configurations of numerical model.
(a)
(b)


### 2.3. Numerical Results

In order to investigate the microscopic mechanism underlying the difference inK
0 of clays at high and low pressure, it is necessary to study the incremental repulsion between particles resulting from the same amount of vertical displacement of clays at both low and high pressure during one-dimensional compression. Table 1 presents such numerical results of both vertical and horizontal interparticle repulsion.Table 1
Changes in repulsion between clay particles due to the change of configuration.
Initial state
Repulsion
Configuration
After rotation
After translation
Low pressure
Δ
F
x
+45.95
−4.91
Δ
F
y
+0.30
+0.38
High pressure
Δ
F
x
+37.37
−26.97
Δ
F
y
+1.51
+3.89It can be seen from Table1 that, experiencing the same amount of vertical compression in the case of low pressure, the required increase in vertical external forces for pure rotation is 0.3 which is comparable with that for pure translation (0.38). However, in the case of high pressure the required increase in vertical external forces for pure rotation (1.51) is much smaller than that for pure translation (3.89). Therefore, in the case of low pressure, the amount of rotating particles in clay specimens is comparable with that of pure translation. But, in the case of high pressure, the amount of rotation should be larger than that of pure translation.Another observation from Table1 is that pure rotation leads to an increase in horizontal interparticle force, but pure translation leads to a decrease in horizontal interparticle force. This holds true for low and high pressure. According to the definition of K
0, it can be deduced that pure rotation would cause an increase in K
0, while pure translation would lead to a decrease in K
0. Apart from the numerical simulations shown above, we carried out a great number of calculations with variables dV1, dH1, and θ, and similar observations have been obtained.


## 2.1. Repulsive Forces between Clay Particles

Since it has been revealed that the calculated relations between void ratio and vertical pressure, based on the double layer repulsive forces between clay particles, agree well with the measured results of both low and high pressure one-dimensional compression tests [12, 13], the noncontact forces between clay particles in this study are limited to the repulsive forces. The repulsion consists of osmotic pressure Π and electrical stress T
e. Π and T at arbitrary positon x in a system of clay-electrolyte are given as follows [14]:(1)
Π
x
=
R
T
∑
i
=
1
N
c
i
x
(2)
T
e
x
=
1
2
ε
E
x
·
E
x
I
-
ε
E
x
E
x
,where x is the position vector, c
i is the concentration of ith ion in the electrolyte, R
=
8.3114 J/(mol·K) is the universal gas constant, T is the absolute temperature, ε is the permittivity of the electrolyte, and E is electric field intensity. The repulsive pressure between two clay particles is(3)
F
=
Π
x
1
-
Π
x
2
-
T
e
x
1
-
T
e
x
2
,where x
2 is placed at any position between two particles and x
1 is an arbitrary position outside the region between these two particles.Bothc
i in (1) and E in (2) are related to the electric potential:(4)
c
i
=
c
0
cosh
⁡
v
e
ψ
k
T
(5)
E
=
-
∇
ψin which c
0 is the concentration of background electrolyte, e is the electronic unit charge, v is the valence of ion in the electrolyte, k
=
1.3806505
×
10
-
23 J/K is Boltzmann’s constant, ∇ denotes the gradient operator, and ψ is the electric potential. It should be noted that (4) only applies to binary monovalent electrolytes which is the common case in the study of saturated clay.In order to calculate the repulsion between clay particles, it is necessary to obtain the potential distribution in the clay-electrolyte system. Dimensionless potentialϕ can be determined by solving the following well-known Poisson-Boltzmann equation:(6)
∂
2
ϕ
∂
ζ
2
+
∂
2
ϕ
∂
η
2
=
sinh
⁡
ϕ
(7)
ϕ
=
v
e
ψ
k
T
(8)
ζ
=
K
x
(9)
η
=
K
y
,where(10)
K
=
2
n
v
e
2
ϵ
k
Tis the reciprocal of double-layer thickness in which n
=
c
0
×
A
v
o
g
a
d
r
o’s number.Combined with the following boundary conditions, the particular solution of (7) in a given domain Ω can be determined. Consider(11)
a
ϕ
+
h
∂
ϕ
∂
n
=
b
on
Γ
1
ϕ
=
g
on
Γ
2in which Γ
=
Γ
1
+
Γ
2 is the boundary enclosing the domain Ω and a, b, g, and h are constant parameters used to determine the boundary conditions.There are few analytical solutions of (6) due to its strong nonlinearity, and it is necessary to solve it using numerical method. Using a partial differential equation solver FreeFem++ which is open finite element software [14], (6) can be solved numerically. The first step to solve (6) in this software is to derive its weak form as follows: (12)
J
ϕ
,
w
=
∬
∂
w
∂
ζ
∂
ϕ
∂
ζ
+
∂
w
∂
η
∂
ϕ
∂
η
d
Ω
+
∫
w
a
ϕ
-
w
b
d
Γ
1
+
∫
w
sinh
⁡
ϕ
d
Γ
=
0
,in which w is an arbitrary weight function which is forced to be zero on the boundary Γ
2. Equation (12) is nonlinear and can be solved using iterative method. Assuming that there is only Dirichlet boundary condition, the following linear iterative equation based on Newton method can be obtained:(13)
J
ϕ
i
,
w
+
∬
∂
w
∂
ζ
∂
δ
ϕ
i
∂
ζ
+
∂
w
∂
η
∂
δ
ϕ
i
∂
η
d
Ω
+
∫
w
δ
ϕ
i
cosh
⁡
ϕ
i
d
Γ
=
0
,
ϕ
i
+
1
=
ϕ
i
+
δ
ϕ
i
,where the superscript i denotes the number of iteration steps.Figure1 presents the finite element mesh used for computation of potential around two inclined charged plates with finite length and the corresponding potential nephogram. Substituting the calculated potential into (4), the repulsion between clay particles can be solved.Figure 1
Finite element mesh and corresponding potential nephogram.


## 2.2. Numerical Model for Analysis

Actual microstructure of clay is complex and the amount of particles in a clay sample is huge as shown schematically in Figure2(a). Because the very fine mesh near the particle surface as shown in Figure 1 is required to accurately calculate the electrical repulsion, it is impractical to simulate actual saturated clay samples. The numerical model shown in Figure 2(b) can be seen as a local point in the clay sample. This particle-scale model consists of four clay particles. Each particle can feel both the vertical and horizontal forces from the neighbor particles. At the same time, the effects of macroscopic vertical compression and lateral constraint on the motion of particles can be considered in this model.Figure 2
Diagram of numerical model.
(a)
(b)The particle is 100 nm long and 1 nm thick and the potential on the surface of particle is 248.2 mV, which are the properties of typical montmorillonite clays [15]. The pore space is occupied by binary monovalent electrolyte solution, and the background electrolyte concentration is 1 mol/m3. According to (7)–(10), the corresponding dimensionless length of clay particle and surface potential are 10 and 9.94, respectively.Two kinds of initial states of the numerical model will be used in this study. The first one represents the clays consolidated at low pressure, and the second is corresponding to the compacted clays at high pressure. According to the reported microscopic structure of clays during one-dimensional compression [3, 11], edge-to-edge distance dH1 and face-to-face distance dV1 between clay particles shown in Figure 2(b) at high pressure are smaller than those at low pressure and the orientation of clay particle θ shown in Figure 2(b) at high pressure is also smaller than that at low pressure. θ corresponding to high and low pressure is assumed to be 10° and 30°, respectively. Based on the following relation between pore size and void ratio [16], dH1 and dV1 can be determined. Consider(14)
t
=
2
e
S
d
s
,where S is the specific surface area, d
s is the relative density for deep clay, t is the average distance between two neighbor parallel clay plates which represents the pore size, and e is void ratio. It is assumed that S
=
1
2
0 m2/g and d
s
=
2.7 g/m3 [13]. Initial void ratios related to high and low pressure are 0.49 and 0.97, respectively, and the calculated initial distances between clay particles using (14) are 3 nm and 6 nm.The basic models of motion of rigid clay particles are translation and rotation. Due to the effects of macroscopic vertical compression and lateral constraint on the motion of particles, the initial configuration of numerical model shown in Figure2(b) can change into the two new configurations. The first one is shown by the dotted lines in Figure 3(a), which results from the rotation of particles. The second as shown by the dotted lines in Figure 3(b) results from the vertical translation of particles.Figure 3
Possible configurations of numerical model.
(a)
(b)


## 2.3. Numerical Results

In order to investigate the microscopic mechanism underlying the difference inK
0 of clays at high and low pressure, it is necessary to study the incremental repulsion between particles resulting from the same amount of vertical displacement of clays at both low and high pressure during one-dimensional compression. Table 1 presents such numerical results of both vertical and horizontal interparticle repulsion.Table 1
Changes in repulsion between clay particles due to the change of configuration.
Initial state
Repulsion
Configuration
After rotation
After translation
Low pressure
Δ
F
x
+45.95
−4.91
Δ
F
y
+0.30
+0.38
High pressure
Δ
F
x
+37.37
−26.97
Δ
F
y
+1.51
+3.89It can be seen from Table1 that, experiencing the same amount of vertical compression in the case of low pressure, the required increase in vertical external forces for pure rotation is 0.3 which is comparable with that for pure translation (0.38). However, in the case of high pressure the required increase in vertical external forces for pure rotation (1.51) is much smaller than that for pure translation (3.89). Therefore, in the case of low pressure, the amount of rotating particles in clay specimens is comparable with that of pure translation. But, in the case of high pressure, the amount of rotation should be larger than that of pure translation.Another observation from Table1 is that pure rotation leads to an increase in horizontal interparticle force, but pure translation leads to a decrease in horizontal interparticle force. This holds true for low and high pressure. According to the definition of K
0, it can be deduced that pure rotation would cause an increase in K
0, while pure translation would lead to a decrease in K
0. Apart from the numerical simulations shown above, we carried out a great number of calculations with variables dV1, dH1, and θ, and similar observations have been obtained.


## 3. Micromechanism Underlying NonlinearK
0 of Clays over a Wide Range of Pressures

According to the numerical results described in Section2.3, rotation and translation would cause an increase and a decrease in K
0, respectively, and in the case of low pressure, the amount of clay particles rotating during one-dimensional compression is comparable with that of translating particles. Accordingly, the increase in K
0 caused by rotation is almost counterbalanced by the decrease in K
0 caused by translation. Therefore, the measured K
0 during low pressure one-dimensional compression of clays largely remains unchanged.However, in the case of the high pressure, the amount of rotating particles is larger than that of translating particles, and so the increase inK
0 caused by rotation is more significant than the decrease in K
0 caused by translation. Consequently, the measured K
0 should increase with the vertical pressure during high pressure one-dimensional compression of clays. This is consistent with the reported test results shown in Figure 4.Figure 4
NonlinearK
0 of reported clays at a wide range of pressures.Apart from the nonlinearK
0 of clays as shown in Figure 4, it has been observed by Min [9] using triaxial testing system that K
0 of deep clay in East China at high pressure increases with the decreasing ratio of height to diameter of clay specimens as shown in Figure 5. This test result can also be interpreted well based on the above analysis.Figure 5
K
0 at various pressures versus h
/
d (min, 2010).Although it was usually thought that the ideal one-dimensional compression tests do not have the problem of end friction, the small lateral deformation of specimens during such tests using triaxial apparatus and so a nonignorable end friction is unavoidable. If there were not such end frictions,K
0 of clay specimens in these tests with different height-diameter ratio would be the same.The end friction will exert a horizontal shear stress on the ends of the clay specimens as shown in Figure6 and induces shear stresses inside the specimens. The shear stresses which would produce a torque on clay particles may make it easier to rotate the particles. Therefore, in contrast to those only subjected to vertical stress, the inclined clay particles subjected to both vertical stress and shear stress are more inclined to rotate themselves. According to the early analysis in this section, the increase in K
0 during one-dimensional compression test with end friction will be more significant than that without end friction. As a result, a higher K
0 was observed in the tests with end restraints.Figure 6
Schematic of the effect of end friction.


## 4. Conclusions

In order to explore the micromechanism underlying the reported nonlinearK
0 of clays at high pressure, a particle-scale numerical model is established and the analysis on the links between particle-scale interaction forces, structure, and macroscopic K
0 is carried out. The main conclusions of this study are as follows.(1)
In the case of low pressure, the amount of purely rotating particles in a clay sample is comparable to that of purely translating particles. But the number of purely rotating particles is larger than that of purely translating particles in the case of high pressure.
(2)
Pure rotation which tends to make clays reach a parallel configuration would cause an increase inK
0 while pure translation which causes smaller face-to-face distance between clay particles would lead to a decrease in K
0.
(3)
Based on the above observations, the links between interparticle forces, microstructure, and macroscopicK
0 have been investigated. The links can be used to well interpret the nonlinear changes in K
0 of clays at high pressure with both vertical pressure and height-diameter ratio of clay specimens.



---

*Source: 540621-2015-10-04.xml*