Analysis of Dynamic Loading and Penetration of Soils

International Conference on Geotechnical Engineering. 2013
Analysis of Dynamic Loading and Penetration of Soils Applications to Site Investigation and Ground Improvement
J. P. Carter, M. Nazem
Australian Research Council Special Research Centre for Geotechnical Science and Engineering, The
University of Newcastle, Callaghan NSW, Australia
and
D.W. Airey
Centre for Geotechnical Research, The University of Sydney, Sydney, NSW, Australia
ABSTRACT: A robust numerical method is suggested for dealing with the complex and difficult
problem of analysing dynamic soil penetration and compaction. The approach is based on the
Arbitrary Lagrangian-Eulerian (ALE) method of analysis, whose main features and challenges are
briefly described. The method is then employed to perform a parametric study of a variety of
penetrometers free-falling into layers of inhomogeneous clay deforming under undrained
conditions. The effect of the mechanical properties of the clay soil on the penetration
characteristics is discussed and some of the predictions are compared with the results of laboratory
experiments. How these results might be used in the interpretation of in situ penetration tests is
then suggested. The numerical method is also adopted to examine the fundamental compaction
behaviour of soils and the trends in predicted results are broadly compared with the results of
physical experiments.
1
INTRODUCTION
Successfully analysing the dynamic compaction and penetration of soil deposits, by either free
falling or propelled objects, is one of the most challenging problems in geomechanics. To date
most methods used to predict the effects of an object impacting and possibly penetrating the ground
surface have been empirical, based on observations in field trials and laboratory experiments.
However, there is now a real need for rigorous theories to describe the fundamental response of the
soil in this situation, theories that are soundly based on the principles of impact mechanics and
theoretical soil mechanics. One aim of this paper is to describe such a theoretical approach.
The use of a falling weight (or “pounder”) is the most common means of applying a dynamic
impact to the soil surface in order to cause compaction of the underlying material, although in
recent times the use of impact rollers has increased in this area of ground treatment and
improvement. The analysis of a free-falling mass will be described in later sections of this paper.
Obviously the theory describing this type of behaviour has much in common with that required to
analyse the impact and subsequent penetration of a soft soil.
Penetrometers are widely used to investigate the mechanical properties of soils and to embed
objects in soil deposits, especially anchors and other items in the sea floor. Free falling
penetrometers have been employed to provide information on the strength of seabed sediments.
Such tests can be envisaged as the dynamic equivalent of a static cone penetration test (CPT).
They can provide useful data such as the total depth and time of penetration and the deceleration
characteristics of the falling penetrometer. Potentially, these data can then be used to deduce
strength parameters for the soil in situ. Other types of projectiles are also used to embed anchors
and other objects deeply in the sea floor. Large scale anchors of this type have been employed for
ship moorings and the mooring of floating (offshore) petroleum production facilities.
To date, dynamic penetrometers have been used for offshore oil and gas industry applications,
such as determining soil strengths for pipeline feasibility studies and anchoring systems, in military
applications for naval mine countermeasures and terminal ballistic studies, in extra-terrestrial
exploration, and they have also been proposed and investigated for deep sea nuclear waste disposal
(Chow 2011). The specifications of a dynamic penetrometer, including its shape, geometry, mass,
and initial velocity, normally depend on its application. For instance, torpedoes used to anchor
flexible risers are typically 12-15 m long, weigh 240-950 kN in air, and are 762-1077 mm in
diameter (Medeiros 2002), while the STING penetrometer, developed by Canada’s Defence
Research Establishment Atlantic, to investigate the upper sea bed consists of a 1 m long 19 mm
diameter steel rod with an enlarged tip at its base with diameters of 25-70 mm (Mulhearn 2002).
Numerical and experimental studies have shown that the penetration characteristics of these
penetrometers depend on their geometry and impact energy as well as the mechanical properties of
the soil layer, particularly its undrained properties. Scott (1970), Dayal and Allen (1973), True
(1975), and Beard (1977) are among the pioneers who presented the earliest applications of free
falling penetrometers in providing information on the upper few metres of seabed sediments. At
the time of these studies powerful numerical methods were not available to solve the corresponding
theoretical problem and researchers had to rely on empirical correction factors to interpret
empirical test data. Scott (1970) also demonstrated that technological limitations prevented reliable
data from being obtained. Over the last decade, improvements in sensor technology have led to
renewed interest in the dynamic penetration resistance of soil with several new penetrometer
systems being developed to investigate soil strength and bearing capacity. One of the factors
preventing the widespread use of this technology is the uncertainty in the method of interpretation.
For example, experimental studies (Mulhearn, 2002, Abelev et al, 2009) with penetrometers of
different diameters have indicated limitations of the widely used empirical methods of
interpretation developed by True (1975), and have demonstrated the need for a better understanding
of the factors affecting the dynamic soil resistance.
For a variety of reasons, the numerical simulation of these types of impact compaction and
penetration problems is probably one of the most difficult problems in theoretical geomechanics.
In the case of poorly draining material the analyst must consider the relative incompressibility of
the soil during the short testing times required for full penetration and the likely inhomogeneity and
rate dependence of the soil deposited on the seabed. For free-draining materials the volume
changes resulting from the passage of stress waves through the soil must also be considered. These
aspects of the soil response require a robust stress-integration algorithm to accurately evaluate the
stresses for a given strain field. Inertia effects may also be important in cases where the penetration
occurs rapidly. Moreover, the penetration of objects into layers of soil usually involves severe
mesh distortion caused by the large deformations. More importantly, the boundary conditions of
the problem do not remain constant during the analysis since the contact surface between the
penetrometer and the soil changes continuously during penetration. Finally, an advanced timeintegration scheme must be employed to accurately predict the highly nonlinear response of the
soil.
This paper describes a robust numerical method for dealing with such complex and difficult
problems. The favoured approach is based on the Arbitrary Lagrangian-Eulerian (ALE) method of
analysis, whose main features and challenges are described briefly in the paper. Application of this
method to the simulation of the dynamic compaction of soils and the dynamic impact compaction
and penetration of instruments into undrained layers of non-uniform soil is discussed, and
comparisons with experimental observations are made.
2
2.1
NUMERICAL APPROACH
Background
The Arbitrary Lagrangian-Eulerian (ALE) method, being well established in fluid mechanics and
solid mechanics, has been shown to be robust and efficient in solving a wide range of static as well
as dynamic geotechnical problems involving large deformations. Nazem et al. (2006) presented a
robust mesh optimisation procedure used with the operator-split technique for solving large
deformation problems in geomechanics. This technique refines a mesh by relocating the nodal
points on all material boundaries followed by a static analysis. Later, Nazem et al. (2008)
presented an ALE formulation for solving elastoplastic consolidation problems involving large
deformations. These authors compared the performance of the ALE method with the Updated
Lagrangian (UL) method and showed the efficiency and robustness of the ALE method by solving
some classical problems such as the consolidation of footings and cavity expansion. Sheng et al.
(2009) addressed a successful application of the ALE method, incorporated with an automatic load
stepping scheme and a smooth contact discretisation technique, in solving geotechnical problems
involving changing boundary conditions, such as penetration problems including the installation of
piles. The application of the ALE method in the dynamic analysis of geotechnical problems was
also addressed by Nazem et al. (2009a).
2.2
Operator-split technique
The Arbitrary Lagrangian-Eulerian (ALE) method has been developed to eliminate the mesh
distortion that usually occurs in the Updated Lagrangian (UL) method by separating the material
and mesh displacements. In a UL formulation, all variables are calculated at the end of the last
equilibrium configuration. This assumption necessitates updating the spatial coordinates of all
material points according to the incremental displacements at the end of each time step. While
capable of solving some problems with relatively large deformations, the UL method often fails to
succeed due to the occurrence of excessive mesh distortion and the entanglement of elements.
Mesh distortion usually starts with elements twisting or otherwise distorting out of their favourable
shapes. For instance, in a domain consisting of triangular elements it is usually preferred that all
triangles remain roughly equilateral during the analysis. However, material displacements,
especially in regions with higher deformation gradients, will gradually cause distortion of elements,
decreasing the accuracy of the analysis or even ultimately resulting in a negative Jacobian. One
way to overcome this shortcoming of the UL method is to separate the material and mesh
displacements, leading to the development of the ALE method.
The equilibrium equation in the ALE method can be written in two different forms. It is
possible to write the governing global system of equations in terms of two sets of unknown mesh
and material displacements, leading to the so-called coupled ALE method. A supplementary set of
equations in terms of the material and mesh displacements needs to be established through a mesh
motion scheme and the two sets of unknown displacements are then solved simultaneously.
Alternatively, in the decoupled ALE method or the operator-split technique, the analysis can be
performed in two separate steps: a UL step followed by an Eulerian step. In the UL step, the
governing equations are only solved for the material displacements in order to fulfil the
requirements for equilibrium. This step usually results in a distorted mesh. In the subsequent
Eulerian step, the main goal is to minimise the mesh distortion by refining the mesh. Mesh
refinement can be achieved by generating a new mesh for the entire domain or by moving current
nodal points into new positions. Regardless of which strategy is adopted for mesh refinement, the
topology of the domain, i.e., the number of nodes, number of elements, and the connectivity of
elements should not be changed. The Eulerian step is particularly important if significant mesh
distortion occurs during the UL step.
After mesh refinement, all state variables at nodal points as well as at Gauss points must be
mapped from the distorted mesh to the new mesh. This remapping is usually performed using a
first order expansion of Taylor’s series as
∂f
fr = f + vi − vir ⋅
∂xi
(
)
(1)
and f denote the time derivatives of an arbitrary function f with respect to the mesh
where f
r
and material coordinates respectively, vi is the material velocity, and vir represents the mesh
velocity. Other important aspects of the ALE method are discussed as follows.
2.3
Updated Lagrangian formulation
During each increment of the operator split technique, the analysis starts with a UL procedure at
time t. The main goal of the UL step is to find the incremental displacements, velocities and
accelerations that satisfy equilibrium at time t+Δt. The matrix form of equilibrium is usually
derived from the principle of virtual work. The weak form of this principle states that for any
virtual displacement δu, equilibrium is achieved provided
" − ∫ σ δε dV − ∫ δ u ρu dV −
ij
ij
k
i
i
k
V
V
∑$$$+ k δ u b dV + δk u f dS
k
# ∫Vk i i k ∫ Sk i i k
+∫
Sc
(t
N
∫
Vk
δ ui cui dVk %
'
''
&
δ g N + tT δ gT )dSc = 0
(2)
where k is the total number of bodies in contact, σ denotes the Cauchy stress tensor, δε is the
variation of strain due to virtual displacement, u, u and u represent material displacements,
velocities and accelerations, respectively, ρ and c are the material density and damping, b is the
body force, q is the surface load acting on area S of volume V, δgN and δgT are the virtual normal
and tangential gap displacements, tN and tT denote the normal and tangential forces at the contact
surface Sc. The solution of Equation (2) requires the discretisation of the domain as well as the
contact surfaces. In this study we adopt the so-called node-to-segment contact discretisation as
shown in Figure 1, in which a node on the slave surface may come into contact with an arbitrary
segment of the master contact surface. With this assumption the last two terms in Equation (2),
representing the virtual work due to normal and tangential contact forces, can be written in the
following form:
Master surface
Slave node
Active master
segment
Figure 1. Node-to-segment contact discretisation
∫
Sc
∫
nc
t N δ g N dSc ≈ ∑ δ uiT FNci
i =1
t δ g N dSc ≈ ∑ δ uiT FNci
Sc N
(3)
nc
i =1
(4)
in which nc is the total number of slave nodes, and FNc and FTc represent the normal and tangential
nodal forces of the contact element, respectively.
The equation of motion for a solid can be obtained by linearisation of the principle of virtual
displacement and is written in the following matrix form (e.g., Wriggers 2006)
t+Δt
t+Δt + Cu t+Δt + R t+Δt = Fext
Mu
(5)
where M and C represent respectively the mass and damping matrices, R is the stress divergence
term, u denotes the displacement vector, and Fext is the time-dependent external force vector. Note
that the right-hand super-script denotes the time when the quantities are measured and a
superimposed dot represents the time derivative of a variable. The solution of the momentum
equation requires a step-by-step integration scheme in the time domain. Explicit and implicit
methods are available for this purpose. In explicit methods the solution at time t+Δt depends only
upon known variables at time t. These methods, such as the central difference method, do not
require a factorization of the effective stiffness matrix and are easy to implement. However,
explicit methods are conditionally stable, i.e., the size of the time steps must be smaller than some
critical time step. Implicit methods, on the other hand, require the solution of a nonlinear equation
at each time step and have to be combined with another procedure such as the Newton-Raphson
method. The main advantage of implicit methods is that they can be formulated to be
unconditionally stable, allowing the analyst to use a bigger time step than used in the explicit
methods.
In this approach we adopt the implicit generalised-α method presented by Chung and Hulbert
(1993) to integrate the momentum equation over the time domain. This method has been shown to
be very efficient in solving dynamic problems in geomechanics (e.g., Kontoe et al. 2008 and
Nazem et al. 2009a). In the generalised-α method, two integration parameters α f and α m , are
introduced into the momentum equation to compute the inertia forces at time t + (1 − α m ) Δt and the
internal and damping forces at time t + (1 − α f ) Δt , respectively. The accelerations and velocities
are approximated by Newmark equations according to
t+Δt =
u
1
1 t 1− 2 β t

