close

Вход

Забыли?

вход по аккаунту

?

j.ijplas.2018.08.007

код для вставкиСкачать
Accepted Manuscript
A discrete element modelling approach for fatigue damage growth in cemented
materials
Nhu H.T. Nguyen, Ha H. Bui, J. Kodikara, S. Arooran, F. Darve
PII:
S0749-6419(18)30246-8
DOI:
10.1016/j.ijplas.2018.08.007
Reference:
INTPLA 2397
To appear in:
International Journal of Plasticity
Received Date: 20 April 2018
Revised Date:
14 August 2018
Accepted Date: 14 August 2018
Please cite this article as: Nguyen, N.H.T., Bui, H.H., Kodikara, J., Arooran, S., Darve, F., A discrete
element modelling approach for fatigue damage growth in cemented materials, International Journal of
Plasticity (2018), doi: 10.1016/j.ijplas.2018.08.007.
This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to
our customers we are providing this early version of the manuscript. The manuscript will undergo
copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please
note that during the production process errors may be discovered which could affect the content, and all
legal disclaimers that apply to the journal pertain.
ACCEPTED MANUSCRIPT
A discrete element modelling approach for fatigue damage
growth in cemented materials
RI
PT
AUTHORS:
AUTHOR’S AFFILIATIONS
M
AN
U
SC
Nhu H.T. Nguyen1, Ha H. Bui1*, J. Kodikara 1, S. Arooran1 and F. Darve2
Department of Civil Engineering, Monash University, Melbourne, Australia
2
Grenoble CNRS, UMR 5521, 3SR, Grenoble Alpes University, Grenoble, France
EP
TE
D
1
AC
C
*CORRESPONDING AUTHOR
Dr. Ha H. Bui, Senior Lecturer
Department of Civil Engineering, Monash University
Room 126, Building 60, Monash University, Victoria 3800, Australia
Email: ha.bui@monash.edu Tel : +61 3 9905 2599
1
ACCEPTED MANUSCRIPT
Abstract
This study proposes a modelling approach capable of capturing both the fatigue response and
localisation of failure in cemented materials. This modelling approach takes advantage of
Discrete Element Method (DEM) to reproduce the heterogeneous microstructure and crack
development in cemented materials. In conjunction with this, a new constitutive model is
RI
PT
developed to characterise the fatigue behaviour of the materials at the grain scale. The model
formulation is based on coupling damage mechanics and plasticity theory and combining
with a fatigue damage evolution law to describe the degrading response of cemented
materials subjected to cyclic loading. The proposed model is employed to govern the explicit
SC
behaviours of DEM bonding contacts representing cement bridges between aggregates in the
physical materials. The macro-behaviour is then obtained in DEM simulations as the
M
AN
U
collective response of all contacts and particles in the material domain. Through numerical
experiments, the proposed modelling approach is shown to capture well the fatigue behaviour
and cracking process in cemented materials subjected to cyclic loading. The microstructural
effects on the fatigue response of the materials are naturally reproduced in simulations thanks
to the discrete nature of DEM. These results demonstrate the capability of the proposed
modelling approach as well as its potential to be a faithfully numerical technique for
TE
D
modelling and investigating fatigue damage in cemented materials.
Keywords: Discrete element method, fatigue damage, cemented materials, couple damage
AC
C
EP
plasticity, fatigue model
2
ACCEPTED MANUSCRIPT
Introduction
Fatigue cracking due to repetitive loads has been recognised as one of the main distress
modes in many structures made by cemented materials such as highway pavements, airport
pavements and concrete foundations of wind power plans (Dattoma et al., 2006; Fathima &
Kishen, 2013). When these structures are subjected to cyclic loading, plastic deformation is
RI
PT
observed to develop even though the applied nominal stress is below the material yielding
limit. If the applied nominal stress is sufficiently small, plastic deformation gradually
stabilises and eventually unchanged after a number of load cycles. This behaviour is called
shakedown (Weichert & Maier, 2014; Zhang et al., 2017), and the magnitude of cyclic loads
SC
under which the plastic deformation can stabilise is the shakedown limit or endurance limit of
the materials. Nevertheless, if the applied stress is greater than the shakedown limit, plastic
damage of cemented structures.
M
AN
U
deformation keeps accelerating with the increase of cyclic loading, leading to the fatigue
At the microstructural scale, fatigue damage in cemented materials is considered as a
progressive process of microstructural changes with increasing the number of load cycles.
This process is associated with the propagation and coalescence of micro-cracks in the
materials (Liang & Zhou, 1997; Fathima & Kishen, 2013). The development of micro-cracks
TE
D
in cemented materials subjected to cyclic loading is mainly attributed to their heterogeneous
microstructure, which consists of aggregates and cement mortar with pre-existing microcracks at the mortar matrix or at the interfaces between mortar and aggregates (Zhang, 2001;
Ray & Kishen, 2010). This heterogeneous feature of cemented materials causes stress to
EP
concentrate in areas near micro-crack tips, making local stresses in these areas to be higher
than the applied nominal stress. As a consequence, the micro-cracks can grow although the
applied nominal stress has not reached the material yielding limit. In this context, the fatigue
AC
C
1.
behaviour of cemented materials can be considered to depend on the rate of micro-cracking
per load cycle and material microstructural features. Modelling of fatigue damage in
cemented materials thus requires appropriate descriptions of the fatigue crack growth rate and
microstructural effects.
Existing modelling approach for fatigue damage in cemented materials can be classified into
two main methods, based on fracture mechanics or damage mechanics. The fracture
mechanics-based method describes the fatigue behaviour of cemented materials as the
relationship between the crack growth rate per load cycle and the fracture mechanics terms
such as stress intensity factor and fracture energy (Ray & Kishen, 2010; Simon & Kishen,
3
ACCEPTED MANUSCRIPT
2017), or based on the cohesive fracture theory (Nguyen et al., 2001; Yang et al., 2001). On
the other hand, the damage mechanics-based method focuses on developing the evolution law
of damage variable with respect to load cycles (Papa & Taliercio, 1996; Alliche, 2004). The
evolution of damage variable, in this sense, represents the loss of material load-carrying
capacity when undergoing cyclic loading. As the fatigue damage of cemented materials is
RI
PT
considered to be an irreversible process (Horii et al., 1992; Fathima & Kishen, 2013), the
theory of plasticity has also been adopted to describe the increment of plastic deformation in
the material fatigue response (Voyiadjis & Abu-Lebdeh, 1994; Ristinmaa, 1995; Moreo et al.,
2007; Al-Rub & Darabi, 2012; Darabi et al., 2012; Carrara & De Lorenzis, 2015; Cheng et al.,
SC
2016). Although a great effort has been dedicated to developing constitutive models for the
fatigue damage of cemented materials, current modelling approaches still face challenges in
capturing both the fatigue response and localised failure in the materials subjected to cyclic
M
AN
U
loading. This can be attributed to the fact that modelling fatigue damage is more complicated
than other types of damage (Xiao et al., 1998), as it requires the appropriate characterisations
of different stages of crack growth (i.e. stable crack growth and unstable crack growth (Gao
& Hsu, 1998)) as well as the material heterogeneous microstructure. The microstructural
effects on the macro-behaviour of cemented materials can be addressed by micromechanics-
TE
D
based models (Krajcinovic & Fanella, 1986; Mohamed & Hansen, 1999; Shen et al., 2012;
Shen & Shao, 2016). Besides, most existing fatigue models have been developed for
continuum modelling and implemented into continuum computational methods such as finite
element method (FEM) for numerical simulations. The continuum methods, nonetheless,
EP
suffer from over-complicated model formulation and ad-hoc treatments to deal with the
heterogeneous response and localised failure occurred in the fatigue damage process of
AC
C
cemented materials. Therefore, it is necessary to have an alternative approach that can
effectively handle the heterogeneous response and localised failure in cemented materials,
hence easing the complexity of constitutive model formulation.
In recent decades, Discrete Element Method (DEM) (Cundall & Strack, 1979) has emerged
as an effective numerical tool to model the damage behaviour of cemented materials and
reproduce the material microstructural responses. Different from continuum methods, which
consider materials as a single continuous body, DEM discretises material domains as an
assemblage of particles, carrying their own physical properties and motions. These particles
interact with each other through contacts and these interactions are governed by contact
models. This discrete nature enables DEM to naturally replicate the heterogeneous
4
ACCEPTED MANUSCRIPT
microstructure of cemented materials. Furthermore, as DEM allows large displacement and
full detachment between particles, this method can handle the localised failure and
discontinuity in the materials with no need for ad-hoc treatments on constitutive contact
models. Such advanced features make DEM a potential method to model the fatigue damage
of cemented materials. In fact, several bonding contact models have been developed to
RI
PT
characterise the damage behaviour of cemented materials under monotonic loadings (Scholtès
& Donzé, 2013; Gui et al., 2015; Gui et al. 2016; Oñate et al., 2015; Nguyen et al., 2017a;
Nguyen et al., 2017c). As the fatigue damage is due to the micro-cracking induced from the
heterogeneous microstructure of the materials, these bonding contact models might be used to
SC
characterise fatigue damage in cemented materials if the DEM particles and contacts are
modelled at the molecular scale and the true microstructure of cemented materials is
accurately replicated. Nonetheless, this modelling approach would remarkably increase the
M
AN
U
overall computational cost as millions of particles will be generated in simulations. For more
computational efficiency, particles in DEM simulations are generally modelled at larger
scales (e.g. grain scale or aggregate scale). At these scales, the above bonding contact models
are unable to describe the fatigue damage growth in cemented materials owing to the lack of
intrinsic fatigue damage mechanisms in the cement bridges. As a result, it is necessary to
TE
D
have bonding contact models for DEM incorporating the grain-scale fatigue damage
mechanisms to capture the fatigue responses of cemented materials. However, to the authors’
knowledge, there is no existing contact model for DEM capable of capturing the fatigue
behaviour of cemented materials. This indeed prevents the application of DEM for modelling
EP
fatigue damage in the materials.
The aim of this paper is to propose a modelling approach based on DEM to characterise
AC
C
fatigue damage in cemented materials. In this approach, DEM plays the role in modelling the
heterogeneous response and localised failure in the materials, while a new fatigue bonding
model is developed to incorporate physical fatigue damage mechanisms at cement bridges in
the materials. The unique feature of this approach is its capability in representing the
microstructural processes of cemented materials and easily capturing crack initiation and
propagation. The new fatigue model remains simple as it is formulated in two principal
directions (i.e. normal and shear) of the DEM contact and only characterises the grain-scale
fatigue behaviour of cement bridges, while different stages of fatigue crack growth exhibited
at the macroscopic scale of cemented materials is captured as the collective response of all
particles and contacts in the material domain. The capability of this proposed modelling
5
ACCEPTED MANUSCRIPT
approach in capturing and predicting the fatigue behaviour of cemented materials is first
compared with the experimental data. Subsequently, its ability in reproducing the
microstructural effect on the fatigue damage response of cemented materials is demonstrated
in further simulations.
This article is organised as follows. The basis of DEM computation is briefly explained in
RI
PT
Section 2. Details of the model formulation and implementation algorithm are then presented
in Section 3. This is followed by the demonstrations of model behaviour at the constitutive
level and parametric studies in Section 4 in order to illustrate the influences of the model
parameters on contact behaviour at the constitutive level. Next, Section 5 presents the
SC
application of the proposed modelling approach to capture the fatigue damage responses and
microstructural behaviour of cemented materials. Finally, Section 6 summarises conclusions
Numerical background
Particle motion equations
The DEM is adopted in this study as the numerical framework to simulate the behaviour of
cemented materials. In DEM, the computational domain is discretised into a number of rigid
second law as follows:
−
=
where
+
EP
−
=
TE
D
particles in different sizes. These particles can move and rotate in accordance with Newton’s
is the translational acceleration;
AC
C
2.
M
AN
U
drawn from this study.
is the moment of inertia of the particle.
(1)
(2)
is the rotational acceleration;
and
is the mass and
are the total force and moment acting on the
particles as illustrated in Fig. 1 and are computed as the vectorial summation of all applied
and contact forces as follows:
=
=
where
+
(3)
× +
is the applied force;
(4)
is the contact force, which is obtained based on the force-
displacement relationship determined in contact models;
particle centre to the contact point and
is the vector connecting the
is the applied moment. On the other hand,
6
and
ACCEPTED MANUSCRIPT
are the damping force and moment, incorporated into Eqs. (1) and (2) to damp out the
particle accelerations. This is essential to setup a quasi-static condition in simulations as the
DEM formulations are fully dynamics. In this study, the damping force and moment are
computed using the global damping formulations proposed in Potyondy & Cundall, (2004).
The global damping force and moment are proportional to the total force and moment acting
on the particle, but imposed in opposite directions as follows:
+
=
RI
PT
=α
(5)
(6)
where α is the global damping coefficient. To achieve a quasi-static condition in DEM
SC
simulations, Potyondy & Cundall, (2004) suggested using α = 0.7 so that particle’s
M
AN
U
acceleration is heavily damped. In addition, the loading rate in simulations should be chosen
small enough to ensure numerical results at the structural level are not affected by changes in
loading rate. In this sense, the effect of kinetic energy on the macro-behaviour in DEM
simulations can be negligible.
Subsequently, the motion equations of particles are integrated in time. In DEM computation,
the standard Leap-Frog algorithm is normally used to integrate the particle motion equations.
TE
D
In particular, once particle accelerations at a time step “ ” are known, the particle velocities
and displacements at the next time step are updated as follows:
!∆ /$
!∆
=
=
%∆ /$
+
+
%∆ /$
+
∆
∆
!∆ /$
=& +
∆
!∆ /$
(7)
(8)
(9)
∆
AC
C
&
!∆
=
EP
!∆ /$
(10)
As the DEM computation uses the explicit time integration scheme, it is necessary to impose
a constraint for the computational time step as below to ensure numerical stability:
∆ ≤∆
()
in which ∆
∆
()
= C,
(11)
()
is the critical time step determined as:
-.
⁄/
(12)
7
where C,
ACCEPTED MANUSCRIPT
-
≤ 1 is the stability coefficient and / is the stiffness parameter. In the modelling
of samples with numerous particles, the critical time step is normally attained as the
minimum of all critical time steps computed on all particles in the assembly.
DEM modelling of intact materials
DEM was conventionally developed to model the non-cohesive behaviour of granular media.
RI
PT
This method has later been applied to simulate the behaviour of intact materials with the
introduction of a bonding contact model (Potyondy & Cundall, 2004). In this model, particles
are allowed to be bonded when they are not touching each other, as long as the distance
between their centres is smaller than a distance tolerance. The bonded particles can interact
even when they are non-overlapped and the bonding contact can resist all tension, shear and
links the contact force to contact stress as follows:
M
AN
U
= 3 2 = 34 2 + 3, 2 SC
compression stresses. The bonding contact is assigned a specific shape and area 2 , which
(13)
where 34 and 3, are the normal and shear stress components acting in the bonding contact
(see Fig. 1), respectively. The normal and shear stresses in the bonding contact are induced
from the relative displacements of the bonded particles based on the bonding contact models.
TE
D
In DEM modelling of intact materials, the bonding contact model plays an important role in
characterising the local behaviour of materials while the macro-behaviour of materials is
AC
C
EP
captured as the combination of all contact and particle responses in numerical samples.
Fig. 1: The total force and moment acting on particles and the normal and shear stress vectors
acting on the bonding contact
3.
A fatigue damage contact model for DEM
8
ACCEPTED MANUSCRIPT
As mentioned earlier, DEM is able to reproduce the heterogeneity of cemented materials.
However, the uses of circular particles to represent aggregates and contact bonds to represent
mortar matrix in the general DEM simulation have significantly simplified the true
microstructural details of the materials (e.g. morphologies of aggregates and cohesive mortar
matrix). Furthermore, existing bonding contact models in DEM are unable to capture the
RI
PT
fatigue damage growth in cemented materials owing to the lack of fatigue crack growth
mechanisms at the grain scale. These issues have indeed limited the applicability of DEM for
modelling fatigue cracking in cemented materials. To overcome these limitations, a new
fatigue bonding contact model for cemented materials is developed in this study. This
SC
bonding contact model phenomenologically takes into account those microstructural details
of cemented materials missing in DEM simulations and more importantly captures the fatigue
miro-cracking mechanism of the materials. The proposed model is formulated on the basis of
M
AN
U
the cohesive damage-plastic model developed by the authors to model the grain-scale damage
behaviour of cemented materials under monotonic loads (Nguyen et al., 2017b). In the
proposed model, new fatigue variables and fatigue evolution law are incorporated into the
original cohesive damage-plastic law to account for the micro-crack growth and irreversible
displacement occurred in the cement bridge when subjected to cyclic loading. Details of the
3.1. Model formulation
TE
D
fatigue bonding contact model formulation are described as below.
The expression of the relative displacement between two particles connected by a bonding
contact in the original model is first redefined to incorporate the fatigue displacement
EP
variable. This variable represents the irreversible displacement at cement bridges caused by
fatigue cracking, which distinguishes itself from the plastic displacement developed when
AC
C
stress reaches the yielding limit. The separation of fatigue and plastic displacements is
necessary to reflect fundamental differences between the underlying mechanisms of these
two irreversible processes. Specifically, fatigue occurs due to the degradation of elastic
modulus produced by micro-cracking. The micro-cracking is an irreversible process
associated with the local micro-mechanisms of failure. Micro-cracks can grow even when the
macro-behaviour of the materials is quasi-static and cyclic stress state is still below the plastic
yield surface of the materials. On the other hand, once micro-cracks grow into macro-cracks,
the plastic displacement occurs as the result of frictional interaction between two sides of the
macro-cracks. At this stage, the stress state of the material lies on its plastic yield surface. As
9
ACCEPTED MANUSCRIPT
a consequence, the relative displacement in a bonding contact is now equal to the sum of
elastic, plastic and fatigue displacement components as follows:
=
5
+
6
+
7
(14)
Accordingly, the contact stress-displacement relationships in normal and shear directions are
84 = /49 1 − : ;<4 − <46 − <4 = − :/49 〈−;<4 − <46 − <4 =〉
7
7
8, = /,9 1 − : ;<, − <, − <, =
7
6
RI
PT
rewritten as:
(15)
(16)
SC
where /49 and /,9 are the initial normal and shear stiffness of the bonding contact; : is the
damage variable. The Macauley brackets embraced in Eq. (15) is to impose the full recovery
of contact normal stiffness in compression modes.
M
AN
U
The expression of damage variable : is also revised to take into account the contribution of
fatigue damage. In this model, it is considered that within a load cycle (i.e. loading and then
unloading), if the stress state in bonding contact is below its yielding limit, only fatigue
damage can develop in the contact. On the other hand, if stress in the bonding contact reaches
the yielding limit, fatigue damage is lumped into the plastic damage due to yielding. In this
sense, the damage variable is expressed as the sum of damage increments in all load cycles,
:=
A
)BC
@: =
A
)BC
D 1 − E) @:6 + E) @:7 F
(17)
EP
(@:7 ):
TE
D
each of which consists of an incremental damage due to either yielding (@:6 ) or fatigue
AC
C
where E) is a shift factor in each load cycle, which equals 0 if the stress state in the contact
reaches its yield surface and 1 if the stress state in the contact is below its yield surface within
a load cycle. It can be seen from Eqs. (15), (16) and (17), the normal and shear stresses in the
bonding contact are influenced by two sets of internal variables, including plastic variables
(<46 , <,6 and :6 ,) and fatigue variables (<4 , <, and :7 ). The plastic variables represent the
7
7
plastic damage state of the bonding contact, while the fatigue variables manifest its fatigue
damage state. The yielding and fatigue processes in the bonding contact are thus
characterised by the increments of these internal variables.
6
6
Increments of plastic variables: @<4 , @<, and @:6
10
ACCEPTED MANUSCRIPT
The increments of plastic variables are calculated based on the original damage-plasticity
model. Details of the original model formulations can be found in Nguyen et al., (2017b).
Herein, only the key formulations and modifications made to this model for the calculations
of incremental plastic displacements and damage variable are presented.
Based on the original damage-plasticity model, the plastic displacements and damage
RI
PT
variable can develop when stress in the bonding contact reaches a yield criterion which takes
the following form:
FH = I $ 8,$ − JK 9 1 − : − 84 L MN$ + JK 9 1 − : − 8 9 1 − : L MN$
$
SC
.2K 9 8 9 L M − 8 9 L M
where I =
K9
(18)
In the above equations, 8 9 is the contact tensile strength, K 9 is the contact cohesion and M is
M
AN
U
the friction angle. When the bonding contact is yielding, the increments of normal and shear
plastic displacements are updated based on the standard plastic flow rule as follows:
@<46 = 2@λ L ψJK 9 1 − : − 84 L ψN
@<, = 2@λ8,
6
(19)
(20)
TE
D
where ψ is the dilatancy angle and @λ is the plastic multiplier. The evolution of plastic
damage variable is then derived from the formulation of damage variable evolution proposed
in the original model, yet taking into account the reformulation of total displacement and
damage variables in Eqs. (14) and (17). This results in the following function of incremental
V
V
ST
ST
%R WU ! WX Y
TU
TX
AC
C
@:6 = χ P1 − Q
EP
plastic damage variable:
Z
(21)
In Eq. (21), χ is the term associated with the deformation state of the bonding contact and is
expressed as follows:
χ=Q
V
]
V
]
∑ STU
!∑ STU ∑ STX !∑ STX
^
%[
!
W
W
TU
TX
(22)
where <4 and <, are the softening parameters as in the original model. Eqs. (21) and (22) can
facilitate coupling between damage mechanics and plasticity theory to fully describe both
cracking and irreversible displacement occurred in the damage process of cemented materials.
In addition, the power law obtained for the incremental plastic damage variable can provide a
11
ACCEPTED MANUSCRIPT
smooth transition between the cohesive and non-cohesive behaviour of a contact bond when
damage variable is close to one (Nguyen et al., 2017b).
Increments of fatigue variables: @<4 , @<, and @:7
7
7
RI
PT
Different from the plastic variables, the fatigue variables are assumed to develop when the
maximum stress in the bonding contact within a load cycle satisfies a fatigue criterion, i.e.
_7 ` 0. In this study, the fatigue criterion is derived from the yield criterion by incorporating
a fatigue threshold parameter, a, as follows:
SC
Fb = I $ 8,$ − JaK 9 1 − : − 84 L MN$ + JaK 9 1 − : − a8 9 1 − : L MN$
(23)
From Eq. (23), it can be seen that a determines the proportion of the fatigue limit to the
M
AN
U
plastic yielding limit of the bonding contact. The value of a has to be chosen between 0 and
1 to ensure the fatigue criterion is always below the yield criterion, hence satisfying the predetermined condition that fatigue damage is able to evolve under the yield criterion. Different
from the yield criterion, the fatigue criterion only plays the role in determining the stress
states in which fatigue damage can evolve in bonding contacts, while the increment of fatigue
TE
D
damage variable is defined independently later on. The shape of the fatigue criterion is
displayed in Fig. 2 along with the yield criterion. The shaded area represents the range of
AC
C
EP
maximum cyclic stress state in the bonding contact wherein fatigue damage can develop.
Fig. 2: The yield and fatigue criteria of the proposed model. The shaded area represents the
range of maximum cyclic stress states wherein fatigue damage is allowed to develop
12
ACCEPTED MANUSCRIPT
Once the above fatigue conditions are met, the fatigue damage variable will increase.
Different from the plastic damage variables whose increments are calculated cycle-by-cycle,
the increment of the fatigue damage variable is only updated once at the onset of unloading
path in a load cycle. This updating method is considered to be more computationally efficient
than the cycle-by-cycle update, especially in the modelling of high-cycle fatigue damage
RI
PT
(Turon et al., 2007). Moreover, it sets a precedence for embedding the “cycle jump”
technique to reduce computational time in the simulations of cyclic tests. As a consequence
of this fatigue damage updating method, the incremental fatigue damage variable is
formulated as the function of the number of load cycles in accordance with Paris’ law (Paris
SC
& Erdogan, 1963). Originally, Paris’s law defined the fatigue crack growth rate per load
cycle as a power function of the stress intensity factor as follows:
e
M
AN
U
@L7
= K ∆d
@c
where L7 is the fatigue crack length; c is the number of load cycles; K and
(24)
are the material
fatigue parameters; and ∆d is the stress intensity factor range in a load cycle. The fatigue
crack length L7 is then linked to the fatigue damage variable based on the geometrical
:=
L
f
TE
D
expression of the damage variable for rectangular DEM contacts as follows (see Fig. 3):
(25)
(25) leads to
(26)
AC
C
@L7 = f@:7
EP
where L is the total crack length. In case the crack growth is due to fatigue, i.e. @L = @L7 , Eq.
Fig. 3: The shape of a bonding contact and presumed crack length
By substituting Eq. (26) into Eq. (24), the incremental fatigue damage variable is obtained as:
13
@:7 K ∆d
=
@c
f
ACCEPTED MANUSCRIPT
e
(27)
In Eq. (27), the stress intensity factor range, ∆d, represents the effect of cyclic stress state on
the fatigue behaviour of cement bonds. According to this equation, ∆d = ∆dg when the
bonding contact is under mode I loading, while ∆d = ∆dgg when the bonding contact is under
RI
PT
mode II loading. The contributions of both modes I and II are also taken into account in a
mixed-mode loading based on an energy release rate analysis (Paris & Sih, 1965). As a result,
the stress intensity factor is considered to take the following form:
∆d = . ∆dg
$
+ ∆dgg
$
(28)
SC
In Eq. (28), mode I and mode II stress intensity factors are calculated based on the fracture
mechanics theory (Irwin, 1957):
k
∆dgg = hgg √jL∆8, = hgg √jL;8,e
− 84e)4 =
k
M
AN
U
∆dg = hg √jL∆84 = hg √jL;84e
− 8,e)4 =
(29)
(30)
where hg and hgg are the geometrical factors calculated by Eq. (31). This equation is derived
by Tada et al., (1973) based on the analysis results of mode I and mode II tests on rectangular
plates with a horizontal centre notch, which is similar to the shape of the DEM contacts.
TE
D
L $
L p
1 − 0.025 m n + 0.06 m n
f
f
hg = hgg = h =
jL
qrs m n
2f
(31)
EP
Subsequently, substituting Eq. (25) into (29), (30) and (31) yields:
∆dg = h√jf:;84e
h=
k
− 84e)4 =
(32)
− 8,e)4 =
AC
C
∆dgg = h√jf:;8,e
k
(33)
1 − 0.025:$ + 0.06:p
qrs mj: n
2
(34)
In Eq. (32) and (33), the value of damage variable has to be greater than 0 so that Paris’ law
(i.e. Eq. (27)) can be valid. To satisfy this condition, the damage variable is assumed to take a
small initial value :)4) in numerical simulations. This assumption is reasonable in the sense
that :)4) represents the pre-existing micro-cracks in cement bonds after the sample generation
process.
14
ACCEPTED MANUSCRIPT
Along with the evolution of the fatigue damage variable, the fatigue displacements in normal
and shear directions increase. The incremental fatigue displacements are related to the
incremental fatigue damage variable by Eq. (35). This equation is derived in the same way as
the formulation of incremental plastic damage variable in Eq. (21):
]
Z
(35)
RI
PT
@:7 = χ P1 − Q
]
ST
ST
%[ WU ! WX ^
TU
TX
This leads to
@:7
@<4 @<,
+
= −t R1 −
Y
<4
<,
χ
7
(36)
SC
7
As Eq. (36) is insufficient to solve for the normal and shear incremental fatigue
M
AN
U
displacements in mixed-mode conditions, an additional constraint between these two
variables is required. In this work, the additional constraint is defined as the ratio between
@<4 and @<, . This constraint represents the proportion of incremental fatigue displacements
7
7
in normal and shear directions in a load cycle and can be derived following the concept of
plastic flow rule used to calculate the normal and shear plastic displacements in Eqs. (19) and
(20) as follows:
7
@<,
7
=
L ψJaK 9 1 − : − 84 L ψN
8,
TE
D
@<4
(37)
Eqs. (36) and (37) make a set of equations with two unknowns, @<4 and @<, . By solving this
7
7
EP
set of equations, the values of @<4 and @<, can be obtained. Accordingly, the plastic and
7
7
fatigue damage variables are calculated and the total damage variable is updated using Eq.
AC
C
(17). During the failure process, the total damage variable increases from :)4) to 1. In parallel
to this, both the yield and fatigue criterions are allowed to shrink in order to represent the loss
of load-carrying capacity of the bonding contact (see Fig. 4). The yield criterion shrinks faster
than the fatigue criterion, and these two criteria meet when the damage variable reaches one,
which corresponds to the final damage stage of the bonding contact.
15
RI
PT
ACCEPTED MANUSCRIPT
3.2. The model implementation
M
AN
U
damage variable
SC
Fig. 4: The evolutions of the yield and fatigue criteria with respect to the increment of
For DEM numerical experiments, the proposed constitutive model is implemented into PFC
5.0 (Itasca, 2014). The implementation algorithm of the proposed fatigue model is illustrated
in Fig. 5. The computational process involves updating the internal variables and updating
stresses. Once all internal variables are updated, the stresses in the bonding contact are then
TE
D
rectified from the stress predictors, 846(5 and 8,6(5 . Specifically, if the bonding contact is
yielding, the stresses are updated using the stress return algorithm given in the original model
(Nguyen et al., 2017b). On the other hand, if the contact is subjected to fatigue damage, the
EP
increments of the fatigue displacements and fatigue damage variable are updated at the onset
of unloading paths. The stresses at this position are then rectified using Eq. (38) and (39) to
account for the change of contact state caused by the increments of fatigue displacements and
AC
C
fatigue damage variable.
84A = 846(5
8,A
=
8,6(5
1 − :A
7
− @<4 1 − :A /49
1 − :A%C
(38)
1 − :A
7
− @<, 1 − :A /,9
A%C
1−:
(39)
where c and c − 1 denotes the current and previous number of load cycles, respectively.
16
EP
TE
D
M
AN
U
SC
RI
PT
ACCEPTED MANUSCRIPT
AC
C
Fig. 5: Flow chart of the implementation algorithm of the proposed fatigue model
3.3. Cycle jump technique
The high-cycle fatigue damage simulations of cemented materials involve thousands of load
cycles, which requires extensive computational time to model cycle by cycle the entire
fatigue damage process. In the simulations of regular cyclic loading scenarios with the same
loading range, a “cycle-jump” technique can be employed to reduce the computational cost.
This technique was first proposed in Chaboche, (1988) and Chaboche & Lesne, (1988) based
on the assumption that damage evolution in the bonding contacts is steady in a small range of
load cycles. In this sense, the simulation can “jump” from the first load cycle to the final load
cycle of the range and the total damage increment in the bonding contact can be obtained by
17
ACCEPTED MANUSCRIPT
extrapolating the damage increment in the first load cycle. This approach can considerably
reduce the actual number of simulated load cycles, hence saving the overall computational
time. In this study, the number of omitted cycles in a “cycle jump”, cuTe6 , is determined
using the following equation:
cuTe6 =
@:uTe6
@:e k
(40)
RI
PT
where @:uTe6 is an estimated damage increment after a cycle jump and @:e
k
is the
maximum damage increment in all bonding contacts at the first load cycle in the cycle jump
range. As the stress range in the bonding contact is considered to be unchanged in a cycle
SC
jump, the total increment of fatigue damage variable after the cycle jump can be derived as
@:7 =
AvwxV
uBC
R
hu √: u
hC √:C
e
Y @:7C
M
AN
U
follows:
(41)
where : u = :C + ∑zBC @:7y , which is the total damage variable of the bonding contact at the
u%C
beginning of load cycle {
|
in the cycle jump; hu is the geometrical factor corresponding to
TE
D
: u and @:7C is the incremental fatigue damage at the first load cycle in the cycle jump.
Details of the derivation of Eq. (41) is presented in Appendix A.
Model behaviour and the influences of fatigue parameters
In this section, the model behaviour at the constitutive level is examined by performing
EP
simulations of two particles connected by a bonding contact. The particle diameter and
contact width are both equal to 0.1 m. The model parameters assigned to the bonding contacts
are given in Table 1. The initial value of damage variable :)4) is assumed to be 1%. The
AC
C
4.
constitutive responses of the model are illustrated in two fundamental stress modes (pure
tension and pure shear) as shown in Fig. 6. For each mode, both monotonic and cyclic tests
are conducted to show the correlations between the damage and fatigue behaviour at the
constitutive level. In the monotonic test, the fatigue damage is deactivated in the model by
setting the fatigue threshold to be greater than one. As a result, the fatigue variables remain
zero and the constitutive model formulation is transformed to the original model formulation
in Nguyen et al., (2017b). In the cyclic tests, a cyclic loading regime with the maximum
stress level of 60 kPa and minimum stress level of 20 kPa is applied to the contact until the
damage variable reaches 0.4 (see Fig. 6). Subsequently, the contact is loaded monotonically
18
ACCEPTED MANUSCRIPT
until failure. As the maximum cyclic stress can reach to 60% of the contact cohesion, the
fatigue threshold has to be chosen less than 0.6 so that fatigue damage can develop in the
bonding contact.
Table 1: The values of the model parameters used in the pure tension and pure shear tests
SC
20.0
10.0
120
100
2.0×10-5
2.0×10-5
10
5
RI
PT
Selected value
M
AN
U
Model parameters
Elasto-plastic parameters
/49 (GPa/m)
/,9 (GPa/m)
8}9 (kPa)
K 9 (kPa)
<4 (m)
<, (m)
ϕ (0)
ψ (0)
Fatigue parameters
K
EP
TE
D
1.0
1.0×10-8
Fig. 6: Illustrations of the constitutive tests in pure tension and pure shear modes and the
AC
C
cyclic loading regime
4.1. Pure tension and pure shear test results
Fig. 7 illustrates the results of the pure tension tests. The normal stress-displacement curve
obtained in the monotonic test follows the same pattern as in the original model (Nguyen et
al., 2017b), including a linear elastic stage followed by a softening stage when the bonding
contact is yielding (see Fig. 7a). In the cyclic tests, the fatigue damage variable is shown to
evolve with increasing the number of load cycles (see Fig. 7b). This is also observed from the
gradual reduction of the slope of loading/unloading paths as shown in Fig. 7a. Along with the
development of fatigue damage, the fatigue displacement in the normal direction is shown to
increase. The increments of normal fatigue displacement and fatigue damage variable are
19
ACCEPTED MANUSCRIPT
induced from the hysteresis loops of loading/unloading paths, representing the micro-crack
growth and deforming mechanisms occurring in the fatigue damage process of cemented
materials. When damage variable in the contact is equal to 0.4, normal stress in the bonding
contact is still under the yielding limit. Thus, under a continuous strain increment, the normal
stress increases linearly until reaching the yield criterion, then softening in the same path as
RI
PT
in the monotonic test. Fig. 8 presents the results obtained in the pure shear test. It shows that
when the damage variable of the bonding contact is equal to 0.4, shear stress in the bonding
contact has just reached the yielding limit. As a result, when the monotonic load is applied,
the plastic damage immediately develops in the bonding contact and the shear stress is
TE
D
M
AN
U
SC
softened without any elastic stress increment as in the pure tension test.
Fig. 7: (a) The stress-displacement curves obtained in the pure tension tests; (b) the evolution
AC
C
EP
of damage variable in the pure tension cyclic test
Fig. 8: (a) The stress-displacement curves obtained in the pure shear tests; (b) the evolution
of damage variable in the pure shear cyclic test
20
ACCEPTED MANUSCRIPT
To examine the effectiveness of the proposed cycle jump technique, the pure tension and
shear cyclic tests are repeated in which the cycle jump technique is employed. The value of
@:uTe6 is chosen equal to 0.005 and 0.01. Fig. 9 shows the stress-displacement curves for
normal and shear tests obtained in the simulations using the cycle jump technique. The actual
number of simulated load cycles is shown to significantly reduce when using the cycle jump
technique. Specifically, the numbers of simulated load cycles until the damage variable
RI
PT
reaches 0.4 in both the pure tension and shear tests reduce by 388.8% when @:uTe6 = 0.005
and by 797.4% when @:uTe6 = 0.01, as compared to the simulation without the cycle jump
technique. Fig. 10 shows a comparison of the fatigue damage evolution results between the
SC
original simulation without the cycle jump technique and those using the cycle jump
technique. The curves of fatigue damage evolution versus the number of load cycles obtained
M
AN
U
in the simulations using the cycle jump technique excellently match with those obtained from
the original simulations. These results demonstrate the effectiveness of the cycle jump
technique for reducing the overall computational time in fatigue simulations while still
AC
C
EP
TE
D
ensuring the accuracy of numerical results.
21
ACCEPTED MANUSCRIPT
Fig. 9: The stress-displacement results in the pure tension and pure shear cyclic tests using
M
AN
U
SC
RI
PT
the cycle jump technique with (a) @:uTe6 = 0.005 and (b) @:uTe6 = 0.01
Fig. 10: The influences of cycle jump technique on the results of fatigue damage evolution in
(a) the pure tension cyclic test and (b) the pure shear cyclic test
4.2. The effects of fatigue parameters on the constitutive behaviour
The influences of fatigue parameters on the constitutive behaviour of the bonding contact are
TE
D
investigated in this section to provide a basis for calibrations of model parameters in
numerical experiments. The fatigue parameters of the proposed model consist of
and K,
which represent the fatigue properties of the bonding contact. Their effects are investigated
by repeating the pure tension and pure shear cyclic tests. However, the values of parameters
of parameter
EP
and K are varied from 0.8 to 1.2 and from 0.5×10-8 to 1.5×10-8, respectively. The influence
on the results of damage evolution against the number of cyclic loads in the
AC
C
pure tension cyclic test is depicted in Fig. 11a. The increment of fatigue damage in the
bonding contact is shown to be sensitive to the value of
. More specifically, it requires 2008
load cycles for the fatigue damage variable to increase from 0.01 to 0.4 when
the number of load cycles reduces to only 50 when
= 0.8, but
= 1.2 for the same damage increment.
This can be explained based on Eq. (27), which shows the increment of fatigue damage
variable is proportional to the stress intensity factor power of
. Therefore, a small change of
would significantly affect the evolution of fatigue damage variable. Fig. 11b displays the
results of fatigue damage evolution with respect to the number of load cycles when varying K.
The increment of fatigue damage variable is shown to be proportional to the value of K,
which is consistent with Eq. (27). Note that as the same loading magnitude is used in both the
22
ACCEPTED MANUSCRIPT
pure tension and pure shear cyclic tests, the results of fatigue damage evolution against the
SC
RI
PT
number of load cycles in these tests exactly match with each other.
M
AN
U
Fig. 11: Influences of fatigue parameters on the damage evolution in the pure tension cyclic
tests. The results of the pure shear cyclic test exactly match with those in the pure tension
cyclic tests
5.
Applications of the proposed modelling approach
The proposed modelling approach is now applied to simulate the experiments of cemented
TE
D
materials. The four-point bending test of cemented-treated samples subjected to both
monotonic and cyclic loads carried out by Sounthararajah et al., (2018) are chosen for the
simulation. This four-point bending test has been used in the pavement industry to develop
practical design models and measure the material fatigue properties for inputs in pavement
EP
designs. The numerical setup of this test is first described, followed by the calibration of
model parameters. Through the calibration process, the capability of the proposed modelling
AC
C
method in capturing the fatigue behaviour of cemented materials is evaluated by comparing
numerical results with experimental data. To further examine the predictive capability of the
proposed approach, the calibrated model parameters are applied to simulate the fatigue test at
different stress levels and numerical results of the S-N curve are then compared with the
experimental results.
5.1. Simulation setup
The experimental and numerical setups of the four-point bending tests are illustrated in Fig.
12. The beam samples are 200 mm in length, 50 mm in height and 50 mm in width, placed on
two supporting rollers at the span of 150 mm and loaded by two loading rollers placed on the
23
ACCEPTED MANUSCRIPT
top surface of the sample at the distance of 50 mm. A numerical sample in the same
dimensions is generated in DEM by an assemblage of poly-dispersed particles with radii
ranging from 1 mm to 19 mm. The particle size distribution of the assemblage is taken
similar to the aggregate size distribution of the physical material (see Fig. 13). Grains with
diameters less than 1 mm are considered to not participate well in the force transmission of
RI
PT
the aggregate skeleton. Thus, they are simulated by particles of 1 mm diameter for more
computational efficiency. All of these particles are then randomly distributed in the sample
domain and are allowed to rearrange until reaching an equilibrium condition. The supporting
and loading rollers are modelled in DEM by wall elements. Only the intersecting areas
SC
between the rollers and samples, instead of the entire rollers, are modelled in order to reduce
the allocated searching domain in PFC 5.0, thereby saving the computational time. The quasistatic condition as described in Section 2 is set for these simulations. The proposed fatigue
M
AN
U
model is subsequently incorporated into all the inter-particle contacts in the sample. On the
other hand, the classical cohesionless model developed by Cundall & Strack, (1979) is
incorporated into the wall-to-particle contacts and frictionless interactions are set for all of
these contacts.
Next, the numerical sample is loaded by applying velocities to the loading walls. The
TE
D
monotonic and cyclic loading patterns in the experiments are reproduced in the simulations
using a wall servo-control technique. The basis of this technique is to control the velocity of
the loading walls in a way that the interactional forces between the walls and samples match
the objective loading patterns. During the loading processes, the values of applied load, •,
and mid-span deflection, @, are continuously recorded. The flexural stress, flexural strain and
equations:
εb =
Eb =
P‚
W„ $
AC
C
FS =
EP
flexural modulus of the material are then calculated from • and @ using the following
(42)
108„@
23‚$
(43)
23• ‚ˆ
108@ W„ ˆ
(44)
where FS is the flexural stress, ‚ is the support span, ‰ is the average beam width, „ is the
average beam height, εb is the flexural strain, Eb is the flexural modulus, • is the maximum
magnitude of applied cyclic load and @ is the sample deflection corresponding to • .
24
ACCEPTED MANUSCRIPT
Similar to the experiments, the initial flexural modulus EbŠ‹Š is determined as the average
value of the flexural modulus results obtained in first 50 load cycles in the fatigue test. The
fatigue cracking process of the sample can be observed at the macroscopic scale through the
curve of Eb /EbŠ‹Š versus the number of load cycles, or so-called the flexural modulus
SC
RI
PT
degradation curve.
M
AN
U
Fig. 12: The experimental and numerical setups of the four-point bending test. Particles with
AC
C
EP
TE
D
diameters more than 1mm are coloured darker than particles with diameters equal to 1mm
Fig. 13: The particle size distributions of experimental and numerical samples
5.2. Model calibration
The proposed constitutive model consists of two sets of parameters: elasto-plastic parameters
and fatigue parameters. These parameters are calibrated based on a trial and error process by
directly comparing the simulation results to those of the experiments. The elasto-plastic
parameters are first calibrated with the monotonic test data. Details of the calibration
procedure can be found in Nguyen et al., (2017b). The calibrated values of elasto-plastic
parameters are then used in the simulation of fatigue test to further calibrate the fatigue
25
ACCEPTED MANUSCRIPT
parameters. The calibration procedure of the fatigue parameters is as follows. The fatigue
threshold is first calibrated. This parameter controls the number of contacts subjected to
fatigue in the sample, hence affecting the sample deformation corresponding to its fatigue
damage evolution. The fatigue threshold is calibrated to fit the deformation of experimental
samples at the final stage, which is represented by the value of Eb /EbŠ‹Š . Parameters
is then
calibrated to fit the flexural modulus degradation curve at the initial fatigue stage. Finally,
RI
PT
parameter K is calibrated to fit the fatigue life of a numerical sample with its experimental
counterpart.
The deformable method proposed in Potyondy & Cundall, (2004) is also adopted in the
calibration process to make the constitutive law independent of the contact size. Accordingly,
SC
the contact normal and shear stiffness are replaced by the effective modulus, Π9 , and normal-
M
AN
U
to-shear stiffness ratio, /49 ⁄/,9 , in the calibration. The same values of the effective modulus
and normal-to-shear stiffness ratio are assigned to all bonding contacts in numerical samples.
The normal and shear stiffnesses of each bonding contact are then back-calculated using the
relationship given in Nguyen et al., (2017b). The calibrated values of the model parameters
for the chosen quasi-static condition (i.e. global damping coefficient of 0.7) are presented in
simulations.
TE
D
Table 2. The initial value of damage variable :)4) is assumed to be equal to 1% in all
Table 2: The values of calibrated model parameters for the flexural tests
AC
C
EP
Model parameters
Elasto-plastic parameters
Π9 (GPa)
/49 ⁄/,9
8}9 (kPa)
K 9 (kPa)
<4 (m)
<, (m)
ϕ (0)
ψ (0)
Fatigue parameters
K
a
Calibrated values
1.55
2.4
74
74
2.0×10-5
2.0×10-5
20
10
1.2
1.8×10-11
0.84
Fig. 14 compares the flexural stress-strain curves obtained in the experiments and simulation
for the monotonic test. It shows that the numerical curve falls within the experimental range.
26
ACCEPTED MANUSCRIPT
During the loading process, the cemented sample undergoes three stages including a linear
elasticity stage, a non-linear hardening stage prior to peak stress and a residual stage after the
peak. All of these observed responses are captured well in the simulation. In addition, the
cracking pattern in the numerical sample agrees with the experimental result as illustrated in
Fig. 15. In both the experiment and simulations, the pattern of major cracks can slightly vary
RI
PT
depending on the arrangement of particles, but the major cracks always initiate at the bottom
of the beam in the area between the two loading rollers and then propagates upwards along
particle interfaces. The localised fracture and crack development are naturally captured in the
simulation without any ad-hoc treatment on the contact model. This demonstrates the
SC
effectiveness of the proposed method in modelling the damage response of cemented
EP
TE
D
M
AN
U
materials under monotonic loads.
AC
C
Fig. 14: The flexural stress-strain curves obtained in the experiment and simulation of the
monotonic test.
Fig. 15: Sample crack patterns in the experiment and simulation of the monotonic test
The numerical fatigue test is then performed with the maximum cyclic stress equal to 1029
kPa, which is 77% of the maximum flexural stress obtained in the monotonic test (FS•Ž• =
1336.4 kPa). This corresponds to the maximum cyclic load of 857.5 N according to Eq. (42).
27
ACCEPTED MANUSCRIPT
The cyclic load in haversine shape as illustrated in Fig. 16 is continuously applied to the
numerical sample until failure. The cycle jump technique is also employed in this simulation
with @:uTe6 = 0.003. This chosen value of @:uTe6 was checked and shown to marginally
affect the fatigue response of the sample at macro-scale.
Fig. 17a shows a comparison of the flexural modulus degradation curves obtained in the
RI
PT
experiment and simulations using the proposed fatigue model and the original model
(Nguyen et al., 2017b). The reduction of flexural modulus in the experiment is observed to go
through three stages including the initial, middle and final stages. The simulation using the
proposed fatigue model can closely capture the gradual reduction of flexural modulus in the
SC
middle stage and the unstable reduction of flexural modulus at the final stage. However, there
stands a deviation in the initial stage in which the flexural modulus in the simulation is
M
AN
U
observed to drop faster than the flexural modulus in the experiment. This can be attributed to
the difference in microstructural configuration between the experimental and numerical
samples. The reduction of flexural modulus during the cyclic loading process is the
consequence of damage evolution in the numerical sample, evidenced by the reduction of the
slope of loading/unloading curves with increasing the number of load cycles as can be seen in
Fig. 17b. Along with the damage evolution, permanent deformation gradually increases in the
TE
D
sample, exhibited by that flexural strain values measured at the end of unloading paths
gradually shifts away from the origin (see Fig. 17b). On the contrary to the simulation result
of the proposed fatigue model, the simulation using the original cohesive model is unable to
capture the reduction of the flexural modulus as well as the sample failure in the cyclic test.
EP
The flexural modulus in the simulation with the original cohesive model remains constant
with increasing the number of load cycles (see Fig. 17a). These results clearly show the
AC
C
capability as well as the necessity of the proposed fatigue model for the DEM simulations of
fatigue damage in cemented materials.
28
ACCEPTED MANUSCRIPT
Fig. 16: The shape of the applied cyclic load used in the fatigue test
Fig. 18 illustrates the damage evolution and crack development in the numerical sample
during the cyclic test. In the initial stage, damage starts to evolve in several places in the area
between two loading rollers but no clear fracture is observed in the numerical sample. During
RI
PT
the middle stage, damage is shown to localise at the bottom mid-span of the sample, then
gradually developing with increasing number of load cycles. At the final stage, a major crack
initiates and rapidly propagates from the bottom to the top surface of the sample. The final
crack pattern in the numerical sample matches well with the crack patterns observed in the
monotonic test and in the experiment. This result indicates the capability of this modelling
TE
D
M
AN
U
SC
approach in reproducing fatigue cracking in cemented samples.
Fig. 17: (a) The flexural modulus degradation curves obtained in the experiment and
EP
simulations with the cyclic stress level equal to 77% of the maximum flexural stress; (b) The
AC
C
flexural stress-strain result obtained in the simulation of the monotonic and fatigue tests
29
SC
RI
PT
ACCEPTED MANUSCRIPT
Fig. 18: Damage evolution in the numerical sample during the fatigue test
M
AN
U
5.3. Numerical prediction
In Section 5.2, the proposed modelling approach has been shown to be able to capture the
fatigue behaviour of cemented materials. In this section, additional simulations of the fatigue
test at different maximum cyclic stress levels (i.e. 95%, 92%, 84%, 81% and 70% of FS•Ž• )
are carried out to examine its predictive capability. The numerical sample and boundary
condition remain the same as in Section 5.1. The same model parameters as given in Table 2
TE
D
are used in all of these simulations. The results of S-N curves obtained in the simulation are
then compared with the S-N curves of the experiments. Similar to the experimental method,
the fatigue life of the numerical sample is defined as the number of load cycles applied to that
AC
C
EbŠ‹Š = 50%.
EP
sample until the flexural modulus is equal to 50% of the initial flexural modulus, i.e. Eb /
30
ACCEPTED MANUSCRIPT
Fig. 19: (a) The S-N curve obtained in the experiment and simulation of the flexural fatigue
test; (b) The flexural modulus degradation curves obtained in the numerical fatigue tests with
different maximum cyclic stress levels
Fig. 19a displays the S-N curves obtained in the experiments and simulations. It shows that,
RI
PT
given a single set of model parameters and a microstructural configuration of the numerical
sample, the proposed modelling approach can predict the decrease of sample fatigue life with
increasing the maximum cyclic stress level. Similar to the experiment, the results of fatigue
life obtained in the simulations with the maximum stress level ranging from 77% to 95% of
SC
FS•Ž• relatively fall to a logarithmic trend line closing to the experimental trend line. The
small deviation is considered due to the differences in the microstructural configurations of
the numerical and experimental samples. On the other hand, the numerical sample subjected
M
AN
U
to the cyclic stress level equal to 70% of FS•Ž• is shown to be self-stabilised after hundreds
of load cycles. This can be attributed to the shakedown effect, in which particles and contacts
in the sample gradually rearrange to a more stable condition when undergoing cyclic loading.
The difference between the experimental and numerical results, in this case, can be attributed
to the randomness of sample microstructures in the experiment and simulation. Fig. 19b
depicts the flexural modulus degradation curves obtained in the numerical fatigue tests with
TE
D
different stress levels. It can be seen that when increasing the applied cyclic stress level, the
sample fails earlier and the value of Eb /EbŠ‹Š at the failure points, where the curves suddenly
drop, is also larger. This is due to the fact that the reduction of Eb /EbŠ‹Š is as a result of fatigue
EP
damage evolution in the sample, while the sample failure occurs when the applied cyclic
stress reaches the yield criterion after a number of load cycles. A higher cyclic stress level
AC
C
applied to numerical samples results in greater stress magnitudes in the bonding contacts.
Hence, the contacts can reach their yield criterion faster while the evolution of fatigue
damage is still small. The sample also fails earlier and the value of Eb /EbŠ‹Š at the failure
points is greater than those of the fatigue tests with lower cyclic stress levels.
5.4. The influence of particle arrangement on the numerical results
In Section 5.2, the model parameters are calibrated with a specific microstructural
configuration of the numerical sample. However, in the experiments of cemented materials,
the microstructure configuration of samples would vary owing to the randomness of
aggregates and cement bonds arrangement in different samples. This is one of the main
causes of highly scattered results observed in the fatigue tests of cemented materials. As the
31
ACCEPTED MANUSCRIPT
proposed modelling approach can directly replicate the heterogeneous microstructure of
cemented materials through the sample generation process, the influence of particle
arrangement on the material behaviour can be naturally reproduced in simulations. To
examine the influence of particle arrangement on numerical results, four additional beam
samples are generated, namely Samples 1 to 4. The sample dimensions and particle size
RI
PT
distribution remain the same as the original sample given in Section 5.1, yet the particle
arrangement is varied differently in each new sample. The same values of model parameters
as given in Table 2 are assigned to all bonding contacts in these new samples. Both the
monotonic and fatigue tests at the maximum stress level of 1029 kPa (i.e. 77% of FS•Ž• of
SC
the original sample) are then performed on the new samples and numerical results are
compared with those of the original sample.
M
AN
U
Fig. 20 compares the flexural stress-strain curves and the flexural modulus degradation
curves obtained in numerical tests of the new samples with those of the original sample. It
can be seen while the results of the monotonic test are relatively consistent and falling well
within the experimental range (Fig.20a), the results of the fatigue test vary substantially
among the numerical sample (Fig.20b). The scatter results of the fatigue tests can be
attributed to the difference in the microstructure of each numerical sample due to the
TE
D
rearrangement of particles. In fact, the variation of the material microstructure can affect the
distribution of local stress in bonding contacts, causing changes of the damage evolution rate
per load cycle in each bonding contact. As the fatigue behaviour of the sample at the
macroscopic scale is the collective result of all the damage evolution per load cycle, the
EP
change of damage evolution in every load cycle would result in a considerable deviation in
the final fatigue life of the samples. These responses are consistent with the experimental
AC
C
observation of cemented materials, which normally shows high scattered results of fatigue
life in different samples (Do et al., 1993; Dyduch et al., 1994).
32
ACCEPTED MANUSCRIPT
Fig. 20: The variation of flexural strength results (a) and fatigue results (b) obtained in the
simulations when varying the particle size distribution of samples
To provide an insight into the scatter results of the fatigue test, the coordination number is
employed to explain the relationship between the particle arrangement in numerical samples
RI
PT
and their mechanical responses. In DEM modelling, the coordination number of a sample is a
determining factor of the sample packing quality, which is obtained by the following equation:
C4 =
2c
c6
(45)
where c and c6 are the total number of contacts and particles in the sample, respectively.
SC
Fig. 21 presents the relationships between the coordination numbers C4 and the results of
M
AN
U
elastic flexural modulus and fatigue life of the samples, respectively. The elastic flexural
modulus and fatigue life may increase with the increments of the coordination number in the
numerical samples. These observations are reasonable in the sense that the increment of the
coordination number represents the increase of the number of contacts to encounter loads in
numerical samples, which improves the sample stiffness and load-carrying capacity. However,
it is acknowledged that the fatigue response of numerical samples also depends on how stress
TE
D
transmits and concentrates in the samples during the cyclic loading process. The magnitude
of stress concentration can vary from one sample to another due to the differences in their
microstructure configuration. As the evolution of fatigue damage relies on the magnitude of
stress in accordance with Paris’ law, the variation of local stress causes the rate of fatigue
EP
damage evolution to be varied in different samples. This microscopic mechanism contributes
to the scattered results observed in the numerical fatigue tests of cemented materials. As a
AC
C
result of such microstructural effect, the fatigue life of samples 2 and 3 is shown to be quite
similar although their coordination numbers are relatively different (see Fig. 21). This
behaviour can be further explained by having a closer look at the stress distribution in these
numerical samples. Fig. 22 illustrated the distribution of tensile stress in the numerical
samples under the application of flexural stress equal to 400 kPa. It shows that even though
the same stress is applied to these samples, the magnitude of tensile stress in bonding contacts
at the location where cracks are initiated (zoomed area) of sample 3 is greater than those in
other samples. This suggests that the rate of fatigue damage evolution in sample 3 is higher,
and thus this sample can reach its fatigue life faster than the others.
33
RI
PT
ACCEPTED MANUSCRIPT
SC
Fig. 21: The relationships between the coordination number and the elastic flexural modulus
and fatigue life of the samples
M
AN
U
Fig. 23 depicts the crack patterns of the experimental and numerical samples obtained in the
fatigue tests. In all experiments and simulations, cracks initiate in the area between the two
loading rollers and develop from the bottom to the top surface of the samples. However, the
location and shape of cracks are varied in each sample, and this can also be attributed to the
difference in the microstructural configuration of each sample due to the rearrangement of
particles. The difference of sample microstructure induces the change of stress transmission
TE
D
and concentration in the sample, thus changing the crack location. Besides, as the crack
propagates along particle interfaces, the direction of crack propagation highly depends on
how particles arrange near the crack tip. The variation of crack patterns can be well captured
AC
C
behaviour.
EP
by the proposed approach thanks to it physical reproduction of the material microstructural
34
ACCEPTED MANUSCRIPT
Fig. 22: The distribution of tensile stress in numerical samples under the application of
M
AN
U
SC
RI
PT
flexural stress of 400 kPa
Fig. 23: Crack patterns of experimental and numerical samples obtained in the fatigue tests
Conclusion
A discrete-based approach is proposed in this study to model the fatigue damage of cemented
TE
D
materials. This modelling approach takes advantages of DEM to reproduce the heterogeneous
microstructure of cemented materials and capturing real crack development. Alongside with
the application of DEM, a new fatigue bonding contact model is developed to characterise the
fatigue response of cemented materials at the grain scale. This model is derived on the basis
EP
of the cohesive damage-plasticity model proposed in Nguyen et al., (2017b), along with
fatigue damage variables and a fatigue evolution law, to describe both cracking and
irreversible displacements in the fatigue process of cemented materials. Fatigue damage is
AC
C
6.
allowed to evolve when stress is inside the yield surface through the introduction of a fatigue
criterion under the yield criterion. The evolution of fatigue damage is then described by a
single fatigue law of a scalar fatigue damage variable (i.e. Paris’ law), while different stages
of fatigue crack growth exhibited in the materials at the macroscopic scale are captured as the
collective response of all bonding contacts and particles in numerical samples.
Through numerical flexural tests, the capability of the proposed modelling approach in
capturing and predicting the fatigue behaviour of cemented materials under flexural boundary
condition is demonstrated. This is evidenced by the good agreements in the flexural modulus
degradation curves, the S-N curves and sample crack patterns between the experiments and
35
ACCEPTED MANUSCRIPT
simulations. In addition, thanks to the capability of DEM in reproducing physical material
microstructure and the explicit characterisation of contact behaviour in the new fatigue model,
the proposed modelling approach is able to reproduce the microstructural effects on the
fatigue behaviour of cemented materials. The variations of stress-strain results, fatigue life
and crack patterns in the samples with changing particle arrangement can be naturally
RI
PT
captured in the DEM simulations. This demonstrates the advantages of the proposed
approach over existing continuum methods for modelling the fatigue response and localised
fracture in cemented materials, while still retaining the simplicity of contact model
formulation. Moreover, the incorporation of fatigue damage mechanisms of cement bridges
SC
into the DEM bonding contact model is shown to be necessary given the existing damageplastic contact model is unable to capture these mechanisms in the DEM simulations at the
M
AN
U
grain scale.
Acknowledgement
Funding supported by Australian Research Council (grant number LP130100884 and DP
160100775) is gratefully acknowledged.
Appendix A: Derivation of the total increment of fatigue damage variable in a cycle
jump
TE
D
∆d is first derived by substituting Eq. (32) and (33) into (28) as:
∆d = h√jf:.∆84$ + ∆8,$
EP
The incremental fatigue damage in a load cycle is then obtained by substituting the above
equation into Eq. (27) as:
K ∆d
@:7 =
f
K;√jf.∆84$ + ∆8,$ =
= ;h√:=
f
e
AC
C
e
e
Hence, the incremental fatigue damage in each load cycle of a cycle jump is calculated as:
K;√jf.∆84$ + ∆8,$ =
@:7C = mhC .:C n
f
e
@:7$
…
= mh
e
K;√jf.∆84$ + ∆8,$ =
:$ n
f
$.
e
e
e
h $ √:$
=R
Y @:7C
C
C
h √:
36
ACCEPTED MANUSCRIPT
u
@:7
K;√jf.∆84$ + ∆8,$ =
u. u
= mh : n
f
e
e
=R
hu √: u
hC √:C
e
Y @:7C
The total increment of fatigue damage variable in a cycle jump is then obtained as follows:
uBC
u
@:7
=
AvwxV
uBC
Reference
hu √: u
e
R
Y @:7C
C
C
h √:
RI
PT
@:7 =
AvwxV
SC
Al-Rub, R. K. A. & Darabi, M. K. (2012). A thermodynamic framework for constitutive
modeling of time-and rate-dependent materials. Part I: Theory. International Journal
of Plasticity, 34, 61-92.
Alliche, A. (2004). Damage model for fatigue loading of concrete. International Journal of
Fatigue, 26, 915-921.
M
AN
U
Carrara, P. & De Lorenzis, L. (2015). A coupled damage-plasticity model for the cyclic
behavior of shear-loaded interfaces. Journal of the Mechanics and Physics of Solids,
85, 33-53.
Chaboche, J.-L. (1988). Continuum damage mechanics: Part II—Damage growth, crack
initiation, and crack growth. Journal of applied mechanics, 55, 65-72.
Chaboche, J. L. & Lesne, P. M. (1988). A non‐linear continuous fatigue damage model.
Fatigue & fracture of engineering materials & structures, 11, 1-17.
TE
D
Cheng, J., Qian, X. & Zhao, T. (2016). Rheological viscoplastic models of asphalt concrete
and rate-dependent numerical implement. International Journal of Plasticity, 81, 209230.
Cundall, P. A. & Strack, O. D. L. (1979). A discrete numerical model for granular assemblies.
Géotechnique, 29, 18.
AC
C
EP
Darabi, M. K., Al-Rub, R. K. A., Masad, E. A., Huang, C.-W. & Little, D. N. (2012). A
modified viscoplastic model to predict the permanent deformation of asphaltic
materials under cyclic-compression loading at high temperatures. International
Journal of Plasticity, 35, 100-134.
Dattoma, V., Giancane, S., Nobile, R. & Panella, F. W. (2006). Fatigue life prediction under
variable loading based on a new non-linear continuum damage mechanics model.
International Journal of Fatigue, 28, 89-95.
Do, M.-T., Chaallal, O. & Aïtcin, P.-C. (1993). Fatigue behavior of high-performance
concrete. Journal of Materials in civil Engineering, 5, 96-111.
Dyduch, K., Szerszeń, M. & Destrebecq, J.-F. (1994). Experimental investigation of the
fatigue strength of plain concrete under high compressive loading. Materials and
Structures, 27, 505-509.
Fathima, K. M. P. & Kishen, J. M. C. (2013). A thermodynamic framework for fatigue crack
growth in concrete. International Journal of Fatigue, 54, 17-24.
Gao, L. & Hsu, C.-T. T. (1998). Fatigue of concrete under uniaxial compression cyclic
loading. Materials Journal, 95, 575-581.
37
ACCEPTED MANUSCRIPT
Gui, Y., Bui, H. H. & Kodikara, J. (2015). An application of a cohesive fracture model
combining compression, tension and shear in soft rocks. Computers and Geotechnics,
66, 142-157.
Gui, Y. L., Zhao, Z. Y., Kodikara, J., Bui, H. H., & Yang, S. Q. (2016). Numerical modelling
of laboratory soil desiccation cracking using UDEC with a mix-mode cohesive
fracture model. Engineering Geology, 202, 14-23.
RI
PT
Horii, H., Shin, H. C. & Pallewatta, T. M. (1992). Mechanism of fatigue crack growth in
concrete. Cement and concrete composites, 14, 83-89.
Irwin, G. R. (1957). Analysis of stresses and strains near the end of a crack traversing a plate.
Journal of applied mechanics, 24, 361-364.
Itasca (2014). PFC 5.0 User Manual. Itasca Consulting group, Inc, Washington DC, United
State, p. 210.
SC
Krajcinovic, D. & Fanella, D. (1986). A micromechanical damage model for concrete.
Engineering Fracture Mechanics, 25, 585-596.
M
AN
U
Liang, R. Y. & Zhou, J. (1997). Prediction of fatigue life of asphalt concrete beams.
International Journal of Fatigue, 19, 117-124.
Mohamed, A. R. & Hansen, W. (1999). Micromechanical modeling of concrete response
under static loading—part 1: model development and validation. Materials Journal,
96, 196-203.
Moreo, P., Garcia-Aznar, J. M. & Doblare, M. (2007). A coupled viscoplastic rate-dependent
damage model for the simulation of fatigue failure of cement–bone interfaces.
International journal of plasticity, 23, 2058-2084.
TE
D
Nguyen, N., Bui, H. H., Nguyen, G. D., Kodikara, J., Arooran, S. & Jitsangiam, P. (2017a). A
thermodynamics-based cohesive model for discrete element modelling of fracture in
cemented materials. International Journal of Solids and Structures, 117, 159-176.
EP
Nguyen, N. H. T., Bui, H. H., Nguyen, G. D. & Kodikara, J. (2017b). A cohesive damageplasticity model for DEM and its application for numerical investigation of soft rock
fracture properties. International Journal of Plasticity, 98, 175-196.
AC
C
Nguyen, T. T., Bui, H. H., Ngo, T. D., & Nguyen, G. D. (2017c). Experimental and
numerical investigation of influence of air-voids on the compressive behaviour of
foamed concrete. Materials & Design, 130, 103-119.
Nguyen, O., Repetto, E. A., Ortiz, M. & Radovitzky, R. A. (2001). A cohesive model of
fatigue crack growth. International Journal of Fracture, 110, 351-369.
Oñate, E., Zárate, F., Miquel, J., Santasusana, M., Celigueta, M., Arrufat, F., Gandikota, R.,
Valiullin, K. & Ring, L. (2015). A local constitutive model for the discrete element
method. Application to geomaterials and concrete. Computational Particle Mechanics,
2, 139-160.
Papa, E. & Taliercio, A. (1996). Anisotropic damage model for the multiaxial static and
fatigue behaviour of plain concrete. Engineering Fracture Mechanics, 55, 163-179.
Paris, P. C. & Erdogan, F. (1963). A critical analysis of crack propagation laws. Journal of
basic engineering, 85, 528-533.
Paris, P. C. & Sih, G. C. (1965). Stress analysis of cracks. Fracture toughness testing and its
applications. ASTM International.
38
ACCEPTED MANUSCRIPT
Potyondy, D. O. & Cundall, P. A. (2004). A bonded-particle model for rock. International
Journal of Rock Mechanics and Mining Sciences, 41, 1329-1364.
Ray, S. & Kishen, J. M. C. (2010). Fatigue crack propagation model for plain concrete–An
analogy with population growth. Engineering Fracture Mechanics, 77, 3418-3433.
Ristinmaa, M. (1995). Cyclic plasticity model using one yield surface only. International
journal of plasticity, 11, 163-181.
RI
PT
Scholtès, L. & Donzé, F.-V. (2013). A DEM model for soft and hard rocks: Role of grain
interlocking on strength. Journal of the Mechanics and Physics of Solids, 61, 352-369.
Shen, W. Q., Shao, J.-F., Kondo, D. & Gatmiri, B. (2012). A micro–macro model for clayey
rocks with a plastic compressible porous matrix. International journal of plasticity, 36,
64-85.
SC
Shen, W. Q. & Shao, J. F. (2016). An incremental micro-macro model for porous
geomaterials with double porosity and inclusion. International Journal of Plasticity,
83, 37-54.
M
AN
U
Simon, K. M. & Kishen, J. M. C. (2017). A multiscale approach for modeling fatigue crack
growth in concrete. International Journal of Fatigue, 98, 1-13.
Sounthararajah, A., Bui, H. H., Nguyen, N. H. T., Jitsangiam, P. & Kodikara, J. (2018).
Early-Age Fatigue Damage Assessment of Cement-Treated Bases Under Repetitive
Heavy Traffic Loading. Journal of materials in civil engineering, (in press).
Tada, H., Paris, P. C. & Irwin, G. R. (1973). The stress analysis of cracks. Handbook, Del
Research Corporation.
TE
D
Turon, A., Costa, J., Camanho, P. P. & Dávila, C. G. (2007). Simulation of delamination in
composites under high-cycle fatigue. Composites Part A: applied science and
manufacturing, 38, 2270-2282.
Voyiadjis, G. Z. & Abu-Lebdeh, T. M. (1994). Plasticity model for concrete using the
bounding surface concept. International Journal of Plasticity, 10, 1-21.
EP
Weichert, D. & Maier, G. (2014). Inelastic behaviour of structures under variable repeated
loads: direct analysis methods, Springer.
Xiao, Y. C., Li, S. & Gao, Z. (1998). A continuum damage mechanics model for high cycle
fatigue. International journal of fatigue, 20, 503-508.
AC
C
Yang, B., Mall, S. & Ravi-Chandar, K. (2001). A cohesive zone model for fatigue crack
growth in quasibrittle materials. International Journal of Solids and Structures, 38,
3927-3944.
Zhang, J. (2001). Size effect on fatigue in bending of concrete. Journal of materials in civil
engineering, Vol. 13.
Zhang, J., Shen, W. Q., Oueslati, A. & De Saxce, G. (2017). Shakedown of porous materials.
International Journal of Plasticity, 95, 123-141.
39
ACCEPTED MANUSCRIPT
Highlights
AC
C
EP
TE
D
M
AN
U
SC
RI
PT
• The fatigue behaviour of cemented material can be effectively captured by DEM combining
with a new fatigue contact model
• The fatigue contact model coupling damage mechanics and plasticity theory can describe well
the degrading response of cemented materials subjected to cyclic loading
• This proposed modelling method facilitates the characterisation of fatigue damage in
cemented materials for better understanding and design of the material behaviour
Документ
Категория
Без категории
Просмотров
1
Размер файла
9 563 Кб
Теги
007, 2018, ijplas
1/--страниц
Пожаловаться на содержимое документа