ut+Δt − ut −
u −
u
2
βΔt
2β
βΔt
u t+Δt =
⎛ α⎞
⎛
α
α ⎞ t

ut+Δt − ut + ⎜ 1− ⎟ u t + ⎜ 1− ⎟ u
βΔt
⎝ β⎠
⎝ 2β ⎠
(
(
)
)
(6)
where α and β are the Newmark integration parameters and Δt represents the time step.
By introducing a tangential stiffness matrix
⎛ ∂R ⎞
K T ( uit +Δt ) = ⎜ t +Δt ⎟
⎝ ∂u ⎠uti +Δt
(7)
we can use the Newton-Raphson method to solve the momentum equation in the following form
for the ith iteration
$
'
α 1− α f
& 1− α m ⋅ M +
⋅ C + 1− α f ⋅ K T (i−1) ) ⋅ Δu(i ) =
&% β ⋅ Δt 2
)(
β ⋅ Δt
(
t+Δt
ext
(1− α ) ⋅ F
f
)
(
)
t
t+Δt
+ α f ⋅ Fext
− 1− α f ⋅ Fint
(i−1)
(
)
$ 1− α
1− α m t '
t+Δt
t
m
⋅
u
−
u
−
⋅ u )
&
2
(i−1)
β ⋅ Δt
& β ⋅ Δt
)
t
−α f ⋅ K T (i−1) ⋅ u − M &
)
* 1− α
- t
m
&− ,
)

−1/⋅ u
)(
.
%& + 2β
(
)
$ α 1− α
* α
- t'
f
t
&
⋅ u(t+Δt
−
u
+
1−
1−
α
⋅ u )
,
f /
i−1)
& β ⋅ Δt
+ β
. )
−C&
)
t
&−Δt *, α −1-/ 1− α ⋅ u
)

f
&%
)(
+ 2β .
(
)
(
)
(
(
)
u(t+Δt
= u(t+Δt
+ Δu(i )
i)
i−1)
where
)
,
u(t+Δt
= ut
0)
(8)
and where Fint , the internal force vector, is obtained according to the Cauchy stress tensor, σ, and
the nodal forces at contact surfaces from
t +Δt
Fint
=∫
V
nc
(
BT ⋅ σ t +Δt dV t +Δt − ∑ FNc i + FTci
t +Δt
i =1
)
(9)
Note that the tangential stiffness matrix in Equation (7) includes the contribution of the material
stiffness, Kep, the stiffness due to geometrical nonlinearity, Knl, and the stiffness due to normal and
tangential contact, KNs and KTs, i.e.,
KT = K ep + K nl + K Ns + KTs
(10)
In the problems studied here the penetrometer is idealised as a rigid body, i.e., its size and shape
remain unchanged during penetration. To achieve this condition all nodal points on the
penetrometer are prescribed to undergo equal vertical displacements while horizontal displacement
is prohibited.
2.4
Definition of contact
To describe the contact at the interface between two bodies, constitutive equations must be
provided for both the tangential and normal directions. Among several strategies available in
contact mechanics we use the penalty method to formulate these constitutive relations. The normal
contact in the penalty method is described by
tN = ε N g N
(11)
where εN is a penalty parameter for the normal contact.
In general, the response in the tangential direction is captured by the so-called stick and slip
actions. In the former, no tangential relative movement occurs between the bodies whilst the latter
represents a relative displacement of bodies in the tangential direction. This assumption facilitates
splitting the relative tangential velocity between the bodies, g T , into a stick part gTst and a slip gTsl
part, as in the following rate form
gT = gTst + gTsl
The stick part can be used to obtain the tangential component of traction by
(12)
tT = ε T gT
(13)
where εT is a penalty parameter for tangential contact. Note that the assumption in (12) is
analogous to the theory of plasticity in which the incremental strains are divided into an elastic and
a plastic part. Following this analogy, we must provide a slip criterion function. This can be
achieved by writing the classical Coulomb friction criterion in the following form
f s ( t N , tT ) = t T − µt N ≤ 0
(14)
where µ is the friction coefficient.
For the numerical examples considered later in this paper, the interface between the
penetrometer and the soil is assumed, for simplicity, to be perfectly smooth. This situation
corresponds to the special case where µ = 0.
2.5
Material behaviour
Undrained behaviour of the soil to be penetrated is represented by an elastoplastic Tresca material
model with an associated flow rule or by the Modified Cam Clay (MCC) stress-strain model, which
is formulated in terms of effective stresses. For cases involving the fully drained response of the
soil the MCC model has been adopted.
In a large deformation analysis, the stress-strain relations must be frame independent. This
requirement, known as the principle of objectivity, can be satisfied by introducing an objective
stress rate, such as the Jaumann stress rate σ ∇J , into the stress-strain relations as
ep
dσ ij∇J = Cijkl
⋅ d ε kl
(15)
ep
where C represents the elastoplastic constitutive matrix and ε is the strain vector. Among several
available options to integrate equation (15), Nazem et al. (2009b) showed that it is slightly more
efficient to apply rigid body corrections during integration of the constitutive equations. This
strategy is adopted in this study.
It is noted that the undrained shear strength of cohesive soils often increases with the rate of
straining. This phenomenon, termed the strain rate effect, usually introduces uncertainties in
predicting the shear strength of soils in dynamic penetration tests. Graham et al. (1983) proposed a
simple equation which prescribes approximately 5-20% increase in shear strength for each order of
magnitude increase in the rate of shear strain. This effect can be expressed by the following widely
used equation
⎡
⎛ γ
su = su,ref ⎢1+ λ log ⎜
⎢⎣
⎝ γ ref
⎞⎤
⎟⎥
⎠ ⎥⎦
(16)
where su is the undrained shear strength of the clay, su,ref is a reference undrained shear strength
measured at a reference strain rate of γ ref , γ denotes the shear strain rate and λ is the rate of
increase of strength per log cycle of time. Typically γ ref = 0.01 per hour for clays (Einav and
Randolph 2006), and this value was adopted in the current study. For simplicity, the shear strength
of a material can be assumed to be constant while integrating the constitutive equations during an
individual time-step of a nonlinear finite element analysis. However, the shear strength parameters
must be updated according to the known shear strain rates at the end of each increment. Note that
the parameter γ appearing in equation (16) should be evaluated using a frame-independent
quantity such as the maximum plastic shear strain.
In cone penetration tests it is found that the penetration resistance in clays increases as the rate
of penetration (v) or the diameter of the penetrometer decreases. In static tests this is mainly due to
partial consolidation occurring in advance of the cone tip, which in turn increases the effective
stresses around the cone followed by an increase of the sleeve friction. Usually a normalised
velocity, V, is used to assess this degree of consolidation according to
V=
vd
cv
(17)
where cv represents the coefficient of consolidation of the soil (e.g., Lehane et al. 2009). Values of
the normalised velocity between 10 (O'Loughlin et al. 2004) and 100 (Brown and Hyde, 2008)
have been recommended to ensure that the soil behaviour is undrained. In the penetration problems
studied here the value of normalised velocity is generally above 10,000 due to the fast penetration
of the object (see Section 4), and the relatively short penetration time (0.1-0.9 s). Therefore, we
assume undrained conditions throughout the penetration, and thus neglect the shear strength
increase due to partial consolidation.
Centre
line
FFP
Soil
(a) Boundary
between soil and
object before
deformation.
(b) Distorted
boundary
(c) Boundary after
nodal relocation
Figure 2. Nodal relocation on the boundary
2.6
Mesh refinement
In the ALE operator split technique the material displacements are obtained at the end of the UL
step and will normally result in a distorted mesh. Refining the mesh at the beginning of each Euler
step is very important since a distorted mesh can lead to inaccurate results. Most mesh refinement
techniques are based on special mesh-generation algorithms, which must consider various
parameters such as the dimensions of the problem, the type of elements to be generated and the
regularity of the domain. Developing such algorithms for any arbitrary domain is usually difficult
and costly. Moreover, these algorithms do not necessarily preserve the number of nodes and the
number of elements in a mesh and they may cause significant changes in the topology. A robust
method for mesh refinement, based on the use of a simple elastic analysis, was presented by Nazem
et al. (2006). The method has been implemented for two-dimensional plane strain problems as well
as axi-symmetric problems. To obtain the mesh displacements, we first re-discretise all the
boundaries of the problem, which include the boundaries of each discrete body, the material
interfaces and the loading boundaries, resulting in prescribed values of the mesh displacements for
the nodes on these boundaries. Each boundary node is then relocated along the boundary as
necessary. With the known total displacements of these nodes on the boundaries, an elastic
analysis is then performed using the prescribed boundary displacements to obtain the mesh
displacements for all the internal nodes and hence the optimal new mesh. An important advantage
of this mesh optimisation method is its independence of element topology and problem dimensions.
The method does not require any mesh generation algorithm, does not change the topology of the
problem, and hence can be easily implemented in existing finite element codes. This method is
potentially a good candidate for three-dimensional mesh optimisation.
The relocation of nodes along the boundaries plays a key role in the mesh refinement scheme
proposed by Nazem et al. (2006).
This relocation demands an efficient mathematical
representation of the boundaries, which strongly depends on the physical problem being simulated.
Nazem et al. (2008) proposed an efficient method for nodal relocation based on a quadratic spline
technique. This method was successfully employed for solving many geotechnical problems such
as consolidation and bearing capacity of soils under a footing. However, it is not suitable for
simulating boundaries in a FFP problem due to the linear shape of the penetrating object, which
consequently introduces abrupt changes in the slope of the straight lines defining the boundary.
Instead, we introduce a simpler, but more efficient, algorithm for nodal relocation along the
boundary of the soil in contact with the penetrometer, which is based on the assumption of linear
boundary segments. This procedure is depicted schematically in Figure 2. The boundary between
the soil and the penetrating object at the beginning of a time step (Figure 2a) may become distorted
after the UL step, as shown in Figure 2b. Such cumulative distortion, if not prevented, can
potentially cause significant distortion of elements in that region ultimately indicated by a negative
Jacobian. Thus, nodal relocation on the boundary can be performed as shown in Figure 2c. Note
that the nodes on a curved boundary are only permitted to move in the tangential direction of the
boundary, i.e., the normal component of the convective velocity on a boundary is zero, but not
necessarily the tangential component. However, on a multi-linear boundary the nodes may be
allowed to move along a line defining a linear segment. The nodal relocation on linear segments is
shown schematically in Figure 3, where node i must be moved into a new position, i′ , along line 1
or line 2 such that l = li . Note that li represents a normalised length of the segment and can be
calculated at the beginning of a time step before the mesh is distorted. The Cartesian coordinates
of the new location can be calculated explicitly. For brevity, the full mathematical details are not
provided here but the procedure for relocating all nodes along a boundary is described in Nazem et
al. (2012).
2.7
Energy absorbing boundaries
One of the well-known issues in computational dynamic analysis of soil-structure-interaction (SSI)
problems is how to simulate an infinite medium. Employing the finite element method, for
instance, with boundaries that are not infinitely distant, one must guarantee that the outgoing waves
from the source (usually a structure) do not reflect back from the finite boundaries toward the
source, since in reality these waves should propagate to infinity and dissipate at a far distance from
the source. If it is allowed to occur, such reflection will most probably affect the accuracy of the
numerical results. To assure no waves are reflected back from truncated computational boundaries,
it is common to use artificial boundaries that absorb the energy of incoming waves.
A simple, but efficient, boundary was developed by Lysmer and Kuhlemeyer (1969), which is
known as the “standard viscous boundary” in the literature. The standard viscous boundary is
probably the most popular artificial boundary since it possesses an acceptable dissipation
characteristic at a low computational cost. Kellezi (1998) suggests that absorbing boundaries must
not be located closer than (1.2 - 1.5)λs (where λs is the length of the shear waves) from the
excitation source. A recent study by Kontoe (2006) showed that the standard viscous boundary is
capable of absorbing dilatational waves (P-waves) as well as shear waves (S-waves) in the analysis
of plane strain and axisymmetric problems. When required, this type of boundary was used in the
problems solved in this study.
y
l=li
i’
i
i+1
i-1
x
Figure 3. Linear nodal relocation
3
PENETRATION INTO NON-UNIFORM SOIL
In previous studies (e.g., Nazem et al. 2012), we considered a smooth free-falling object
penetrating into a uniform soil layer, in which the material properties do not change in any
direction. However, in many natural settings the shear strength of soils actually increases with
depth. This variation can be expressed by a linear equation as
su ,ref ( y ) = su ,ref ( 0 ) + ks y
(18)
where su,ref(y) represents the reference undrained shear strength of soil at depth y, ks indicates the
increase in shear strength per unit depth, and su,ref(0) is the reference shear strength of soil at the
ground surface.
Table 1. Parameter values adopted in Section 3
Parameter
Unit
Value
d
M
vo
su,ref(0)
G/su,ref
L
ksd/su,ref(0)
m
kg
m/sec
kPa
0.04, 0.06,0.08, 0.1
0.5, 1.0
5, 10
1.0
67
0.2
0, 0.25, 0.5, 1.0
-
The effect of the shear strength increase parameter ks, on soil penetration was studied by
analyzing more than 50 different problems with parameters summarized in Table 1. Note that in all
cases the ratio between the shear modulus of the soil, G, and its undrained shear strength, su,ref, is
assumed to be 67 at all integration points, and the rate parameter, λ, is 0.2. Depending on the site,
the rate of strength increase with depth ks, usually varies between 1-3 kPa/m. Adopting these
typical values for a free-falling penetrometer (FFP), the effect of ks on the penetration
characteristics will not be significant, mainly due to the small diameter (size) of the penetrating
object. Alternatively, a normalized parameter, ksd/su,ref(0), is adopted in this study, which quantifies
the significance of the rate of shear strength increase on the penetration.
The finite element mesh and boundary conditions of this axi-symmetric problem are depicted in
Figure 4. The finite element mesh includes 10,252 nodal points and 4,988 6-node triangular
elements each containing 6 integration points. For simplicity, the material damping and friction
between the soil and penetrometer are assumed to be zero. Note that the right-hand and bottom
boundaries are able to dissipate the energy of oncoming waves.
m, v0
Smooth energy absorbing boundary
10d
16d
Smooth energy absorbing boundary
7.5d
Figure 4. Finite element mesh of free falling penetrometer
In a previous study, Nazem et al. (2012) showed that an appropriate strategy for studying the
dynamic penetration problem is to plot the initial impact velocity of the penetrometer, 0.5mv2,
normalised by 0.25πd3su,ref(0), versus the total depth of penetration, normalized by d. Adopting this
strategy, the normalised impact energy of the penetrometer is plotted versus the normalised depth
of penetration in Figure 5. According to Figure 5, the depth of penetration of objects with an
identical impact energy decreases as the normalised parameter ksd/su,ref(0) increases. Figure 5 also
shows that there is a linear relation between the normalised depth of penetration and the normalised
Normalised impact energy
impact energy, but the slope of such lines depends on the value of ksd/su,ref(0). This dependency,
for the specimens analysed in this example, can be eliminated by plotting the initial impact energy,
normalised by 0.25πd3su,ref(d), versus the normalised depth of penetration, as represented in Figure
6. Note that su,ref(d) indicates the reference shear strength of the soil at depth d. Similar
independence for other material properties and penetrometer parameters is yet to be investigated.
ksd/su,ref = 0.0
ksd/su,ref = 0.25
ksd/su,ref = 0.5
ksd/su,ref = 1.0
Penetration in diameters
Normalised impact energy
Figure 5. Initial impact velocity of the penetrometer, normalised by 0.25πd3su,ref(0), versus the total depth
of penetration, normalised by d
Penetration in diameters
Figure 6. Initial impact velocity of the penetrometer, normalised by 0.25πd3su,ref(d), versus the total depth
of penetration, normalised by d
To study the deceleration characteristics of the penetrometer, its variation of velocity versus
penetration is plotted in Figure 7 for specimens where m = 0.5 kg, d = 0.04 m, v0 = 10 m/sec, and
ksd/su,ref(0) = 0.0, 0.25, 0.5, and 1.0. The plots in Figure 7 indicate that the velocity of the object
reduces approximately quadratically as it decelerates and penetrates further into the non-uniform
layer of soil. Note that the quadratic reduction of velocity with penetration has also been reported
in experimental studies as well as numerical analysis results of penetration into a uniform clay
layer (Nazem et al. 2012). It is also notable that for all cases studied here the total time of
penetration is less than 0.3 sec.
Velocity (m/sec)
Penetration in diameters
ksd/su,ref = 0.0
ksd/su,ref = 0.25
ksd/su,ref = 0.5
ksd/su,ref = 1.0
Figure 7. The velocity of the penetrometer versus the penetration in diameters
4
TORPEDO ANCHORS
An interesting application of dynamic penetration problems is the numerical analysis of the
installation of torpedo anchors, which have proved to be efficient for deep-water anchoring systems
due to their relatively easy and economic installation. Torpedo anchors are usually 10-15 m long,
0.75-1.2 m in diameter, and may weigh between 250-1,000 kN (Randolph & Gourvenec 2011).
Assuming that the installation phase is relatively fast, we consider only undrained soil behaviour
here. For simplicity, material damping and the friction between the soil and anchor are ignored.
The impact velocity of the anchor is assumed to be less than its terminal velocity in a fluid, i.e., the
torpedo may accelerate after impacting the soil due to the gravitational force.
To investigate the soil behaviour during torpedo penetration, 27 specimens were analysed with
the properties summarised in Table 2. The finite element mesh and the boundary conditions are the
same as those used in the previous example, and are shown in Figure 4.
Table 2. Values of parameters used in torpedo anchor analysis
Parameter
Unit
Value
d
m
vo
su,ref (0)
G/su,ref
λ
ksd/su,ref (0)
m
kg
m/sec
kPa
-
1.0
1,000
5, 10
1.0
67
0.25, 0.5
0, 1.0, 2.5
The numerical results for all specimens, in terms of the total depth of penetration, p, and the
total time of penetration, tp, are summarised in Table 3. According to Table 3 the total penetration
time of a torpedo anchor is comparatively higher than the penetration time of a miniature freefalling penetrometer (FFP). In addition, the strength rate and shear strength increase with depth
parameters, λ and ks, significantly influence the total depth and time of penetration.
Table 3. Summary of numerical results for torpedo anchors
No
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
ks (kPa/m)
0
0
1
1
1
1
1
1
2.5
2.5
2.5
2.5
2.5
2.5
0
0
1
1
1
1
1
2.5
2.5
2.5
2.5
2.5
2.5
λ
0.25
0.25
0.25
0.25
0.25
0.25
0.25
0.25
0.25
0.25
0.25
0.25
0.25
0.25
0.5
0.5
0.5
0.5
0.5
0.5
0.5
0.5
0.5
0.5
0.5
0.5
0.5
m (kg)
2,000
2,000
2,000
2,000
5,000
5,000
10,000
10,000
2,000
2,000
5,000
5,000
10,000
10,000
2,000
2,000
2,000
5,000
5,000
10,000
10,000
2,000
2,000
5,000
5,000
10,000
10,000
v0 (m/sec)
5
10
5
10
5
10
5
10
5
10
5
10
5
10
5
10
10
5
10
5
10
5
10
5
10
5
10
p (m)
5.34
14.34
1.96
3.19
4.15
6.07
7.67
10.24
1.49
2.31
2.69
3.87
4.24
5.87
2.21
6.93
2.46
3.01
4.71
5.56
7.93
1.22
1.85
2.06
3.11
3.23
4.66
tp (sec)
1.85
2.89
0.53
0.50
0.93
0.83
1.43
1.26
0.38
0.35
0.63
0.54
0.88
0.77
0.71
1.41
0.39
0.74
0.66
1.16
1.03
0.31
0.28
0.50
0.44
0.72
0.62
To study the penetration properties of a torpedo anchor in more detail, the velocity variations of
Specimens 17 and 23 are plotted versus time as well as the normalised penetration in Figures 8 and
9, respectively. For these two specimens all input parameters, except ks, are identical. According
to Table 3, the value of ks in Specimens 17 and 23 are 1.0 and 2.5 kPa/m, respectively. Figures 8
and 9 indicate that the torpedo accelerates at the early stages of penetration. This behaviour was
observed in all specimens, indicating that, unlike miniature FFPs, the acceleration of the torpedo
changes approximately linearly with time.
To study the soil behaviour due to torpedo penetration we consider two different cases where
ks = 0 (uniform soil) and ks > 0 (non-uniform). For a uniform soil, Nazem et al. (2012) showed that
the dynamic soil resistance due to penetration of a FFP can be estimated according to:
qdyn = Ndp su ,ref ( 0)
(19)
where
⎛ G
N dp = 2.321 − 50.488λ + (1.690 + 34.546λ ) ln ⎜
⎜s
⎝ u ,ref
⎛ p⎞
⎜ ⎟
⎝d⎠
−1
( −2.698 − 148.36λ + 410.27λ
⎛
( 2.399 + 46.52λ − 103.91λ ) ln ⎜⎜ s G
2
2
⎝
u , ref
⎞
⎟⎟ −
⎠
+
⎞⎞
⎟⎟ ⎟
⎟
⎠⎠
(20)
Velocity (m/sec)
Time (sec)
Figure 8. The velocity of a torpedo anchor versus time
Velocity (m/sec)
Penetration in diameters
Figure 9. The velocity of torpedo anchor versus the penetration normalised by diameter
Using the values in Table 3, the dynamic penetration factors, Ndp, for specimens 2 and 16 are
estimated as 31.75 and 53.01, respectively. On the other hand, the dynamic soil resistance versus
the normalised penetration predicted by numerical analysis is plotted in Figure 10 for the same
specimens. According to Figure 10, the predicted normalised dynamic soil resistance is in good
agreement with the value predicted by Equation (19) for Specimen 2, and is slightly lower than
53.01 for Specimen 16. This suggests that the validation of Equations (19) and (20) for torpedo
anchors penetration into a uniform soil layer requires further investigation.
Unlike uniform soils, the dynamic resistance of a non-uniform layer of soil is not likely to
converge to a constant value as the normalised penetration, p/d, approaches infinity. For
Dynamic soil resistance / su,ref (0)
Specimens 8 and 14, where ks is respectively 1.0 and 2.5 kPa/m, the normalised dynamic soil
resistance is plotted versus the normalised penetration of the torpedo in Figure 11. It is evident that
the soil resistance due to penetration increases approximately linearly with penetration. This
observation also indicates that the force acting on the anchor due to soil resistance changes with
time as the object penetrates into soil, and hence the deceleration of torpedo is not constant during
penetration.
Penetration in diameters
Dynamic soil resistance / su,ref (0)
Figure 10. Normalised dynamic soil resistance versus normalised penetration of torpedo into a nonuniform soil layer
Penetration in diameters
Figure 11. Normalised dynamic soil resistance versus normalised penetration of torpedo into a uniform
soil layer
5
VALIDATION OF PENETRATION PREDICTIONS
Previously, Carter et al. (2010) reported the validation of the finite element results predicted by the
method of penetration analysis by comparing the numerical predictions with experimental test
results undertaken at the University of Sydney, Australia (Chow 2012). These tests involved model
penetrometers being dropped into pots of normally consolidated clay with constant undrained shear
strength. An illustration of the experimental apparatus used in these tests is presented in Figure 12.
Some of the results of the testing program conducted to validate the numerical approach are
presented in Table 4. This table also provides the results obtained from the ALE finite element
analysis, which were all conducted assuming λ = 0.2. Table 4 shows that there is generally good
agreement between the experimental results and the numerical predictions.
Adjustable
release
mechanism
Wall
400
2450
Penetrometer
315
Kaolin
pot
Note : All dimensions in mm
Drawing not to scale
(a)
(b)
Figure 12. (a) Side elevation of laboratory free falling test setup; (b) A completed free falling
penetrometer drop into the kaolin pot
6
USE OF PENETRATION PREDICTIONS TO INTERPRET IN SITU TESTS
Some examples of the predictions of the final penetration of an object free falling into a uniform
clay deposit obtained using the dynamic ALE method were presented in the previous sections. For
relatively uniform cohesive soil deposits deforming under undrained conditions, four material
parameters are required to define completely the response of the soil to penetration. These are G,
su,ref, γ ̇ref, and λ.
Table 4. Comparison between experimental results and ALE predictions
Test
su,ref
(kPa)
d
(m)
m
(kg)
vo
(m/s)
p/d
(Measured)
p/d
(ALE)
Difference ALE and
tests (%)
1
5.15
0.02
0.26
4.77
4.33
3.99
7.7
2
5.15
0.02
0.35
4.77
4.80
5.22
8.7
3
5.15
0.02
0.45
4.75
5.20
6.50
25.0
4
5.15
0.02
0.54
4.74
7.05
7.80
10.6
5
5.15
0.02
0.63
4.76
7.94
9.36
17.9
6
7.46
0.04
0.71
4.77
1.14
1.39
22.2
7
6.91
0.04
0.71
4.75
1.39
1.44
3.7
8
6.91
0.03
0.74
4.75
2.75
2.71
1.4
9
4.45
0.02
0.26
4.75
4.45
4.38
1.6
For undrained deformations Poisson’s ratio of the soil takes the value 0.5 because the soil
deforms at constant volume. It was found from this study that the final penetration depth was
relatively insensitive to G, but was very dependent on the uniform strength su,ref, somewhat
dependent on the rate parameter λ and less dependent on the rigidity index, G/su,ref. Provided
reasonable assumptions can be made for the values of the parameters G/su,ref, γ ̇ref and λ, the final
penetration depth can be used to estimate the uniform value of the rate-independent undrained
strength, su,ref. Hence the normalisation of the numerical predictions shown, for example, in Figure
6 proves a useful means of interpreting that value of undrained shear strength from the measured
final penetration of the object into the soil deposit. This process is illustrated schematically in
Figure 13.
Figure 13. Schematic illustraion of the interpretation method for the falling penetrometer test
7
ANALYSIS OF IMPACT COMPACTION
In this example numerical analysis of the dynamic compaction of a soil layer under a free-falling
weight is considered. The problem domain, axi-symmetric finite element mesh, boundary
conditions, and material properties are shown in Figure 14. The water table is located at the ground
surface. A fully coupled consolidation analysis under dynamic loading conditions has been carried
out. The Modified Cam-Clay (MCC) material model, which facilitates plastic volumetric changes,
was used to model soil stress-strain behaviour. The material parameters shown in Figure 14 are
defined as follows: λ = the slope of the normal compression line (NCL) in the space of the
logarithmic mean stress lnp′ versus the void ratio e, κ = the slope of the unloading-reloading line
(URL) in the lnp' - e space, e0 = the initial void ratio, OCR = the over-consolidation ratio of the
soil, φ = the friction angle of soil, γ = the unit weight of the soil, k=the hydraulic conductivity of the
soil, and ν = Poisson’s ratio of the soil.
Smooth / impermeable
8R
Smooth / permeable
Smooth / impermeable
R=1 m
γ=17 kN/m3
ν=0.3
φ=250
λ=0.25
κ=0.05
e0=1.8
k=10-9 m/s
OCR=1
Smooth / impermeable
8R
Figure 14. Dynamic compaction problem
The finite element mesh includes 1447 nodal points and 690 6-node triangular elements. The
analysis includes three stages. In the first stage, a body force loading due to the self-weight of the
soil was used to establish an initial stress field over a long period of time, allowing all excess pore
pressure to dissipate. After generating the nonzero stress field, the location of the yield surface at
each integration point was adjusted according to the initial effective stresses and the designated
value of the overconsolidation ratio (OCR). The second stage of analysis included applying an
overburden pressure of p0 = 5 kPa. Finally, the weight was dropped onto the soil layer at an initial
velocity, and removed as soon as it comes to rest. The process of releasing the falling weight and
allowing it to impact the soil surface was repeated in quick succession, a further four times, thus
applying 5 dynamic impact blows to the surface of the soil layer. Four different values of the
impact velocity, viz., 1, 2, 3 and 4 m/s, were used in the numerical simulations. To avoid further
complexities, the object was assumed to be rough so that tangential displacement of the soil at the
falling weight-soil interface was not permitted.
It should be noted that the value selected for the hydraulic conductivity of the soil (10-9 m/sec)
is typical of a fine-grained soil such as clay. Each successive blow from the falling weight occurs
almost immediately after the weight comes to rest during the previous blow. Furthermore, the
period of the each impact is only a fraction of a second, so that the overall period to effect 5 blows
is very short, on the order of 1 sec. Thus, in this particular example the behaviour of the soil is
effectively undrained, i.e., it will occur under the constraint of constant volume deformation.
Some results of these dynamic consolidation analyses are presented in Figures 15 to 18. Figure
Time (sec)
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
Velocity (m/sec)
0
-0.5
Impact
velocity = 1
m/sec
-1
-1.5
Impact
velocity = 2
m/sec
-2
Figure 15. Velocity versus time
Settlement / R
0
0.1
0.2
0.3
0.4
0.5
Velocity (m/sec)
0
-0.5
Impact
velocity = 1
m/sec
-1
-1.5
Impact
velocity = 2
m/sec
-2
Figure 16. Velocity versus normalized settlement
15 indicates the velocity of the falling weight from the time it first impacts the soil surface until it
comes to rest in the soil. After each blow it is raised virtually instantaneously and allowed to
impact the soil again with the same velocity as the previous blow. Results are included in Figure
15 for initial impact velocities of 1 and 2 m/sec. The same results are presented in Figure 16 in the
form of the velocity of the falling weight versus the settlement it causes in the soil immediately in
Dynamic soil resistance (kPa)
70
Impact velocity = 1 m/sec
60
Impact velocity = 2 m/sec
50
40
30
20
10
0
0
0.1
0.2
0.3
0.4
Settlement / R
Figure 17. Dynamic soil resistance versus normalised settlement
contact with the falling weight. Both Figures 15 and 16 indicate the cumulative effects of each
successive impact. The final settlement after 5 blows is dependent on the initial impact velocity,
and for the two cases considered here is either almost 20% or 40% of the radius of the impact
weight.
The variation of the soil resistance with penetration of the falling weight into the soil, for each
impact blow, is illustrated in Figure 17. The resistance is plotted as the average contact pressure
between the soil and the falling weight. It is of interest to observe that the peak in this resistance
increases with each blow. It would also be of interest to predict whether this increase continues
with blows beyond the mere 5 studied here. However, such calculations have not yet been
conducted and are the subject of an on-going research work.
Because the impact loading occurs relatively rapidly on this saturated clay-like soil, it is also of
interest to predict the excess pore water pressures generated in the soil and how they vary as a
result of repeated impact blows. Figure 18 indicates contours of total pore water pressure (static
plus excess) initially, before any impact loading, and immediately after each impact blow from the
falling weight, for the case of an initial impact velocity of 2 m/sec. The steady accumulation of
excess pore water pressure with each blow is indicated by the regular and increasing distortion of
the pattern of initial static pore water pressure. The increasing penetration of the falling weight
with successive blows is also evident from the deformed shapes plotted in Figure 18. It would also
be of interest to predict the further deformation of the soil as these excess pore water pore water
pressures generated by the repeated impact dissipate and the soil consolidates. Such predictions are
also the subject of on-going work.
8
EXPERIMENTAL RESULTS FOR IMPACT COMPACTION
To enable comparison with the numerical results model scale dynamic compaction tests have been
conducted. These have been performed in a 2-D configuration in a tank with a Perspex face,
shown schematically in Figure 19, allowing observation of the mechanism of dynamic compaction.
These tests have been limited to dry and partially saturated sands and sand-silt mixtures and thus
are not directly comparable to the numerical analyses of the impacts on clay.
The tests have made use of high speed photography and digital image correlation (DIC)
techniques to investigate the deformation patterns, calculate soil strains and observe strain
localisation. Additional instrumentation has included accelerometers on the plunger and buried in
the soil, and load cells. The testing apparatus, shown in Figure 19, consists of a fabricated
container fixed to a steel frame open at the top with a clear acrylic (Perspex) front. The metal
a.Before impact
c. After 2nd drop
e. After 4th drop
b. After 1st drop.
d. After 3rd drop
f. After 5th drop
Figure 18. Contours of pore water pressure (kPa) – Impact velocity = 2 m/sec
frame extends to a height of 2.8 m above the floor with a fixed pulley system resting on an axle
positioned at the top of the frame. The internal dimensions of the soil container are 350 mm x
150 mm x 750 mm. A 25 mm thick Perspex sheet is bolted to the front of the frame. A high speed
camera was positioned 1.2 m from the testing apparatus. The field of view of the 1024 x 1024
pixel image is 403 mm x 403 mm giving an object space pixel size of 0.39 mm. Photographs have
been taken at rates between 1000 frames per second (fps) and 3000 fps. Further details of the
photographic technique are provided in Nazhat et al (2011).
Figure 19. Schematic of test configuration
The guiding trolley and attached plunger are able to free fall up to 1 m, although in the tests
discussed here drop heights ranged between 150 mm and 300 mm (impact velocities of 1.7 to
2.4 m/s). The weight of the free falling assembly ranged from 50 N to 180 N, depending on the
width of the tamper which varied between 35 mm and 110 mm. The tamper length (150 mm) was
equal to that of the container to ensure a 2-D configuration.
A typical test has involved 5 to 8 drops of the tamper, with high speed photography taken over a
time interval of 3.2 seconds for each drop. The test programme has included tests of dry Sydney
sand at different relative densities, tests to investigate the effect of adding non-plastic fines to the
sand, tests to investigate the effect of degree of saturation, and tests to investigate the effect of
tamper geometry.
The soils used have been Sydney sand, which has the properties d50 = 0.3 mm, coefficient of
uniformity cu = 3, maximum and minimum voids ratios emax = 0.8, emin = 0.58, and non-plastic
feldspar fines which has been mixed with the sand to form a 2:1 sand:silt mixture with properties
d50 = 0.2 mm, cu = 22, emax = 0.89, emin = 0.51.
The results of high speed photography combined with DIC have been used to investigate the
pattern of deformation in the soil. Figure 20(a) shows a typical displacement pattern after an
impact on a dry loose sand specimen. It is noticeable that the impact mobilises a general bearing
capacity mechanism with much of the sand displaced sideways and upwards as the plunger
penetrates the soil. It can also be observed that there is a region beneath the plunger where the soil
is displaced downwards, and this is the region of most interest for dynamic compaction in practice.
If the region beneath the plunger is looked at in more detail the strain pattern shown in Figure 20(b)
is obtained. This figure shows the cumulative shear strains 0.01 sec after the impact and indicates a
sequence of shear (compaction) bands that have been shown (Nazhat and Airey, 2012) to propagate
down from the plunger. Volumetric strains in this region are only significant within these shear
strain bands. After a compaction band has passed any location the shear and volumetric strains
return almost to zero. However, at the end of the process the compaction bands are locked in the
soil. In a series of dynamic impacts, each impact causes additional shear (compaction) strain bands
to propagate and accumulate in the soil so that after several impacts the cumulative strain contours
appear as in Figure 21a. The pattern associated with these strain bands depends on the soil type
and tamper geometry (Nazhat and Airey, 2012). As an example, the effect of 6 impacts on a
compressible sand-silt mixture is shown in Figure 21(b). In this case the general bearing capacity
mechanism is absent and the energy from the impact is used more effectively to compact the soil
beneath the plunger.
Figure 20. (a) Displacement vectors produced by an impact on loose sand, and (b) Cumulative shear
strain contours 0.01 s after impact beneath the plunger
Figure 21. Cumulative shear strains after 6 impacts on (a) loose sand, (b) loose sand-silt mixture
In addition to the photography accelerometers and load cells have been placed on the plunger
and buried in the soil. For comparison with the numerical results for clay a typical result showing
the average stresses mobilised by the plunger during an impact and during slow static loading is
shown in Figure 22. This figure shows that very similar responses are observed, at least during the
early stage of soil penetration. It is expected that these particular results will be more alike than
would be the case for clay, because rate effects are less significant in sand. It has also been
observed that the mechanism for sands is similar in both static and dynamic loading (Nazhat,
2013), the primary difference is the presence of the compaction bands at depth, but these represent
only a small part of the total energy dissipation.
120
100
Average applied stress (kPa)
Dynamic Test
Static Test
80
60
40
20
0
0
10
20
30
40
50
60
70
80
90
100
Displacement (mm)
Figure 22. Soil resistance versus depth of penetration for static and dynamic loading
The existence of shear/compaction bands observed in these experiments presents a significant
challenge for numerical modelling, and is unlikely to be predicted by conventional dynamic
analyses adopting the more common, continuum-based constitutive models used regularly to
represent the static behaviour of soils. Further research is under way to understand these shear
bands and to develop numerical analyses capable of reproducing this behaviour.
9
CONCLUSIONS
Finite element solutions of problems involving impact compaction and dynamic penetration of soil
deposits, by either free falling or propelled objects, have been studied. A robust numerical method
was described for dealing with such complex and difficult numerical problems.
The method was employed in a parametric study of a variety of penetrometers free-falling into
layers of inhomogeneous clay deforming under undrained conditions. The effect of the mechanical
properties of the clay soil on the penetration characteristics was discussed and some of the
predictions were compared with the results of laboratory experiments. It was found that these
smaller penetrometers decelerated at approximately uniform rate soon after initial impact with the
soil.
A penetrometer of much larger physical dimensions, viz., a torpedo anchor typical of those used
in offshore moorings, was also investigated numerically. The key difference in response for these
more massive penetrometers was the prediction of a non-uniform rate of deceleration during
penetration of the soil. As for the smaller penetrometers the depth of final penetration was a
function of the kinetic energy with which the penetrometer first impacts the soil.
How the numerical results for essentially cylindrical penetrometers with a conical tip might be
used in the interpretation of in situ penetration tests was then suggested. It was found that the
known initial impact energy and final penetration depth could be used to estimate the undrained
shear strength of the soil.
Finally, the numerical method was also adopted to examine the fundamental compaction
behaviour of soils and the trends in predicted results were compared, at least broadly, with
observations made in physical experiments. It was noted that further research is required to
understand the deformation mechanisms observed in the compaction experiments and to develop
numerical analyses capable of reproducing such behaviour. It appears that current methods based
on continuum models of soil behaviour are unable to reproduce the behaviour observed in the
physical experiments.
10 REFERENCES
Abelev, A., Simeonov, J. & Valent, P. 2009, Numerical investigation of dynamic free-fall penetrometers in
soft cohesive marine sediments using a finite element approach, Oceans 2009, IEEE/MTS conference,
Biloxi.
Beard, R.M. 1977, Expendable doppler penetrometer, Technical Report, Naval Construction Battalion
Center, Port Hueneme, California.
Brown, M.J. & Hyde, A.F.L. 2008, High penetration rate CPT to determine damping parameters for rapid
load pile testing. Proc. 3rd International Conference of Geotechnical and Geophysical Site
Characterisation, Taiwan. Pub. Balkema. 657-663.
Carter J.P., Nazem M., Airey D.W. and Chow S.W. (2010) Dynamic analysis of free-falling penetrometers in
soil deposits. Geotechnical Special Pub. 199. GeoFlorida 2010 - Advances in Analysis, Modeling and
design. American Society of Civil Engineers, West Palm Beach, Florida, 53-68.
Chow, S.H. 2012. Rate effects in free-falling penetrometer tests in clay. PhD Thesis, The University of
Sydney, Australia.
Chung, J. & Hulbert G.M. 1993, A time integration algorithm for structural dynamics with improved
numerical dissipation: the generalized-a method. Journal of Applied Mechanics, 60, 371–5.
Dayal, U. & Allen, J.H. 1973, Instrumental Impact Cone Penetrometer. Canadian Geotechnical Journal, 10,
3, 397-409.
Einav, I. & Randolph, M. 2006, Effect of strain rate on mobilised strength and thickness of curved shear
bands” Géotechnique, 56, 501–504.
Graham, J., Crooks, J.H.A. & Bell, A.L. 1983, Time effects on the stress–strain behaviour of natural soft
clays. Géotechnique, 33, 327–340.
Hu, Y. & Randolph. M.F. 1998, A practical numerical approach for large deformation problems in soils.
International Journal for Numerical and Analytical Methods in Geomechanics, 22, 327-350.
Kellezi, L. 1998, Dynamic soil-structure interaction. Transmitting boundary for transient analysis. PhD
Thesis, Series R No. 50, Department of Structural Engineering & Materials, DTU, Denmark.
Kontoe, S., Zdrakovic, L. & Potts, D.M. 2008, An assessment of time integration schemes for dynamic
geotechnical problems”, Computers and Geotechnics, DOI: 10.1016/j.compgeo.2007.05.001.
Kontoe, S. 2006, Development of time integration schemes and advanced boundary conditions for dynamic
geotechnical analysis, PhD thesis, Imperial College, University of London.
Lehane, B.M., O'Loughlin, C.D., Gaudin, C. & Randolph, M.F. 2009, Rate effects on penetrometer resistance
in kaolin. Géotechnique, 59, 41–52.
Lysmer J.R. & Kuhlemeyer L. 1969, ‘Finite dynamic model for infinite media. Journal of the Engineering
Mechanics Division, ASCE, 95(EM4), 859-77.
Medeiros, C.J. 2002, Low cost anchor system for flexible risers in deep waters. Proceedings of the 34th
Annual Offshore Technology Conference, Houston, Texas, Paper OTC 14151.
Mulhearn, P.J. 2002, Influence of penetrometer probe tip geometry on bearing strength estimation for mine
burial prediction, DSTO, TR 1285.
Nazem M., Sheng D. & Carter J.P. 2006, Stress integration and mesh refinement in numerical solutions to
large deformations in geomechanics. International Journal for Numerical Methods in Engineering, 65,
1002-1027.
Nazem, M., Sheng, D., Carter, J.P. & Sloan, S.W. 2008, Arbitrary-Lagrangian-Eulerian method for largedeformation consolidation problems in geomechanics. International Journal for Analytical and
Numerical Methods in Geomechanics, 32, 1023-1050.
Nazem, M., Carter, J.P. & Airey, D. 2009a, Arbitrary Lagrangian-Eulerian Method for dynamic analysis of
Geotechnical Problems. Computers and Geotechnics, 36 (4), 549-557.
Nazem, M., Sheng, D., Carter, J.P. & Sloan, S.W. 2009b, Alternative stress-integration schemes for large
deformation problems of solid mechanics. Finite Elements in Analysis and Design, 45, 934-943.
Nazem, M., Carter, J.P., Airey, D.W. & Chow, S.H. 2012, Dynamic analysis of a smooth penetrometer freefalling into uniform clay. Géotechnique, [http://dx.doi.org/10.1680 /geot.10.P.055].
Nazhat, Y. (2013) The deformation mechanism in dynamic compaction, forthcoming PhD thesis, University
of Sydney.
Nazhat, Y. and Airey D.W. (2012) The Effect of Different Tamper Geometries on the Dynamic Compaction
of Sandy Soils, Int. Symp. Ground Improvement, IS-GI, Brussels, May.
Nazhat Y and Airey DW (2011) Applications of high speed photography and X-ray computerised
tomography (Micro CT) in dynamic compaction tests, International Symposium on Deformation
Characteristics of Geomaterials, September 1~3, 2011, Seoul, Korea.
O'Loughlin, C.D., Randolph, M.F. & Richardson, M.D. 2004, Experimental and theoretical studies of deep
penetrating anchors, Proceedings of Offshore Technology Conference, Houston, TX, USA.
Randolph, M. & Gourvence, S. 2011. Offshore Geotechnical Engineering, New York, Spon Press.
Scott, R.F. 1970, In-place ocean soil strength by accelerometer, Proc ASCE, J. Soil Mechanics and
Foundations, 96, 1, 419-444.
Sheng, D., Nazem, M. & Carter, J.P. 2009, Some computational aspects for solving deep penetration
problems in geomechanics. Computational Mechanics, 44, 549-561.
True D.G. 1975, Penetration of Projectiles into Sea Bed Floor Soils, Civil Engineering Laboratory, Technical
Report R-822, Port Hueneme, CA.
Wriggers, P. 2006, Computational Contact Mechanics, Springer-Verlag, Austria.
11 ACKNOWLEDGEMENTS
The authors would like to acknowledge the support of the Australian Research Council in funding
some of the work described in this paper.

Download Report