Pile Soil Interaction

ST A TIC PILE-SOIL-PILE INTERACTION IN
OFFSHORE PILE GROUPS

Brasenose College
by
Anjum Rashid Chaudhry
A thesis submitted for the Degree
of Doctor of Philosophy at the
University of Oxford
Michaelmas Term 1994
ABSTRACT
STATIC PILE-SOIL-PILE INTERACTION IN OFFSHORE PILE
GROUPS
A. R. Chaudhry
Brasenose College, University of Oxford
A thesis submitted for the Degree of Doctor of Philosophy
Michaelmas Term 1994
This thesis is a theoretical study, using both finite element and boundary element methods,
of the behaviour of single-piles and pile groups under vertical and lateral loading. It offers
an improved understanding of the soil-structure interaction that occurs in pile groups,
particularly closely spaced piles subjected to lateral loads. The potential of a two-dimensional
idealisation of what is a three-dimensional problem is demonstrated by achieving real insight
into the complex nature of pile-soil and pile-soil-pile interaction in pile groups.
A new load transfer mechanism is presented for a rigid, axially loaded vertical pile. From this
an improvement is then derived to the analytical solution for pile head settlement given by
Randolph and Wroth (1978). The improved mechanism has the further merit that it can be
applied also to solutions for flexible piles and pile groups.
The improved analytical solution is further adapted in the development of two correcting
layers specifically for venically loaded piles to model infinite boundaries in the finite element
model. The correcting layers help in establishing superiority of the finite element method over
the boundary element method.
To model pile-soil interaction, a purely cohesive interface element is developed and then
validated by performing various two-dimensional test problems, including stability analysis
of flat surface footings. Footing-soil interface tension is successfully modelled in this way -
an outcome that entails a significant modification to the Hansen (1970) bearing capacity
solution. Stability analysis is also carried out of conical footings using a three-dimensional
finite element model: the results help to explain the applicability of the existing bearing
capacity theories to conical footings.
The ultimate lateral soil reaction is determined and various pile loading stages are investigated
through parametric studies. Study of the stage immediately following pile installation (i.e. the
ii
consolidation stage) highlights the need to develop an effective stress analysis for laterally
loaded piles. Pile-soil interaction is studied using the cohesive interface element presented
earlier, which proves to be quite successful in smoothing out the stress discontinuities around
the pile.
A new material model for frictional soils is presented. and validated by using it to model an
extension test: it captures well post-peak behaviour and takes care of the effects of dilation
on the response of laterally loaded piles.
Finally, mechanisms of interaction in closely spaced pile groups are studied. Simple analytical
expressions are derived which quantify the effects of interaction. A new method of analysis
is presented for single-piles and pile groups which offers a considerable degree of reliability
without having to do either impossibly expensive full scale field tests or prohibitively
expensive full three-dimensional analysis using the currently available computers.
References
Hansen, J.B. (1970)
"A revised and extended formula for bearing capacity",
Danish Geotechnical Institute, Copenhagen, Bull. No. 28: pp 38-46
Randolph, M.F. and Wroth, C.P. (1978)
"Analysis of deformation of vertically loaded piles",
Proc. ASCE, Journal of the Geotechnical Engineering Division, 104 (GT12): pp 1465-
1488
iii

Many thanks and praise are due first of all and above all to my Creator who made this thesis
possible despite all the hurdles and constraints in the way of its completion: 'If Allah helps
you, none can overcome you' (Qura'n, 3.160).
I acknowledge with gratitude the support I have received as an ICI-scholar: my thanks to ICI
(UK), the Foreign and Commonwealth Office and Oxford University, who jointly
administered this award. I am grateful also to the High Commission for Pakistan in London
for making a grant from the Dr Wali Mohammad Trust; to the Mott MacDonald Group
(Consulting Engineers) for their travel grant; to the Sir Ernest Trust (Mountbatten Memorial
Grant), the Sir Richard Stapley Educational Trust; and to the Professional Classes Aid Council
for their very valuable (and much needed) financial assistance during the final year of my
studies.
Sincere thanks go to the late Professor C.P. Wroth for his encouragement and advice in the
early stages of this research: the privilege of studying under him will always remain a
memorable part of my education.
I am grateful to Professor G.T. Houlsby for his assistance and inspiration to complete this
work. I acknowledge especially his allowing me to use the following: OXMESH - finite
element mesh generation and 2CAN - data contouring programs.
Special thanks to my supervisors, Dr H.J. Burd and Dr M.S. Williams, for their constant help,
and patience throughout my study years at Oxford, especially while I got used to a new
country and a new language. I am indebted to Dr H.l. Burd for allowing me to use his finite
element program OXFEM and arranging for me to use the commercial finite element package
ABAQUS. I would also like to thank him for his financial help through paid work for
"Edmund Hambly Limited".
My moral tutor Dr R.W. Daniel and Dr G.C. Sills helped me with consistently positive,
constructive criticism of my work: my warmest thanks.
Because it is not possible to mention everyone individually, may I express a 'global'
appreciation to all members of the Civil Engineering Group for providing a lively
environment to work in, and for moral support when needed.
Last but not least, I should like to thank Dr Z.J. Qureshi, Dr AL. Bello and Mr M. Appabai:
none of them are engineers; their help was as friends who made themselves available
whenever needed. I cannot express my appreciation better than in the Arabic of
jaziikum Alliihu khairan, being certain they will understand perfectly.
iv
CONTENTS
Abstract
Acknowledgements
Contents
Nomenclature
CHAPTER 1 Introduction
1.1
1.2
1.3
1.4
Introduction
Description of the problem
Objective of this study
Structure of the study
CHAPTER 2 Analytical and Numerical Methods:
2.1
2.2
2.3
2.4
2.5
2.6
A review of the relevant literature
Analytical methods
Numerical methods (linear-elastic soil behaviour)
Numerical methods (for non-linear elastic soil behaviour)
Numerical methods (for elasto-plastic soil behaviour)
Plane strain solutions for laterally loaded piles
Summary
CHAPTER 3 Elastic Analysis of a Single-Pile under
Axial Loading
3.1
3.2
3.3
3.4
Introduction
Axi-symmetric quasi constant
boundary element formulation
Two-dimensional elasto-plastic
boundary element method
Analytical methods for the analysis of
rigid single-piles
V
11
IV
V
IX
1-1
1-1
1-2
1-4
1-5
2-1
2-1
2-7
2-12
2-17
2-18
2-21
3-1
3-1
3-3
3-10
3-13
3.4.1
3.4.1.1
3.4.1.2
3.4.2
3.4.3
3.5
3.5.1
3.5.2
Simplified analysis of a rigid single-pile
(after Randolph and Wroth (1978
Effect of lower layer of soil on upper layer
Modelling the base of the pile
Simplified analysis of a rigid single-pile:
(proposed new method)
Comparison of results
Finite element analysis
using a correcting layer approach
Development of the correcting layer
around the pile shaft
Correcting layer under the pile base
3.5.3 Finite element analysis results
3.6 Conclusions
Appendix 3A Determination of influence coefficients
CHAPTER 4 Interface Analysis
4.1
4.2
4.3
4.3.1
4.3.2
4.3.3
4.3.4
4.4
4.4.1
4.4.2
Introduction
Background
Six-noded isoparametric interface element formulation
Kinematic equations
Constitutive equations
Finite element equations
Stress update calculations
Verification of the plane strain interface element
Two jointed solid blocks
Plane strain footing subjected to
moment and vertical loading
4.4.2.1 Analytical solution (loss of footing contact)
4.4.2.2 Finite element analysis
Appendix 4A
CHAPTER 5 Analysis of Single Piles under Lateral Loading
5.1
5.2
5.3
5.4
5.5
5.6
5.7
Introduction
Choice of material constitutive behaviour
Choice of the finite element
Finite element analysis of a single-pile
Effect of mesh density on solution quality
Effect of varying soil parameters
Analysis of a single-pile assuming initial
cavity expansion to model installation
vi
3-14
3-17
3-18
3-20
3-23
3-25
3-26
3-28
3-30
3-33
3A-l
4-1
4-1
4-2
4-4
4-5
4-6
4-10
4-12
4-17
4-17
4-19
4-20
4-23
4A-1
5-1
5-1
5-3
5-6
5-9
5-19
5-21
5-22
5.8
5.9
5.10
Analysis of a single-pile assuming radial
variation of soil strength around the pile
Interface analysis
Effect of no-tension behaviour of pile-soil interface
CHAPTER 6 Plane Strain Analysis of Pile Groups
6.1 Introduction
6.2 Plane strain analysis of 2-pile groups
6.3 Plane strain analysis of a closely spaced 2x2 pile group
6.4 Determination of pile-soil-pile interaction factors
6.4.1 Approximate prediction of load-displacement curve
for circular pile groups
6.5 Plane strain analysis of multiple pile groups
6.5.1 Plane strain analysis of 5-pile group
6.5.2 Plane strain analysis of lO-pile group
6.5.3 Back analysis of field tests
Appendix 6A
CHAPTER 7 Analysis of Laterally Loaded Piles in
Frictional Soils
7.1
7.2
7.3
7.3.1
In troduc tion
Background
Material stiffness matrix
The Matsuoka yield surface - constant friction
and dilation angle (model 1)
7.3.1.1 Finite element analysis of a single-pile
7.3.2 The Matsuoka model - variable friction
and dilation angle (model 2)
7.4
7.4.1
7.5
7.5.1
7.5.2
7.6
7.7
Stress update calculations
Calculation of yield surface intersection
Verification of the Matsuoka model (2)
Plane strain analysis of a unit square soil element
Laterally loaded single-pile
(after Lane and Griffiths (1988
Laterally loaded single-pile - the effect of
initial relative density J DO
Conclusions
CHAPTER 8 Bearing Capacity of a Jackup Spudcan Footing
8.1 Introduction
VB
5-26
5-31
5-32
6-1
6-1
6-4
6-11
6-13
6-17
6-19
6-20
6-23
6-23
6A-l
7-1
7-1
7-2
7-6
7-7
7-10
7-13
7-14
7-15
7-17
7-17
7-19
7-20
7-21
8-1
8-1
8.2 Background
8.3 Description of the Problem
8.4 Choice of material constitutive behaviour
8.5 Finite element Model
8.6 Finite element Analyses
8.6.1 Vertical bearing capacity of a spudcan footing
(vertical preload to cause yield)
8.6.1.1 Vertical bearing capacity solution
8.6.1.2 Finite element analyses results (load path 1)
8.6.2 Reduction of vertical load (vertical unloading)
with increase of horizontal load to cause yield
8.6.2.1 Bearing capacity solutions
8.6.2.2 Finite element analysis result (load path 2)
8.6.2.3 A further reduction of vertical load followed by
further increase in horizontal and vertical loads
(load path 3)
8.7 Comparison of results
8.8 Conclusions
CHAPTER 9 Conclusions
9.1
9.2
9.3
9.4
9.5
Concluding remarks
Axial load capacity and axial response of piles
Interface analysis
Lateral load carrying capacity and response of
single piles and pile groups
Future research
BIBLIOGRAPHY
viii
8-5
8-7
8-8
8-11
8-12
8-13
8-13
8-16
8-20
8-20
8-22
8-27
8-27
8-31
9-1
9-1
9-3
9-3
9-4
9-4
NOMENCLATURE
Additional symbols and occasional alternative usage of the symbols below are defined where
they occur in the text.
A Cross-sectional area of the pile
Ac Circumferential area of a pile element
A Plan contact area of a partially penetrated spudcan
p
[B] Matrix relating nodal velocities
B Footing width
B' Effective footing width
B Penetrated diameter of a spudcan
p
C
L
Centre line
c Cohesion at the interface
d Pile diameter
D p Penetration depth of a spudcan
[D] Material stiffness matrix
[D]E Elastic material stiffness matrix
[DY Plastic material stiffness matrix
[D]EP Elasto-plastic material stiffness matrix
EI Bending stiffness
EjEs Flexibility ratio
ix
e Void ratio
Soil flexibility matrix
/<cr)
Yield function
G Shear modulus
H Horizontal/lateral load
H Peak lateral resistance
p
g(cr) Plastic potential
Relative density
Stress invariants
1
r
Rigidity index
Relati ve dilatancy index
[J] Jacobian matrix
K Bulk modulus
[K] Stiffness matrix
K
p
Coefficient of passive earth pressure
k Coefficient of subgrade reaction
Coefficient earth pressure at rest
Interface element shear stiffness
kll
Interface element normal stiffness
I Length of the pile
Critical length
I
e
Length of a pile element
l.
I
Length of an interface element
x
M
[N]
N
c
Pull
p

P
P
p ..
q
R
s ..
s
t
u
u
V,.
w
Bending moment
Shape function matrix
Bearing capacity factor due to cohesion
Bearing capacity factor due to self-weight
Ultimate bearing capacity of a plane strain footing
Axial load at depth z
Vector of applied loads
Volumetric stress
Ultimate lateral soil resistance
Deviatoric stress
External radius
Magical radius
Pile radius
Undrained shear strength of clay
Centre-lo-centre pile spacing
Load transfer
Pile head deflection
Vector of nodal velocities
Relative tangential velocity
Relative normal velocity
Vertical load
Ultimate vertical load capacity of the spudcan
Pile head settlement
xi
x
y
z
a
E
y
<l> cv
v
Vector of global nodal co-ordinates
Lateral deflection
Depth below ground surface
Isoparametric reference co-ordinate
Interaction Factor between piles kJ
Angle of applied loading
Principal strains
Volumetric strain
Soil effective unit weight
Specific volume at a reference pressure P a
Slope of the Critical State Line
Displacement of a single isolated pile under unit horizontal load
Lateral displacement of pile k
Depth factor
Load spread angle through the soil under the pile base
Critical state friction angle
Mohr-Coulomb friction angle
The angle of peak friction
Triaxial compression friction angle
Dilation angle
The angle of maximum dilation
Triaxial compression dilation angle
Poisson's ratio
xu
cr, 't
Cone angle
Normal stress, Shear stress
Principal stresses
Xlll
Chapter 1 Introduction
CHAPTER 1
INTRODUCTION
1.1 Introduction
For just over three decades now, with rapidly advancing technology and ever-increasing
energy demands, oil and gas exploration and extraction activity have been moving from
shallow water to deep offshore waters. Most of the offshore oil extraction platforms are pile-
supported steel structures known as 'template' or 'jacket' structures. These consist of a space
frame with legs through which the piles are driven into the sea-bed. Once driven to the
required depth, the piles are cut off at the head of the template and a prefabricated deck
section is placed on the piles at an elevation above the crests of anticipated storm waves. For
deep water applications, the template or jacket structure is modified. In the modified form,
called a 'tower structure', the deck res ts on four large-diameter columns, each in turn
supported by several closely spaced piles in a circular group. Such pile groups or deep
foundations are subject to severe environmental loads (wind, wave action and current loading)
in addition to gravity loads from the weight of the structure acting vertically.
Other offshore structures include concrete gravity platforms and mobile jackup rigs. Gravity
structures have flat circular bases directly resting on the seabed. These structures achieve
stability through their own weight Mobile jackup rigs have either flat circular footings or
conical shaped with a protruding tip at the centre (commonly known as spudcan footings).
These footings are referred to as shallow foundations.
Environmental loads act laterally and produce horizontal and moment loads. A significant
1-1
proportion of these loads is cyclic in nature (see Figure 1.1). The pile-soil interaction in deep
foundations and footing-soil interaction in shallow foundations, arising due to cyclic loads is
an intricate phenomenon and is beyond the scope of this thesis. We will be dealing with static
loads only.
Waves
---------
Earthquake
~ ~
Figure 1.1: Environmental loads on a template platform (after SeInes (1982
1.2 Description of the problem
Generally, onshore pile groups are subjected to vertical loads but often there are situations,
when a significant lateral load is also experienced by these pile groups e.g. bridge abutment
piles located adjacent to surcharge loads on alluvial soils (Springman (1984),(1989)). The
magnitude of vertical and lateral loads carried by onshore pile groups is substantially lower
compared to that supported by offshore pile groups. It is essential therefore to calculate the
lateral response of piles, as well as their axial response, in order to design for the stability and
safety of the overall structure. Similarly, the behaviour of shallow foundations under storm
1-2
loading conditions is also critical for the safety of a platfonn/rig. In this thesis, the axial and
lateral behaviour of closely spaced pile groups, and the behaviour of flat and spudcan footings
under combined (horizontal, vertical and moment loading) is discussed both in soft clays and
in sands. An undrained material behaviour or the 'total stress approach' is adopted for clays,
and the fully drained material behaviour or the 'effective stress approach' is considered for
sands.
The behaviour of a pile in a group is affected by the spacing between the neighbouring piles.
If the piles are wide apart (in the case of onshore pile groups), the displacement of a pile may
have an insignificant influence of the movement of another pile in its neighbourhood. On the
contrary, if the piles are closely spaced (as in the case of offshore pile groups), neighbouring
piles interact with each other through the soil medium. In this situation, interactive effects will
cause an increase in displacement of a pile in such a group. Before studying the behaviour
of a pile affected by interaction, we shall study what happens to the soil when a single
isolated pile is loaded vertically or displaced laterally.
The resistance to deflection which a laterally loaded single-pile develops is a non-linear
function of the relative displacement between it and the soil. There is an increase in lateral
stresses in the soil near the pile which decay rapidly in magnitude with radial distance. There
is an annular region or a zone of high stresses close to the pile in which the soil gets yielded.
For closely spaced piles within groups, these yielded zones of the soil around individual piles
overlap, fonning larger yielded zones in the soil surrounding the pile group (see Figure 1.2).
The interaction arising due to overlapping of these yielded zones results in a larger deflection
for the group before a level of lateral resistance is reached equivalent to that for a single-pile.
In order to develop a method for the analysis of pile groups, it is first necessary to obtain
estimates of the reduction in lateral resistance of the soil which result from the overlapping
of the yielded zones in the soil around the group. Once a reasonable estimate is obtained, then
1-3
Chapter I Introduction
lones of plastic now
a) wide spacing
b) close spacing
Figure 1.2: Schematic illustration of pile group under lateral loading (after Bogard
and Matlock (1983
the analysis of pile groups can be performed using the conventional methods of analysis
developed for single-piles, provided the structural restraints are identical for each pile in the
group, i.e. all piles in the group respond equally to the applied load.
1.3 Objective of this study
The phenomenon of overlapping yielded zones (see Figure 1.2 (b)) in offshore pile groups is
very complex. It is difficult to model it accurately unless we have a suitable soil model
capable of representing real non-linear soil behaviour. The objective of this study is to explore
the most appropriate method of analysis equipped with a suitable, simple and realistic soil
model, i.e. one able to predict the response of a pile group in a way that mimics actual
response in the field. The larger goal is to develop a method of predicting lateral response
which is reliable enough to be applied in the design of routine pile foundations. In the case
of shallow foundations, the objective is to examine how the soil-structure interaction affects
1-4
the stability of these foundations.
1.4 Structure of the study
Previous research related to this problem has been based on various analysis procedures and
techniques which, for one reason or another, are unsatisfactory. A critical review of the
methods currently in practice with reference to the relevant literature is presented in Chapter
2.
Chapter 3 discusses analytical methods and compares two numerical methods, namely the
integral equation method and the finite element method for investigating the response of a
single vertical pile under axial loading. The problem is deliberately chosen to show how its
axial symmetry can be exploited by simplifying it to two dimensions. The discussion also
shows that, in the design of an offshore pile foundation, the precise determination of axial
capacity and response is as important as the lateral load capacity and lateral response, and that
both need to be considered together.
The comparison of the numerical methods clearly demonstrates the superiority of the finite
element method over the integral equation/ quasi boundary element method for problems that
can be idealised into two-dimensional problems and when elaslO-plastic soil behaviour is an
important consideration. The preferability of the finite element method is reinforced by
presenting a computationally economical approach for carrying out the finite element
calculations. The finite element method is therefore adopted in all the subsequent chapters.
Chapters 4, 5, 6 and 7 are devoted to study of the lateral loading of single-piles and pile
groups.
1-5
Chapter 4 is concerned with the modelling of incompatibility in the displacements which
occur at the pile-soil interface where two different materials meet at the common boundary.
For this purpose, a purely cohesive interface element is developed. Its ability to model the
interface behaviour is demonstrated by using it in two two-dimensional test problems. One
of these problems is concerned with the calculation of the ultimate bearing capacity and
stability analysis of a flat surface footing.
It is important to understand the problem of a single-pile under lateral loading before moving
on to the more complex problem of a pile group. Therefore, Chapter 5 provides this necessary
preliminary investigation before the study of pile group behaviour in Chapter 6. It
demonstrates that a two-dimensional idealisation for laterally loaded piles is quite appropriate
to study their behaviour. Each of the three soil loading stages, namely pile installation,
subsequent consolidation of soil around the pile, and the application of foundation loads, is
studied. The lateral capacity of a single-pile is found to be enhanced by almost 35% due to
consolidation of soil around the pile. It seems appropriate to consider the effect of
consolidation in the analysis of laterally loaded piles as previously suggested by Burland
(1973) for axially loaded vertical piles. Finally, the interface element is shown to be quite
successful in smoothing out the stress discontinuities around the pile and to be quite effective,
therefore, in modelling the pile-soil interaction.
Chapter 6 considers the behaviour of pile groups subjected to lateral loading and, in the
absence of a complete three-dimensional finite element analysis, presents a new approach for
the analysis of laterally loaded pile groups. It combines a two-dimensional finite element
1-6
analysis with a standard p-y curve analysis for predicting the pile head response.
Following the study of pile-soil interaction in cohesive soils presented in Chapter 5, a new
material model to represent the constitutive behaviour of frictional soils is presented in
Chapter 7. The validity of the model is tested by predicting the results of an extension test.
It captures the post-peak behaviour successfully and models variation of dilation as the
analysis proceeds. In this way, effects of dilation on the ultimate response of laterally loaded
piles are duly taken into account.
Following the brief discussion, in Chapter 4, of the ultimate bearing capacity of a flat footing,
Chapter 8 considers the ultimatt bearing capacity of the conical or spudcan footing now
widely in use in offshore situations. The importance of being able to predict the load such
footings can support, and the need (in the context of this thesis) to see how three-dimensional
finite element analysis can be applied to carry out the required calculations, are both self-
evident.
Chapter 9 presents an overview and an assessment of the arguments and discussions in the
preceding chapters. Conclusions are drawn and immediate relevant further research is
proposed.
1-7
Chapter 2 Literature Review
CHAPTER 2
ANALYTICAL AND NUMERICAL METHODS:
A REVIEW OF THE RELEVANT LITERATURE
In this chapter, a review of the existing analytical and numerical methods for the prediction
of both axial and lateral response of single-piles and pile groups is presented. The discussion
first covers analytical methods (Section 2.1), then numerical methods for linear elastic models
of soil behaviour (Section 2.2), and, finally, numerical methods for non-linear elastic (Section
2.3) and elasto-plastic constitutive soil behaviour (Sections 2.4 and 2.5).
2.1 Analytical methods
Vertical loads
The problem of vertical displacement (settlement) of a single-pile under axial loading is
statically indeterminate: no exact analytical solution is available as the actual stress
distribution along the pile is unknown. However, Randolph and Wroth (1978) produced an
approximate analytical solution for a rigid pile by considering separately the load sustained
by the pile shaft and the load sustained by the pile base. They assumed a uniform shear stress
distribution along the pile shaft, and the normal stresses beneath the pile base as those of
arising from a rigid punch at the surface of an elastic half space. In the case of a rigid pile,
the pile head load-settlement ratio is determined by a simple analytical expression without it
being necessary to evaluate the exact distribution of shear stress along the pile. The analytical
2-1
Chapter 2 Literature Rel'jew
expression involves an empirically determined parameter' r m' - the radial distance at which
the vertical displacements are zero. (This parameter will be discussed in detail in Chapter 3
and a new analysis proposed which is thought to be an improvement on the Randolph and
Wroth solution.)
In the same paper, Randolph and Wroth also produced modified solutions to take account of
pile compressibility and radial and vertical inhomogeneity of the soil. In these solutions, a
careful choice of the values of soil shear moduli is very important to get sensible results. In
another paper, Randolph and Wroth (1979a) applied their analytical approach in the analysis
of two-pile groups and presented simple expressions for the elastic interaction between the
piles. The interaction is quantified as interaction factors. These interaction factors, which
depend on the spacing between the piles, were then used in the analysis of groups of piles
symmetrically placed around a pitch circle, where each pile responded identically. They also
suggested modifications in the analysis of pile groups when piles are not subjected to the
same load, for example under a rigid pile cap. (The computer program PIGLET developed
by Randolph (1980), one frequently used by practising engineers, is based on the approximate
analytical solution for single-piles and the simple expressions for interaction between piles.
The soil is modelled as an elastic material with a linear variation of stiffness with depth. Piles
may be raked in any direction, and effect of not only axial but also lateral and torsional
loading can be analysed.) The solutions given by Randolph and Wroth (1978, 1979a) are
excellent for idealised hnear-elastic soil behaviour (i.e. it does not consider soil non-linearity,
soil yielding and soil plasticity); it works very well for most practical situations where piles
are installed at spacings wide enough for the interactions among the piles to remain in an
elastic regime and the principle of superposition to remain valid. However, in the case of
closely spaced piles, especially in an offshore environment, this principle is not applicable.
Lateral loads
In the case of lateral loading of piles, the problem is not as simple as for vertical loading in
2-2
that no simple analysis has yet proved possible. The horizontal pressure distribution along the
pile is unknown as the behaviour of a pile is flexible under lateral loading. (In the case of a
vertically loaded pile. it is the shear stress distribution which is unknown.) However, for
purely cohesive materials, a simplified distribution of soil pressure was suggested by Broms
(1964a) for short, intermediate and long piles. According to Brams, the soil pressure is zero
from ground surface to a depth of 1.5d and, below this depth, a constant 9S" regardless of
the shape of the pile and the nature of the pile-soil interface. Where S" is the undrained shear
strength of the soil. The value of 9S" was largely empirical and no theoretical justification
was attempted. Broms (1964a) also pointed out that the ultimate lateral soil pressure which
develops along a laterally loaded pile in a cohesive soil can be calculated by means of the
theory of plasticity below a depth of approximately 3d, assuming that plane strain conditions
hold, that the soil is weightless, and that failure of soil takes place along plane or spiral slip
surfaces.
Reese (1975) presented two simple expressions for the determination of the ultimate lateral
soil pressure Pus C at any depth z below the surface:
(
'Vz 0.5ZJ
p c=pjd= 3 +_1 +_ S
us S d "
u
for shallow depths
and
for greater depths
where
P" :: the ultimate lateral soil resistance,
d = diameter of the pile and y :: soil unit weight (bulk).
The expressions given above seems to suggest that the ultimate lateral soil resistance, p u' in
a purely cohesive soil increases from a value of 3S"d at the ground surface to 9S
u
d down to
2-3
a depth of about 3d and remains constant for greater depths.
Later, Randolph and Houlsby (1984) showed that the ultimate lateral soil resistance is very
much dependent upon the pile-soil interface. Using the method of characteristics, they
produced exact analytical solutions for the ultimate lateral soil resistance which develops in
front of a laterally loaded rigid circular disc - 11.94S"d for a completely rough pile-soil
interface, and 9.145 d for a smooth interface.
"
Similarly, for sands, Broms (1964b) suggested a triangular distribution of passive pressure
along the front of a laterally loaded pile ignoring the active pressure at the back of the pile.
The passive pressure being equal to 3 times the Rankine passive pressure:
where
y' = effective soil unit weight
<1>' = angle of internal friction (effective stress)
For the more general case, Hansen (1961 a) suggested the following simplified expression for
the determination of the ultimate lateral soil pressure p us
c
-4> at any depth z below the surface:
where
pC-. :: yzK + 5 K
us q" c
z
Kq = factors that are functions of <p and
d
!jJ = friction angle
For a purely cohesive soil, the factor K is zero and the overburden pressure has no effect
q
2-4
on P
il

For any assumed pressure distribution, such as a uniform or a linear variation, the
combination of the ultimate horizontal load and moment to cause failure of a relatively
slender short floating pile with negligible base resistance, can be calculated using a simple
statical approach given by Poulos and Davis (1980). For a general distribution of soil
pressure, a trial and error procedure is recommended by Hansen (1961 a). The analysis
procedure was extended by Poulos (1976) for the case of a pile near a slope. This type of
analysis neither provides any information about the ultimate lateral pile head displacements
at failure nor can be applied to pile groups.
A pile can also be analysed as a transversely loaded beam carrying some axial load with
different end restraints. The surrounding soil is treated as a series of elastic springs attached
to the pile at a number of discrete points which divide the pile into a series of inter-connected
beam-column elements. This soil representation is after Winkler (1867) (see Scott (1981 and
the analogy derives from the analysis of a beam on an elastic foundation carried out by
Hetenyi (1946), who produced closed form solutions for various loading conditions and end
restraints. The stiffness of the soil springs is referred to as the coefficient of sub grade reaction
k of the soil. It is computed by assuming it to be equal to that of a strip footing founded on
the surface of a semi-infinite, ideal elastic medium. Since the pile is surrounded on all sides
by the elastic medium, this assumption gives a lower value for the coefficient of subgrade
reaction than the actual one. It may be kept constant or variable along the length of the pile.
The lateral deflections calculated using this analysis procedure therefore yield only a
conservative estimate, as was pointed out by Broms (1964a). Further, the soil springs are
attached only to the front of the pile, neglecting completely soil behaviour at the back of the
pile. In any case, modelling the soil as a spring cannot truly represent a soil continuum -
indeed, it completely ignores soil continuity. As Schofield and Wroth (1968) pointed out, a
three-dimensional soil continuum cannot be represented by three-orthogonal coil springs. The
basic Winkler soil model was refined by Vlasov and Leontiev (1966) who added a tension
2-5
Chapter 2
Literature Review
layer at the surface of the spring series, causing the springs to interact with each other and
behave like an elastic continuum. In the same paper, Vlasov and Leontiev also gave consistent
closed form solutions, but these solutions cannot be extended to pile groups.
For long flexible piles, Broms (l964a) and, later, Randolph (1981), showed that the sub grade
reaction analysis of a long flexible pile of bending rigidity EI embedded in soil with
coefficient of subgrade reaction k, yields a critical length of pile beyond which the pile
behaves as if infmitely long and then the length of the pile does not influence the deformation
1
at the loaded end. This critical length Ie may be expressed as Ie ... 4[EI/k]4. It follows from
this, that the whole length of the pile does not need to be considered in the subgrade reaction
analysis. For piles that are longer than their critical length, the solutions of Hetenyi (1946),
for a pile loaded by a lateral load H and bending moment M, yield expressions for the
deflection u and rotation 6 at the loaded end given by:
In some cases analytical methods are not very attractive - not only because there is no
information on lateral deflections (Hansen (1961a, but also because of their cumbersome
solution procedures (Vlasov and Leontiev (1966. Moreover, these methods cannot be
extended to pile groups. In this thesis, therefore, preference is given to numerical methods.
However, the numerical solutions are validated against the analytical solutions, wherever
possible.
2-6
2.2 Numerical methods (linear-elastic soil behaviour)
Vertical loads
In the absence of adequate analytical solutions, various theoretical/numerical methods have
been developed over the years and are currently in use to predict the axial response of a
single-pile. These methods can be classified into two groups:
a) Methods which model the soil as a series of non-linear springs.
b) Methods that treat the soil as a continuum.
The methods which belong to group (a) are presented in Section 2.3. However, an
approximate analysis based on Winkler (1867) soil model is presented by Scott (1981). It is
assumed that a series of elastic soil springs are connected to the pile shaft along the length
of the pile and a spring below its tip. The stiffness of these springs is carefully chosen such
that the results obtained from the approximate analysis match with the elastic continuum
methods given by Poulos and Davis (1968) for rigid piles and Mattes and Poulos (1969) for
flexible piles.
The elastic continuum methods of group (b) are: the integral equation method, presented
initially by Poulos and Davis (1968) and later by Randolph (1977); the boundary element
method proposed by Butterfield and Banerjee (1971); more recently, the infinite layer method
of Guo et al. (1987), in which a rigorous approach is presented to model a layered soil
medium for the analysis of a single-pile and Cheung et al. (1988) for pile groups; and
finally, the finite element method adopted by, among others, Ottaviani (1975), Randolph
(1977) and Baguelin and Frank (1980).
All of these methods may be used to obtain solutions for linear elastic soil behaviour.
However, all these methods except the finite element method (see below), suffer from the
2-7
Chapter 2 Literature R e ~ ' i e w
disadvantage that they cannot be readily extended to more complicated constitutive models.
Although some recent refinements of these methods take account of soil non-linearity and soil
plasticity, they have serious limitations when it comes to dealing with soil inhomogeneity
(vertical and horizontal), and the further drawback of slowness of the numerical calculations.
(Elastic continuum methods are explored in Chapter 3 and a new approach is there presented
to make the finite element analyses more economical.)
It is difficult in the finite element method to model an infinite domain. Usually the boundaries
in a finite element mesh are set at sufficiently large distances for the boundary effects to
become negligible, at least beyond the theoretical extent of the collapse mechanism, as studied
by Bell (1991). However, Brocklehurst (1993), following this approach in his study of a two-
layer soil system, observed some sizeable errors in the finite element results due to the close
proximity of the mesh boundaries, and had to shift them further away to obtain acceptable
results. This difficulty of fixing boundaries in a finite element model had been dealt with in
the past either by the use of special elements (infinite elements) described by Bettess (1977),
Medina (1981) and more recently by Chen and Poulos (1993), or by combining the finite
element formulation with a boundary element formulation (Zienkiewicz et al. (1977). (Chapter
3 takes a different approach, in which a correcting layer is used to model the infinite
boundaries. This approach is similar to the one described by Burd and Houlsby (1990) for a
large strain cavity expansion problem.)
Lateral loads
Various attempts have been made in the past to calculate the lateral deflection of piles under
working loads. The methods used to predict the lateral response can also be classified into
two groups, as described above. One of the earliest analysis procedures belongs to group type
(a) and is known as the subgrade reaction method (see above, Section 2.1). The analytical
solutions given by Hetenyi (1946) can also be obtained using numerical methods. One such
2-8
method, proposed by Poulos and Davis (1980), employs finite difference numerical technique
to solve the beam column differential equation. It is easy to incorporate a spring coefficient
which varies linearly with depth. Again, the method is not applicable to analysis of groups
of piles.
Poulos (1971 a) also proposed an elastic solution, known as the Integral Equation Method
(IEM), based on Mindlin's solution for a point load in an homogeneous, isotropic elastic half
space. This method belongs to group type (b) and employs a semi-analytical approach in
which the pile is assumed to be a thin rectangular strip of the same rigidity as the actual pile.
The pile response is calculated by integrating Mindlin's equation over the corresponding
rectangular area of the soil. The method can be extended to study group behaviour (Poulos
1971 b but the amount of computational work involved is substantially reduced by the
introduction of interaction factors. These factors take into account the additional movement
of an individual pile attributable to the movements of neighbouring piles in its group.
However, subsequent extension of the method to deal with a layered soil profile introduces
a degree of approximation into the solution because of the choice of an equivalent shear
modulus. The accuracy of the solution is further affected by the manner in which the pile is
discretised. Evangelista and Viggiani (1976) reported that the use of varying size elements
along the pile gives improved accuracy over the solutions given by Poulos (1971 a). For piles
of large flexibility, the solutions are strongly dependent upon the length of the elements near
the top of the pile. (A similar trend in the use of the integral equation method for the
settlement analysis of a vertical pile, is reported in Chapter 3 below.) However, the cylindrical
pile idealisation proposed by Randolph (1977), rather than the rectangular strip idealisation
by Poulos (1971a) showed close conformity with finite element analysis results.
Randolph (1981) suggested algebraic expressions which allow the behaviour of the flexible
piles under lateral loading to be calculated in terms of fundamental soil properties. The
expressions are based on the results of finite element studies for the response of a laterally
loaded cylindrical pile embedded in elastic soil with stiffness varying linearly with depth.
Randolph performed a special type of finite element analyses in which he applied non-
2-9
symmetrical loading to the axi-symmetric geometry of the circular pile. This type of analysis
is computationally economical compared to a three-dimensional finite element analysis: it
remains to be further explored using more refined soil models, and for the analysis of pile
groups. Earlier, in a similar study, Kuhlemeyer (1979) presented his finite element results
fitted with simple power law relationships. Randolph (1981), in addition to the algebraic
expressions, also suggested an approximate analysis procedure (without the aid of a computer)
for a flexible pile under lateral loading embedded in a soil with a linearly varying soil
modulus with depth, and having zero value at ground level. Due to the very low soil stiffness
around the upper part of the pile shaft, he suggested that limiting soil pressure will be
reached at large pile deformations. Therefore, the deformation of the pile in this region may
be calculated by treating the pile as a simple cantilever. For the lower part of the pile shaft,
i.e. below this failed region, a single pair of algebraic expressions gives the deformations at
the transition point between the failed and the un-failed regions. Both free-headed or fixed-
headed piles can be analysed in this way. Later, Randolph and Houlsby (1984) presented a
more rational approach: they suggested expressions for the limiting soil reaction, for the two
modes of failure i.e., wedge type of failure near the surface, and the failure at depth which
takes place due to the flow of material round the pile. The analytical expressions for the latter
type of failure give limiting soil reaction per unit length of pile, assuming plane strain
conditions. But it is difficult to obtain such expressions for pile groups. Also, the lateral
deflections still remain unknown. (A numerical method, the finite element method, is used
in Chapters 5, 6 and 7 to determine lateral deflection of single-piles and pile groups.)
An extension of the integral equation method mentioned above is the boundary element
method presented by Butterfield and Banerjee (1971), in which soil continuum and pile are
treated as two separate domains whose boundaries are discretised into a finite number of
elements. A set of fictitious tractions are assumed at the pile-soil interface. Butterfield and
Banerjee (1971) showed that these fictitious tractions are identical to the real ones for pile
slenderness ratio greater than five. Each boundary element is associated with unknown
2-10
tractions and displacements. Some of these are known over parts of the boundary, the rest of
them are computed using point-force singular solutions as proposed by Kelvin (Sir W.
Thomson, Cambridge and Dublin Math. J. (1848 or Mindlin (1936). Once these boundary
values are obtained, the displacements and tractions at any point inside the domain, satisfying
the boundary conditions at the soil-pile interface, are computed using the boundary integral
equation, which is the basis of the boundary element method. This method has also been used
to evaluate pile group behaviour. It has the advantage over finite element analysis that only
the boundaries need to be discretised and the singular solutions used in the analysis
automatically satisfy the boundary conditions at the free surface and at infinity. The number
of equations needed to be solved is thereby reduced, making the method computationally
efficient. Various computer programs based on the boundary element method are currently
in use. The Department of Transport's computer program PGROUP was developed by
Banerjee and Driscoll (1976). It is capable of handling pile groups of any geometrical shape
subjected to combined horizontal, vertical and moment static loading. It has since been
upgraded to include a soil model with a linearly increasing modulus following the work of
Banerjee and Davies (1978). Banerjee and Davies (1978) described a non-linear method of
analysis in which volume cells were introduced into the soil domain to handle soil yielding.
In this analysis procedure the boundary element method loses its inherent computational
efficiency over the finite element method i.e. loses the merit of a smaller system of algebraic
equations than is needed in the finite element method. Another computer program DEFPIG,
less rigorous than PGROUP, written by Poulos (1980), is based upon a simplified boundary
element approach for a single-pile analysis and the calculation of the interaction factors for
two equally loaded identical piles. Soil non-linearity is modelled by limiting the stresses at
the pile-soil interface, while soil inhomogeneity is approximated by an averaging procedure
using the point load solutions of Mindlin (1936). The PIGLET computer program (mentioned
earlier, Section 2.1, and developed by Randolph (1980) is based on analytical solutions that
are either derived theoretically or fitted to the finite element results to give the response of
single-piles. Pile-soil-pile interaction is based on interaction factors determined from
2-11
expressions fitted to the results of the finite element analyses. The interaction-factor approach
used in DEFPIG and PIGLET for pile-soil-pile interaction gives loads and bending moments
only at the pile heads, not their distributions along the piles. Limitations of the elastic
approach to pile interaction under lateral loading have been revealed from model pile tests
carried out on a centrifuge by Barton (1982). These tests show that the leading piles of groups
embedded in sand carry a higher proportion of overall applied load than the trailing piles. For
pairs of piles at close spacing (s/d=2) ,where s is the centre-to-centre pile spacing, Barton
found that applied load was shared 60% to the front piles and 40% to the rear piles. In
contrast Matlock et al. (1980) reported results of tests on pile groups in clay, where a uniform
sharing of the applied load was found. Clearly, the adequacy of the elastic approach will
depend on a number of factors such as the load level, soil type, and most importantly pile
spacing. Despite the limitations, these two programs are widely used since (apart from there
being no viable alternative) they provide conservative solutions to most in-land in-service
foundations. Springman (1989) developed a computer program SINPILE (some parts of the
underlying theory are adapted from PIGLET) to predict the bending moment and deflection
of piles subjected to active lateral loading by soil. The program was further upgraded into an
interactive version SIMPLE. This program analyses a group of two vertical piles which are
connected with a stiff cap i.e. under similar lateral loading conditions. More recently a
spreadsheet program SLAP has been presented by Randolph and Springman (1991) for the
analysis of pile response due to external loads and soil movements.
2.3 Numerical methods (for non-linear elastic soil behaviour)
Vertical loads
As real soils seldom obey Hooke's law, the methods that fall within group type (a) (see
above, Section 2.2) take account of non-linear soil behaviour. The earliest, the load transfer
2-12
method, given by Coyle and Reese (1966), involves an iterative analysis procedure. The
procedure relies on curves which relate the ratio of load transfer to soil shear strength versus
pile displacement along the length of the pile. The curves are constructed from data obtained
from the field or laboratory-instrumented piles. The pile is divided into a number of small
segments. An initial guess is made for the movement of the segment at the bottom. Then,
with the help of the appropriate load transfer curve, an iterative procedure is followed to
obtain the load acting on the segment. From this information about the load and the
movement of the bottom segment, the analysis then proceeds upwards, segment by segment,
to the top of the pile. Thus, various initial movements at the bottom pile segment are used
to construct the load displacement curve at the pile head. Clearly, such load transfer methods
ignore continuity of the soil medium and cannot be applied to pile groups. Moreover, in order
to obtain load transfer curves on site, considerably more instrumentation is required than for
a normal pile load test. Also, as soil conditions vary, horizontally as well as vertically,
extrapolation of test data from one site to another may not always be entirely satisfactory. The
load transfer method was refined by Kraft et al. (1981) who suggested the use of theoretical
expressions for the data obtained from field tests on instrumented piles to determine curves
relating load transfer t to vertical pile movement z at various points along the pile. These t-z
curves were then used as non-linear supporting soil springs in an analysis of the load
settlement behaviour of the pile using the finite difference numerical technique. Randolph
(1981) developed a computer program RA TZ, based on the load transfer method.
A hybrid approach, which is a combination of load transfer method and the integral equation
method (mentioned earlier Section 2.2), is used to predict the settlement of vertically loaded
pile groups. In this approach, the soil behaviour of individual piles is modelled using load
transfer curves and the pile-soil-pile interaction is determined using Mindlin's (1936) solution
which is the basis for the integral equation method. Among others, Chow (1986) has
demonstrated the success of this approach by presenting a close match between the predicted
results and the measured results from the field load tests.
2-13
Chapter 2 Literature R e ~ ' i e w
Lateral loads
The subgrade reaction method for the analysis of laterally loaded piles is improved by
introducing non-linear reaction-deflection p-y curves in the analysis. These curves are derived
from actual field measurements which represent non-linear stiffness of soil springs. The
method is generally known as the p-y curve method, developed by Reese and co-workers
(1974,1975), who obtained a finite difference solution to the following equilibrium equation
for a beam on an elastic foundation:
EI d
4
y + P d
2
y P = 0
dz 4 z dz 2
where
y == lateral movement
EI = bending rigidity
z depth
p = soil reaction per unit length
P
z
axial load at depth Z
The finite difference form of the equation and a full description of the resulting equations
after applying the boundary conditions were presented by Reese (1977) who also described
a computer program that solves these equations. The solution requires input of a series of p-y
curves at various depths along the pile. Essentially, an iterative solution procedure is adopted
2-14
as the curves are non-linear: for each iteration the values of the lateral soil reaction are
updated from the corresponding p-y curves by considering the calculated deflections in the
previous iteration until convergence is achieved. The iterative procedure is also needed in
order to have compatible deflections at various levels along the length of the pile.
Design procedures for constructing p-y curves from the results of field measurements on full
size instrumented piles have been developed by Matlock (1970) for soft clays, by Reese et
al. (1974) for sands, and again by Reese et ai. (1975) for stiff clays. To construct these
curves, the ultimate lateral soil reaction per unit length of the pile is required in order to work
out soil stiffness before and after the yielding of soil. There is no consideration of the soil
flow rule and plasticity. In a unified approach, presented by Sullivan et al. (1980), the linear
and non-linear parts of the curves are given by two different empirical relations, slightly
different from those proposed by Matlock (1970). These involve the selection of two
empirical parameters based on intuition and engineering judgement.
The limitations of the p-y curve method are:
1) Difficulties exist in choosing appropriate p-y curves for a given combination of
pile size and soil type.
2) Replacement of the soil continuum by discrete springs restricts the extension of the
analysis to pile groups since interaction between the neighbouring piles may not
be taken into account, when constructing the p-y curves.
Although a fairly successful guess of the initial displacement field is important to obtain a
convergent solution. However, algorithms can be derived which give a convergent solution.
An improvement was suggested by Georgiadis and Butterfield (1982) which introduces a
2-15
Chapter 2 Literature Redew
shear coupling between the springs. The coupled spring model, known as a Pasternak
foundation model (Pasternak (1954)), was originally developed for single-piles. Attempts to
extend it to the analysis of pile groups have not enjoyed any significant success.
O'Neill et ai. (1977) proposed a hybrid model in which the response of a single-pile is
evaluated using p-y curves, and the interaction between the piles through the soil is found
using Mindlin's (1936) solution. An iterative solution procedure was adopted to respond to
the non-linearity of the p-y curves. Later, a similar approach for the analysis of offshore
battered pile groups was presented by Poulos (1980). He considered piles of equal length and
uniform diameter. Leung and Chow (1987) studied pile interaction effects using a hybrid
approach. Beam finite elements for piles and Mindlin's equation for an elastic continuum
were used to model pile-soil-pile interactions. They demonstrated that their method compares
well with the commonly used interaction factor approach for pile groups in a homogeneous
elastic half space. They also modelled field studies made on two offshore pile groups (five-
pile and ten-pile) by Matlock et al. (1980). However, for the ten-pile group, their method
tended to over-predict the load carried by the piles when compared with field test results -
perhaps because of the overlapping of plastic zones around the individual piles. Closely
spaced group interaction cannot therefore be modelled using Leung and Chow's method. More
work is needed to investigate the over-predictive effect further.
A computer program is given in the textbook by Smith and Griffiths (1988) implementing
Leung and Chow's method. Two other widely used computer programs are SPLICE developed
by Clausen et al. (1984) and the W. S. Atkins program SPLINTER. A third, more recent
program, FUGROUP developed by Aldridge and Erbrich (l990), eliminates the inconsistency
in using a constant shear modulus across a layer boundary as computations are made
independent of the shear modulus.
Muqtadir and Desai (1986) described a finite element analysis of a group of battered piles in
which they used a hyperbolic soil model (Duncan and Chang (1970 for sands to incorporate
2-16
material non-linearity. They adopted a non-linear incremental method for the solution of
equilibrium equations. The load is divided into small increments and tangent stiffness defined
by a hyperbolic fonnula. computed after each increment, and used for the succeeding
increment. Although a hyperbolic soil model does not reflect the failure behaviour of the
soils, it does model pre-peak behaviour of both loose and dense sands and also of clays.
2.4 Numerical methods (for eJasto-plastic soil behaviour)
The methods so far described cover linear and non-linear elastic soil behaviour. Kooijman and
Venneer (1988) presented a quasi three-dimensional pile foundation analysis in which, instead
of using p-y curves, they employed a two-dimensional finite element analysis for an elasto-
plastic soil behaviour to calculate lateral soil stiffnesses i.e. theoretical p-y curves using finite
element analysis. These p-y curves take account of interaction between neighbouring piles.
This analysis is based on the assumption that the shear strains in the soil along the vertical
plane of the pile-soil model are zero, whereas the horizontal layers of the soil are modelled
as plane stress slices. This assumption has not been checked and needs numerical
experimentation to test its validity. The analysis was later extended to cover pile group
response by Kooijman (1989).
A three-dimensional finite element analysis for pile group response is computationally
expensive. Lane and Griffiths (1988) carried out such an analysis economically for a laterally
loaded single-pile using a finite element hamlonic analysis. The analysis utilises the advantage
of axial symmetry to perform two-dimensional modelling of the three-dimensional problem
with non-symmetric loading applied in the third dimension in the form of Fourier harmonic
tenns. (This analysis is similar to the analyses done by Kuhlemeyer (1979) and Randolph
(1981) except that the soil is modelled as an elasto-visco-plastic material.) The technique is
an efficient way of handling large three-dimensional problems by reducing computer storage
requirements and running time. Brown and Shie (1990a) performed a realistic three-
2-17
Chapter 2 Literature Reriew
dimensional analysis using a commercial finite element-code, ABAQUS. They modelled the
gap formation and slippage between the soil and pile near the pile head using interface
elements. Two elasto-plastic soil models, the von Mises yield criterion (Calladine (1985)) for
undrained clays and the extended Drucker-Prager plasticity model (Drucker and Prager
(1952)) for sands, were used. Brown and Shie (1990a) reponed numerical problems in using
the latter soil model for purely frictional soils. They managed to eliminate the difficulty by
specifying a small value of cohesion in the model. The finite element model boundaries were
chosen to pass through the pile centre-line and through the soil midway between the piles.
The results were processed to construct p-y curves and a close match with the experimental
curves was shown. The analysis does not strictly correspond to that for a single-pile as the
mesh boundaries are chosen in slich a way that the presence of neighbouring piles are bound
to contaminate the finite element results. A similar type of model for the behaviour of a pile
group has also been presented by Brown and Shie (l990b).
Use of the elasto-plastic boundary element method to study pile group behaviour has not been
actively considered by researchers working in this area. The one reported study, by Banerjee
and Davies (1980), contains very little detail on the actual methods used to incorporate plastic
behaviour in the boundary element method. It is thought (for the reason stated earlier, Section
2.2) that this is not a worthwhile method to pursue.
2.5 Plane strain solutions for laterally loaded piles
In addition to the numerical solutions described above, plane strain analytical solutions for
laterally loaded piles have been derived which provide more quantitative information about
pile and soil behaviour. Baguelin et al. (1977) presented an analytical solution for a circular
disc (two dimensional plane strain model) displaced laterally in an elastic soil. In their
analysis, they justified plane strain conditions to hold near the pile because the displacements
of the various sections of the pile are not arbitrarily independent but vary in a highly
2-18
continuous manner, since the pile's radius of curvature is always very large with respect to
its diameter.
The displacement of the pile section caused by a lateral loading, given by Baguelin et al.
(1977) solution, depends on the choice of outer radius R of a two-dimensional model. The
solution has an ambiguity that for a realistic value of R to be infinity, the displacement of
the pile section in tum tends to infinity, which is unrealistic. It was therefore necessary for
them to determine the value of R for the two-dimensional model used in their theoretical
study, such that the analytical solution would give realistic displacement values. They carried
out a peculiar analysis in which deformation fields from the two-dimensional analysis and a
simplified three dimensional analysis were matched. The analysis is based on the assumption
that there is an annular zone of radius ra in which deformation field follows the law of plane
strain and outside this zone the behaviour of soil is three-dimensional. For the annular zone,
strains were calculated by differentiating the analytical solution with respect to r. For the
outside the zone, they considered a simplified three-dimensional problem of a vertical axis,
loaded horizontally in a semi-infinite elastic solid. The three-dimensional strains were
calculated using Mindlin's (1936) solution for a point source. The distribution of load along
the axis was calculated from the theory of subgrade reaction discussed earlier in Section 2.1.
For short piles, a constant soil stiffness with depth and for long piles, variable stiffness along
the length of the pile was assumed. Using a trial and error procedure, several load
distributions were tried and the value of R was chosen in such a way that the displacement
at r", as calculated by the two- and three-dimensional models was the same. The influence
of the pile section (square pile), disturbance of the soil around the pile due to pile installation,
and plastic yielding of the soil in undrained behaviour, based on Tresca's yield criterion
(Calladine (1985, were studied using finite element method. Graphs were produced to
calculate the corresponding corrective terms to be applied to the analytical solution. Baguelin
et al. (1977) show that such influences are confined to within two pile diameters of the
circular disc. The analytical solution raises doubts about the question of modelling soil
2-19
behaviour as linear elastic, under lateral loading in a two-dimensional model.
Randolph, Carter and Wroth (1979) studied the effects of pile installation and subsequent
consolidation on the performance of driven piles in clay. They modelled pile installation as
an undrained expansion of a cylindrical cavity. The assumption is valid for deep penetrations
when the soil predominantly displaces radially. They considered a work hardening elasto-
plastic material behaviour based on the Modified Cam-Clay model (Roscoe and Burland
(1968)). The excess pore pressures generated in this process were subsequently assumed to
dissipate by means of outward radial flow of pore water. Analyses of the stress changes due
to the expansion of the cylindrical cavity and the subsequent consolidation of the soil with
the dissipation of pore water pressure were presented by Carter, Randolph and Wroth (1979).
Randolph, Carter and Wroth (1979) carried out a parametric study and gave the effects of past
stress history of the soil on the stress changes due to pile installation. In another paper,
Randolph and Wroth (1979b), produced an analytical solution for the consolidation of the soil
around a driven pile. They considered an elastic soil behaviour under plane strain conditions
in the derivation of their solution. The solution gives pore pressure and stress changes in the
soil around an impermeable rigid pile, starting from an initial pore pressure distribution due
to the expansion of a cylindrical cavity in an ideally elastic, perfectly plastic materiaL The
resulting logarithmic variation of excess pore pressure with radius is considered to be close
to that generated around a pile as a result of driving.
These useful studies are related to axial loading of single-piles. No such study has been done
for the case of lateral loading of piles. That is why the main emphasis in this thesis is on
laterally loaded piles and pile groups (Chapters 4, 5, 6 and 7).
When a pile is loaded laterally, two types of soil behaviour are expected, one near the surface
close to the pile head where a wedge type of failure is assumed, and the other at some depth
from the surface where the soil is assumed to flow around the pile at failure. In an offshore
pile foundation, the latter type of failure is more likely to occur than the former. As described
earlier, Randolph and Houlsby (1984) carried out a plane strain analysis of the problem using
the method of characteristics. The soil was modelled as a perfectly plastic cohesive material.
2-20
They presented an exact solution (i.e., one in which both upper and lower bound closed form
solutions match) for the ultimate lateral resistance of the soil. It was found that if the load
was non-dimensionalised with respect to the soil strength, the load factor varied between
(6+1t) for a perfectly smooth pile and (4{i + 2IT) for a perfectly rough pile. The solution
provides information about the size of the expected plastic zone in the soil. A limitation of
the method, however, is that no information is obtained on the displacements needed to reach
the limit load. An understanding of pile-soil-pile interaction phenomena is greatly helped by
determining the ultimate lateral resistance of a soil in a pile group subjected to an external
load under undrained conditions. (A detailed study done on these lines is presented in
Chapters 5 and 6.)
2.6 Summary
Elastic methods are widely used to model the behaviour of inland/onshore pile foundations.
Non-linear methods, which mainly rely on empirical curves based on limited field data, are
used for offshore pile foundations. For closely spaced pile groups in soft clay, the use of these
methods requires great care as it is difficult to predict the effects of interaction on pile group
response. The interactions are highly non-linear due to the overlapping plastic zones around
the individual piles. However, for the case of vertically loaded pile groups, the plastic zones
around the individual piles are not as large in size as in the case of laterally loaded piles.
There is less likelihood that these zones will overlap. The hybrid method presented by Chow
(1986) based on elastic interaction amongst the piles seems quite suitable to predict the
response of vertically loaded pile groups. The use of the boundary element and finite element
methods enables the effects of various soil models to be studied. Application of the boundary
element method to study the elasto-plastic soil-structure interaction has not so far been very
well explored. Current finite element approaches to the problem are computationally
uneconomical, especially for group behaviour. For that reason, a simplified approach for the
2-21
Chapter 2 Literature Rcriew
analysis of laterally loaded piles and pile groups similar to that developed by Kooijman and
Vermeer (1988) is proposed in Chapter 6, in which the ultimate capacity of a single-pile or
a pile group is determined using plane strain finite element analysis as a first step. This
predicts ultimate behaviour more accurately, more or less eliminating the degree of
uncertainty. It therefore helps in optimising foundation design. Then, as a second step,
theoretical p-y curves are determined from a single plane strain finite element analysis of
single-piles or pile groups along their lengths and used as input to non-linear beam-column
equations (see above, Section 2.3) to determine the pile-head displacements. Although the
method is not ideally successful at modelling field tests (see Chapter 6), because of its
simplification of the three-dimensional nature of the problem, it is nevertheless, in the absence
of a complete three-dimensional analysis, and for other reasons to be explained below, the
best method available.
2-22
Chapter 3 Elastic analysis of a single-pile under axial loading
CHAPTER 3
ELASTIC ANALYSIS OF A SINGLE-PILE UNDER AXIAL LOADING
3.1 Introduction
This chapter proposes a new analytical method for calculating the pile head settlement of a
rigid pile. An approximate analytical solution is derived. Its validity is checked using the
integral equation method ClEM). The new analytical method is further adapted in the
development of a correcting layer approach to model the infinite boundaries in the finite
element analysis of a single-pile. The correcting layer approach is computationally economical
to perform the finite element calculations. The response of a single flexible pile under static
axial loading is analysed using the finite element method (FEM) and the IEM and the results
obtained are compared. The ability of the correcting layer approach to model the infinite
boundaries is demonstrated by achieving a close agreement in the results from both of these
methods.
Figure 3.1 shows a pile driven into a linear elastic homogeneous soil. The following
conditions are assumed: there is a perfect bond at the pile-soil boundary (i.e. no relative slip);
and there are no shear stresses acting at the pile base-soil interface (i.e. there is a smooth pile
base).
The aim of this simplification of the initial conditions before the analyses are undertaken is
to focus exclusively on the problem of determining, from the methods currently in use, the
numerical method best suited to the analysis of single-piles. Two out of the various numerical
3-1
p
free surface
flexible pile (E. A )
1
linear elastic homogeneous soil (G.v)
Figure 3.1: A single-pile under axial loading
methods described in the preceding chapter, namely the IEM (also known as quasi-boundary
element method) and the FEM, are used in the analysis following.
A computer program was developed to perform the integral equation analysis for a single
vertical pile under axial loading in an ideal linear elastic homogeneous soil medium. The
numerical formulation of the analysis, presented in Section 3.2 is in some part adapted from
an early work by Randolph and Wroth (1978). In this paper they carried out both the finite
element and the integral equation analyses for an underreamed pile, and compared the results.
Since the FEM requires fixed boundaries to analyse a problem unlike the IEM which does
not require any, they artificially introduced a fixed boundary in their integral equation
formulation at the same depth below the pile base where they had a fixed boundary in their
finite element analysis in order to make a satisfactory comparison, adopting the rigid layer
model described in Butterfield and Banerjee (1971). In this model fictitious stresses are
introduced over the layer to give zero overall vertical displacements at that depth. Randolph
and Wroth (1978) used this model to test their finite element formulation.
The principal aim here is to model infinite boundaries in the finite element analyses while
3-2
keeping the integral equation analysis as simple as possible. However, Section 3.3 goes on
to describe the possibility of extending this analysis into a general purpose two-dimensional
elasto-plastic boundary element method, and the reasons why it is not pursued further in this
thesis.
A new analytical method is proposed in Section 3.4 which provides an improved solution for
the pile head settlement of a rigid single-pile in an ideally elastic soil medium. The method
adopts a load transfer mechanism such that the solution does not have to rely on a fictitious
radius, r
m
, as in Randolph and Wroth's (1978) solution. The proposed method suggests
further a correcting layer approach to model the infinite boundaries in the finite element
analyses. The finite element analyses were carried out using the commercial package
ABAQUS, and the results obtained are compared with the integral equation analysis described
in Section 3.2. The computer program, performing the integral equation analysis, is
significantly adapted from Randolph and Wroth (1978) in order to serve the particular
conditions of this study. The effectiveness of the correcting layer approach is demonstrated
by the close agreement achieved in the results from these two numerical methods (i.e. IEM
and FEM).
The correcting layer approach also has the merit of achieving a very significant saving of
computational effort in the finite element analyses by making meshes almost half the size they
would otherwise have to be (Randolph and Wroth (1978); Baguelin and Frank (1980. The
saving naturally encourages the use of FEM to study the response of not only single-piles but
also closely spaced offshore pile groups under vertical loading, whereas !EM is incapable of
modelling the non-linear pile-soil and pile-soil-pile interactions.
3.2 Axi-symmetric quasi constant boundary element formulation
It is important to appreciate how the integral equation analysis presented here differs from
that of Randolph and Wroth (1978) (see Figure 3.2 below). The pile shown in the figure is
3-3
Chapler3 Elastic analysis oj a single-pile under axial loading
axis of
~
symmetry
SOre
s.oilmass
J: __ J ~
free sanding
1 ~ J I I ~
pile
1"""---1
1
I "
J L _ - 1 ~ ).
it '
1
1/ i
, __ r -, -:"l
rn .i*Hnn
v;
M
A
(a) (b)
- - ~ ___ L-=-rw,

.
I
r - - - ..,
", poinu at which

I
pile I
di'pl.",menu are
'---""'I

,
e&leul.ted
, ,
;- - -,

,
l ___ J

.
,
fheel boundary
.:.=a::J-. _c:::::::'
(e)
(d)
Figure 3.2: Finite element grid and schematic diagram of integral equation method
(after Randolph and Wroth (1978)
underreamed. Various stages involved in the integral equation analysis of the pile are drawn
next to the finite element mesh to model the same pile in the finite element analysis.
Randolph and Wroth (1978) were obliged to consider a rigid layer at some depth below the
pile base in order to compare satisfactorily the results with the FEM. In the present analysis
we need not worry about the rigid boundaries as these are taken care of by adopting the
correcting layer approach being presented in Section 3.5 for the finite element analyses. It is
worth remembering that use of the IEM is reliable only under linear elastic conditions and
is not suited to more complex soil models which are required for detailed analysis. It is for
this reason that the present integral equation formulation is kept simple, and a straight (Le.
not underreamed) pile is analysed, and the effect of a rigid layer is not considered.
The analysis is done using an axi-symmetric quasi-constant boundary element computer
program specifically designed for elastic analysis of a single-pile under axial loading. The
numerical formulation is based on the solution for a point load inside a semi-infinite
homogeneous elastic half space derived by Mindlin (1936). The solution gives stresses and
3-4
Chapter 3
plane z = 0
h
Elastic analysis of a single-pile under axinlloading
mirror image
:r
(O,O,-c)
free surface
;It'P\\R
2
__
load \ (x,y,z)
I R '-. i I
i field point
Z t
x
Figure 3.3: Force in the interior of a semi-infinite elastic half space (after Mindlin
(1936)
displacements at a field point (x,y,z) due to a point load at coordinates (0,0, +c) as shown in
Figure 3.3.
The analysis is based on the following formulation. Figure 3.4(a) shows a cylindrical column
of soil consisting of n elements. These are described as shaft elements. A set of shear stresses
are imagined to be acting at the boundary of this column and vertical normal stress at its base.
The shear stresses are assumed to be constant over each element of pile-soil interface which
is discretised into a number of constant strain boundary elements. The greater the number
of these elements, i.e. the finer the discretisation of the interface, the more accurate the
solution will be. The other assumptions made in the analysis (as mentioned at the beginning
of the chapter) are: a perfectly rough pile which makes a perfect bond with the soil, and
which has a perfectly smooth base. The compatibility of radial displacements is not taken into
account. Poulos and Davis (1968) reported that the radial stresses on the piles are only of the
order of 0.5% of the shear stresses on the pile. Mattes (1969) found that the distribution of
shear stresses and the settlement of the pile, considering both vertical and radial displacement
3-5
Chapter 3 Elastic analysis oj a single-pile under axial loading
applied load (P)
free surface
applied load (P)
ele-l eJe-l t
w
I
1 1
J J
t
soil continuum
I I I
1 1
c1e-n
t
(J ('i
b b
(a) (b) (c)
Figure 3.4: Schematic diagram of the integral equation method
compatibility, is almost identical to that obtained by considering vertical displacement
compatibility only. Butterfield and Banerjee (1971) also found that the additional radial
compatibility condition in the boundary element fom1ulation did not affect the resulting
vertical displacements for a given load.
Each shaft element exerts a constant shear stress on the soil. The stress exerted by j Ih
element is 1), while the base element exerts a normal stress (Jb on the soil underneath it.
Making use of Mindlin's solution, the vertical displacement may be calculated at any point
in the soil continuum by integrating the effects of all the stresses acting on the shaft elements.
In particular, the displacement of the mid-point of each shaft element can be calculated.
Though the discretisation suggests a linear variation of displacement within the element, but
for the purpose of this part of the calculation, the value of displacement at the mid-point is
assumed to apply.
For the pile in Figure 3.4(a), the displacement wJ in a shaft element i due to the traction
at element j is given by:
3-6
Chapter 3 Elastic analysis of a single-pile under axitliloading
(3.1)
where f is an influence coefficient for the shaft element i, determined by integrating
lJ
Mindlin's equation over the element j .
The total displacement Wi S (at any depth) is the sum of individual displacements caused by
the shear stresses acting on each of the n shaft elements and the normal stresses (5b at the
base elements. This can be expressed as:
n
s
Wi L fJ 't
J
+ fb (5b
(3.2)
;01
where the subscript b represents the base element. The influence coefficient for shaft element
i, determined by integrating Mindlin's equation over the base element is fb'
Similarly, the settlement Wb$ at the base is:
n
wb
s
= L fbj 't
j
+ fbb (5/J
i=1
(3.3)
where fbj and fbb are the influence coefficients for the base element integrating Mindlin's
equation over the shaft element j and the base element itself, respectively.
Combining equations (3.2) and (3.3), we can write the n + 1 equations in matrix form as:
3-7
w" = [f"}, "( (3.4)
where W
S
and 't are the vectors of displacements; shear and normal stresses in the shaft and
base elements respectively as shown in Figure 3.4(a), The displacements and the stresses are
linked together through the soil flexibility matrix [f"].
The pile is discretised into exactly the same nllmber of pile elements as the shaft elements.
We will now establish a link between the displacements in the pile and the shear tractions by
satisfying the compatibility condition that, at the pile-soil interface, displacement in the pile
and the soil element will be the same. Figure 3.4(b) shows a free body diagram of the pile
carrying a vertical load P and the fictitious soil tractions acting in the opposite direction. The
applied loading will produce elastic shortening of the pile which is determined from the
theory of elasticity, considering each pile element as a spring of axial stiffness K=EA/l.,
where A is the cross-sectional area, E is the elastic modulus, and ie is the length of the
element. Each of these stiffnesses of the pile elements are assembled in a stiffness matrix [K].
The equilibrium of the pile can then be represented by the equation:
P - 'tA = [K]. w
P
c
(3.5)
where P , 't and w
P
are the vectors of applied forces and the vector of pile displacements
respectively, as shown in Figure 3.4(b). The circumferential area of a pile element is Ac'
Substituting equation (3.4) into equation (3.5), noting that w S = W P, and then rearranging,
3-8
Chapter 3 Elastic analysis of a single-pile under axinlloading
gives equation (3.6) for the n shear tractions and a normal base pressure:
(3.6)
where the matrix [F] is the inverse of the matrix [Kfs+/].
Shear tractions obtained from equation (3.6) are then substituted back into equation (3.4) in
order to find vertical displacements along the length of the pile. Figure 3.4(c) shows the pile
head settlement of the actual pile under an axial load, installed in place of the imaginary pile
initially assumed in the soil medium.
The main part of the program is the evaluation of the influence coefficients for the soil
flexibility matrix [f5] which relate the displacements at points on the eventual pile-soil
interface to the stresses acting over the interface (details can be found in Poulos and Davis
(1980) and are given for completeness in the appendix at the end of the chapter). The non-
singular integration for the off-diagonal terms of [f'] is carried out analytically in the vertical
direction (shaft elements) and radial direction (base elements) in the manner described by
Poulos and Davis (1968). The circumferential integration is done using Simpson's rule with
an angular step size of (n/50) radians. The integration, both analytical and numerical, required
for the evaluation of the leading diagonal coefficients of the flexibility matrix [f S], needs
special attention and has been carried out as suggested in Randolph (1977). These coefficients
relate the displacement of the same point at which the load is assumed to be acting i.e. when i =j
in Figure 3.4. The analytical integration is carried out in two parts leaving a tiny gap
equivalent in size to a very small fraction (.01 %) of the pile radius, at the centre of the
element. The numerical integration is performed at intervals of (n/60) radians over the
3-9
Chapter 3 Elastic analysis of a single-pile under axinlloading
majority of the pile circumference but with a much finer step length in the vicinity of the
singularity. This scheme excludes from the integration a square region of the side, 0.000 1
times the pile radius, thus eliminating the point singularity (see Figure 3.5). For the numerical
integration of the base elements, intervals of (n/50) are used.
I
e
I- 2'
I ~ e = 0.00(11 a
I - - - - - - - . . . . . . J - j ~ half of the square region representing a singularity
fmer step length
i;..
2
8
L to 1l:
2 2
Figure 3.5: Typical pile element showing square region excluded from the
numerical integration
The results obtained from the integral equation analysis are shown in Figures 3.6 and 3.7.
These figures show the distribution of shear stresses along the length of the pile both for the
case of a flexible pile and a rigid pile. Results from the integral equation analyses done by
Mattes and Poulos (1969) for a flexible pile and Poulos and Davis (1968) for a rigid pile are
also plotted in these figures. The results match well with the results of the present integral
equation analysis. (E/Es is the relative pile-soil stiffness or flexibility ratio.)
3.3 Two-dimensional elasto-plastic boundary element method
There is a body of opinion that holds that the boundary integral equation method can be
satisfactorily extended to enable better soil modelling. This would mean that Randolph and
3-10

3.5
0.5

II
0,2 0.4 0,6
depth Z
\1aucs and Poulos (19(9)

d
v = 0,5
E
P
-=50
E
0,8
!EM
Figure 3.6: Shear stress distribution along the length of a flexible pile
Q,9

0,\
.!-.-= JOO
d
v =0.5
E
P
E
s
10000
0,2 0.3 0,4 0,5 0,6 0.7 O,g 0.9
nondimensionalised depth !:....
I
JEM Poulos and Davis (1968) ..
Figure 3.7: Shear stress distribution along the length of a rigid pile
1.0
Wroth's (1978) integral equation formulation is not limited to ideal linear elastic soil
3-11
conditions but adaptable to elasto-plastic soil modelling by extending it to the boundary
element method. But, as we shall see, the adaptation of this method is computationally
expensive and yields results no better than the finite element method.
The boundary element method forms a flexibility matrix which has to be inverted and the
resulting stiffness matrix is fully populated. Therefore the numerical solution procedures are
applied to the whole matrix. The finite element method is based on direct assembly of the
stiffness matrix which is symmetric and banded in nature. The advantage of having a banded
matrix allows the use of efficient Gaussian elimination solution schemes. The boundary
element method is based on singular solutions which require superposition and work well with
elasticity. These solutions are irrelevant to non-linearity. The finite element method can be
readily extended to non-linearity. It is difficult to model infinite boundaries in the finite
element method whereas, in the boundary element method, the condition of infinite
boundaries is automatically satisfied. (A correcting layer approach to model infinite
boundaries in the finite element method is presented in the following Section 3.5.)
Swedlow and Cruse (1971) presented a modified boundary element formulation in which they
introduced fictitious body forces in such a way as to satisfy the yield and flow equations of
plasticity. The resulting new equation, no longer a boundary integral equation, was presented
in an incremental form to take account of material non-linearity. The algorithm required
surface integrals as well as area integrals due to the discretisation of expected plastic zones
in the body of the soil layer into a number of cell elements. Later, Telles et at. (1981)
outlined a procedure to incorporate plasticity, which is the same as the finite element method,
but entailing rather complex computer coding. The numerical procedure for the solution of
the resulting system of algebraic equations cannot be made any more efficient because the
stiffness matrix obtained is very sparse and non-symmetric. The fact that the size of the
stiffness matrix will be much smaller than the one arising in the finite element method does
not make any substantial difference to the problem under consideration here. In both methods,
an almost similar mesh is required since, in the elasto-plastic boundary element method, in
addition to the surface discretisation, domain discretisation of the expected plastic zones
3-12
around the pile is also necessary. By contrast, a carefully graded finite element mesh, having
large finite elements outside the expected high stress concentration areas, results in nearly the
same size stiffness matrix as that which arises in the boundary element method which is also
non-symmetric and sparse. The disadvantage of having to handle a larger stiffness matrix in
the finite element method, compared to the boundary element method, for most geotechnical
problems, is outweighed by the efficient numerical solution procedure of the resulting
algebraic equations, especially when dealing with non-linear problems. Moreover, numerical
evaluation of the terms of the stiffness matrix in the boundary element method is not
straightforward. For example, in order to eliminate the singularity in a triangular element, the
number of integration points near the singularity is increased and each element is subdivided
into a number of sub-elements. Although these sub-elements do not create any additional
nodes in the final system of equations, this procedure does make the computer coding quite
complex and cumbersome. There are, besides, other problems like anisotropies, radial and
vertical inhomogeneities in the soil medium, which remain to be resolved. All of these
reservations combine to discourage the pursuit of the elasto-plastic boundary element method.
3.4 Analytical methods for the analysis of rigid single-piles
Study of the incorporation of an elasto-plastic soil model in the boundary element method
proves that the method to be less efficient as compared to the finite element method. There
is no significant advantage in using an elasto-plastic boundary element method, not only
because of the inefficient performance of the numerical computations involved in the method,
but also because the element discretisation of the problem is more or less the same as that
required in a finite element model. On the other hand, in the finite element method, the
condition of infinite boundaries is not automatically satisfied as is the case in the boundary
element method.
In this section, we describe a new analytical method for the analysis of a rigid single-pile
under axial loading. The method is a development of an early work of Randolph and Wroth
3-13
Chapter 3 Elastic analysis of a single-pUe under axial loading
(1978).
The analytical solution given by Randolph and Wroth (1978) for the analysis of a rigid single-
pile under axial loading is based on the idea of a magical radius r m where the displacements
become zero. This parameter was artificially introduced to ensure analytical tractability. The
new method, described in Section 3.4.2 adopts a load transfer mechanism which eliminates
the need to have a magical radius r m and provides an improved analytical solution. Also, we
will see later that the method is further adapted in the development of two correcting layers
around the finite element mesh of a single-pile to model infinite boundaries in the finite
element analysis. The success of modelling infinite boundaries using correcting layers is
demonstrated for a single-pile but it can equally be applied to closely spaced pile groups
assuming a pile group to be of an equivalent large diameter single-pile. In order to explain
the work involved in the development of the new analytical method, it is necessary to
describe some of the early work on the analysis of a rigid single-pile, and, more particularly,
the work of Randolph and Wroth (1978).
3.4.1 Simplified analysis of a rigid single-pile (after Randolph and Wroth (1978
The analysis is based on the separation of loads carried by the pile shaft and the base. The
approach to the analysis is shown in Figure 3.8(a), where the soil has been divided into an
upper and a lower layer by a horizontal plane AB, at the level of the base of the pile. It is
assumed that the upper layer of the soil will be deformed exclusively by the load shed by the
pile shaft and that the lower layer of the soil will be deformed exclusively by the pile base
load. Figure 3.8(b) shows the separate deformation patterns anticipated. The plane AB has
been "exploded" to AIBI and A
z
B
2
The deformation patterns along AIBl and Aj32 will not
be compatible and this will lead to some interaction between the upper and lower soil layers.
The deformations in the soil around the shaft are modelled as shearing of concentric soil
cylinders.
3-14
A. ------- --
(a) Upper tJ1d lower soil layers
P
t
p

piie
layer
of sail
-------- --- B
lower layer
of soil
(b) Separa':.e deformation panerns of upper and lower layers
(0)
(j
r
e
-.'
01.
(d)
Figure 3.8: Analysis of pile by separation of loads carried by shaft and base (after
Randolph and Wroth (1978
When the pile is loaded, the increase in shear stress, 't, in the vicinity of the pile shaft, will
be much greater than the increase in vertical stress, O'z' Therefore, the vertical equilibrium of
an element of soil shown in Figure 3.8(d) can be expressed as:
... 0
ar
(3.7a)
ignoring the ()O'jaz) term. If the shear stress 'to at the pile face r = r
o
' the above expression
can be integrated to show that the shear stresses around the pile decrease in inverse proportion
to the radius, such that 'tr ... 'toro is approximately constant. (This relationship is valid only
under the assumption that the soil around the pile shaft will deform independently from the
rest of the soil below the level of the pile base i.e. the deformations in the soil, around the
shaft and below the base, will remain confined within distinct horizontal boundaries. This
seems to be an unrealistic assumption.) Using the elasticity theory, the shear strain will be:
3-15
Chapter 3 Elastic analysis oj a single-pile under axial loading
1: au
Y = - =
G az
aw
+-
dr
in which u and ware the radial and venical displacements of the soil respectively, and G is
the soil shear modulus. The primary displacement will be vertical, and thus, ignoring duraz,
integration of the vertical strains results in the following approximate expression for the
displacement of the pile shaft
w '" '"Coro f dr
s G r
(3.7b)
TO
where, Ws is the settlement of the pile shaft; 1: is the shear stress at any radius f from the
axis of the pile; to is the constant shear stress at the pile shaft; and fo is the radius of the
pile.
Equation 3.7b gives infinite displacements for a rigid pile which is physically untrue. This
ambiguity is a direct result of the improper assumption (the soil around the pile being
separated into two distinct layers which deform independent of each other) made in the
analysis. To avoid this ambiguity Randolph and Wroth (1978) suggested that there is some
magical radius r m at which the displacements become negligible, and that, therefore, the
settlement of the pile shaft, w
s
' is given by:
(3.8)
or
w =r 'taro
S "::> G
3-16
where
The parameter, r m' is determined empirically and as we shall see in the next section,
theoretically r m does not have a unique value but varies with the depth of the pile. Randolph
and Wroth (1978) adopted a constant value of rmand suggested an empirical relationship to
determine this value by matching the analytical solution with the results of their integral
equation analysis.
---
-.
'. "'"
,
.
.
,
.
.
,

free surface
f-ot--pile
. '
, '
. '
. '
\"" ~ ", .. /'
~ ~ ~ , ~ - ~ ~ ~ ~ #
.... _-.---.... , .. _- ......... - ......
./"', ..
,
,
.
,
,
.
.
-
Figure 3.9: Hypothetical variation of limit of influence of pile (after Randolph and
Wroth (1978))
3.4.1.1 Effect of lower layer of soil on upper layer
The deformations in the lower layer are assumed to be those arising from a rigid surface
punch (the solution was given by Boussinesq (1876) for a smooth base with non-uniform
pressure distribution at the base-soil interface (see Poulos and Davis (1974)), At some large
distance from the pile the base will appear as a point load, and thus the deformations ofA:fl2
3-17
will decrease approximately inversely with r. The defonnation of Aj32 (Figure 3.8(b will
decrease more rapidly than the logarithmic variation implied in equation (3.7b). Thus, the
lower layer will act as a restraint to the defonnation of the upper layer, causing vertical stress<J.
increments in it, and contradicting to the assumption made in deriving equation (3.8).
Undoubtedly these stress increments die out at the free surface but, in order to satisfy the
equilibrium equation, an additional non-zero Cacrjaz) tenn has to be included. This means
that the shear stresses decrease, more rapidly than initially assumed, with radial distance along
the depth of the pile; similarly, the radius r m' will not be constant but also decrease with
depth, as shown in Figure 3.9.
The assumption that the deformations in the soil below the level of the pile base will be the
same as arising in a soil beneath a rigid punch seems to suggest that the value of r m will
depend not only on the length of the pile l but also on Poisson's ratio v of the soil. This
assumption brings another ambiguity in equation (3.8) that we cannot have a unique solution
for the settlement of a rigid pile. However, to avoid this problem Randolph and Wroth (1978)
ignored the variation of rm with depth, and suggested an average value of 2.5I(l-v).
3.4.1.2 Modelling the base of the pile
Randolph and Wroth (1978) treated the base of a pile as a rigid punch. Since the available
analytical solution corresponds to a rigid punch at the surface of an elastic half space, they
introduced a depth factor, Tl in order to take account of stiffening effect of the soil above the
level of the loaded area. The load settlement behaviour was described as (Timoshenko and
Goodier (1970)):
p /,(l-v)
Wb = 4r G Tl
o
3-18
Chapter 3 Elastic analysis of a single-pile under axio.lloading
where we is the displacement of the pile base.
The depth factor is the ratio of the stiffness of a rigid punch at the surface of an elastic half
space to its stiffness at some depth below the surface.
Pile base can be modelled as either of these two possible ways. It can be treated as a rigid
circular footing at the bottom of an open borehole or a loaded rigid circular area embedded
in a soil mass at a depth, h. The appropriate model to consider is a loaded rigid circular area
embedded at a depth h. A study conducted by Butterfield and Banerjee (1971) based on this
model shows that the limiting value for the depth factor is 0.6. Randolph and Wroth (1978)
suggested that this factor would be greater than 0.85 or close to unity for a straight (Le., not
underreamed) pile. They argued that there will be no stiffening effect of the soil above the
level of the loaded area because the upper layer of soil is already being deformed by the
action of shear stresses along the pile shaft. The pile is not capable of absorbing any of the
load at its base (as would be a column of soil) since it is the means by which the load is
transmitted to the base. In other words they considered the base of the pile as a loaded rigid
circular area at the bottom of an open borehole. This seems to be a conservative approach.
The load-settlement ratio was finally given by (adding together the contribution of the pile
shaft and the base):
=
P
s
4
+ ----- = ------
Gr W ll(1-v)
o s
where P, is the sum of the contributions from the pile base P b and the pile shaft P
s
'
The difficulty in the choice of a proper model for the pile base only arises when we consider
two distinct soil layers, around the pile, deforming independently and interacting at a common
boundary. In the new analytical method, explained in the following section, we do not have
to worry about the choice of a model for the pile base.
3-19
3.4.2 Simplified analysis of a rigid single-pile: (proposed new method)
In this section a new analytical method is described. The soil is not considered to deform
separately in two distinct horizontal layers. A pyramidal load transfer is assumed below the
pile base which eliminates the need to consider the so-called magical radius rm (at which the
displacement are assumed to be zero). The outer boundaries of the soil are taken out to
infinity. This load transfer mechanism is depicted in Figure 3.10.
free surface
linear elastic
homogeneous soil
t
P
-pile
r
00
" ~ theoretical boundaries ~ "
, within which soil defonns in bearini',.:
Figure 3.10: Load transfer mechanism adopted in the development of new analytical
method
The load transfer line marks the boundary between the soil deforming due to shearing of
concentric cylinders of varying depth (unlike constant depth as in Randolph and Wroth's
(1978) analytical approach) around the pile shaft and the soil deforming in bearing below the
level of the pile shaft. This boundary is not fixed but varies with the assumed load spread
angle. Consider a vertical pile embedded in a linear elastic soil medium as shown in Figure
3.1 L Due to axial symmetry, only half of the pile and the surrounding soil needs to be
shown. Zone A will deform in shear and if we assume that loading of the pile generates
uniform shear stresses along the boundary XY at a radial distance r) from the centre line, and
the soil continues to an infinite radius to the right, beyond XY, where the shear stresses will
3-20
Chapter 3 Elastic analysis oj a single-pile under axwlloading
<t r
p!

f
oc
'(
1
'X
:. /"" fO
--
.. pile
G ,v
0

G ,v
A A
E
,
,
,
,
z D
,
,Y
00
,
,
,
,
G)
,
,
,
EB '
v
,
,
B
,
,
,
,
Figure 3.11: Zones of soil deforming in shear and bearing around a pile
be zero, it is possible to estimate vertical deformation at X as follows:
From the geometry of the problem, integration of the approximate equilibrium equation (3.7a)
gives the shear stress 't at radius r, to be:
't =
where 0 is the load transfer angle.
t/I(i+(rl-rO)tanO)
r(l +(r-r o)tanO)
Applying the condition of elasticity, the vertical displacement at point X can be expressed
in this integral form:
3-21
Chapter 3 Elastic analysis oj a single-pile under axia/loading
Evaluating the integral analytically, the following expression is obtained for the vertical
settlement of point X:
(3.9)
where equation (3.9) represents an elaborated form of the new improved solution for the
calculation of settlement of a rigid pile shaft. If the point X (see Figure 3.11) is located on
the pile shaft i.e. r
1
:: ro and tl := 'to in equation (3.9), then once again we obtain equation
(3.8) for the vertical displacement of a pile shaft with the following improvement:
for 8 = 45.
[-r r
o 0
The improvement suggests that the vertical displacement of a pile shaft depends on the
logarithm of the slenderness ratio of the pile and is not dependent on the empirically
determined magical radius r In'
Now consider a pyramidal load transfer beneath the pile base as shown in Figure 3.11. ZoneR
under the pile base will defonn in bearing. If we assume a uniform normal stress a, at a
depth ZI' below the pile base, then, considering the geometry of the problem and satisfying
the equilibrium condition, the venical normal stress at any depth z, can be approximately
wrinen as:
Gr
2
G(z) '" ____ 1_--:::-
(r
l
+(z-Zt)cot8)2
3-22
(It is important to note that the above expression does not require any choice of model for
the pile base as was needed in the earlier work presented by Randolph and Wroth (1978).)
The settlement of point D, assuming soil of infinite extent, can be estimated by evaluating
the following integral expression:
where E is the stiffness (Young's modulus) of the soil underneath the pile base.
The analytical evaluation of this integral gives an approximate stiffness of the zone B:
(3.1 0)
The vertical settlement of the pile shaft (equation (3.9) and the pile base (equation (3.10
will be the same for the conditions assumed above that the pile is rigid. point X lies on the
pile-soil interface and D lies at the centre of the pile base. If we keep to our assumption that
r
l
=ro and cr=a
o
then adding the contribution of load carried by the pile shaft and the base
from equations (3.9) and (3.10), the pile head load-settlement ratio can be written as given
below:
27tl
= + 21t(1 +v)
~ r o
3.4.3 Comparison of results
Load-settlement versus slenderness ratios for rigid and flexible piles are plotted in Figures
3-23
Chapler 3 Elastic analysis oj a single-pile under axial loading
3.12 and 3.13, along with the results obtained from Randolph and Wroth's (1978) analytical
solution and Poulos and Davis's (1980) integral equation method.
2P
WdG
150
100
50
E
-p = 10000
E
v = 0.3
o
o 10 20 30
o
__ propooedmothod(6 =45 )
o
......... propooed method (6 = 47 )
40 50
21
d
60 70
___ Randolph and Wro\h (1978)
_IEM
80 90 100 110
Figure 3.12: ratios for rigid piles
Figures 3.12 and 3.13 show that the assumed load spread angle of 45, in the present
analysis, provides a stiffer pile response compared to the one obtained from Randolph and
Wroth's (1978) solution, for rigid and flexible piles. As mentioned in Section 3.4.1.2, perhaps
Randolph and Wroth's solution (1978) is conservative due to the model they have adopted
for the pile base. (However, the assumed load spread angle can be adjusted to give a less stiff
response. A load spread angle of 47 provides an exact match with the Randolph and Wroth's
(1978) solution for the rigid piles. It is therefore seen that Randolph and Wroth's (1978)
approach is a special case of the above analysis.) Nevertheless, for the case of flexible piles,
the results from the present analysis match well with Poulos and Davis (1980) results.
In the following Sections 3.5.1 and 3.5.2, the development of two correcting layers to be used
in the finite element analyses, one around the pile shaft and the other below the pile base, is
explained. The properties of the two correcting layers are derived using the above analytical
3-24
Chapter 3
120
100
80
2P
--60
wdG
40
20
0
0 10 20 30
o
__ p,_d melhod ( e = 45 )
~ Pow .. md Davi, (1980)
40
Elastic analysis of a single-pile under axial loading
50
21
d
60 70
-.- R8Ildolpb and W,olh (1978)
_IEM
80 90 100 110
Figure 3.13: Load-settlement ratios for flexible piles
approach. Section 3.5.3 shows how successfully the correcting layers model the infinite
boundaries in the finite element analyses.
3.5 Finite element analysis using correcting layer approach
In this section, the results obtained from the integral equation analysis are compared with the
finite element analysis results. The finite element analysis has a problem of fixed boundaries
which effect the results. We can reduce the effect of boundaries by having softer layers
(correcting layers) adjacent to the boundaries. Although in this case it is difficult to derive
exact expressions for material properties in the softer layers but it is possible to use the
proposed analytical method to obtain approximate values of reduced stiffness which are
thought to reduce the effect of rigid boundaries.
The comparison of results given in Section 3.5.3 shows that two correcting layers, one around
the pile shaft and the other below the pile base, successfully model the infinite boundaries in
3-25
the finite element analyses.
3.5.1 Development of the correcting layer around the pile shaft
Let us consider again the vertical pile embedded in a linear elastic soil medium as shown in
Figure 3.14. Due to axial symmetry, only half of the pile and the surrounding soil needs to
be shown. We know that Zone A of an infinite extent will deform in shear and the vertical
settlement of point X can be approximately written as:
1 ~
Pi .... --------<ool1
x
G.v 0
- ' ~ : J 9 - GA' vA
~ fixed boundary
Icz
l
D
,
E "" ...
" y
,
,
o EB " ~ ~ " "
rigid base
,
Figure 3.14: Shaft and base correcting layers
W ...
x
t{)(l-rotane+r
t
tane)
G
tane
lcote-r
o
3-26
r
1
+lcote-r
o
]
.In ___ _
r)
(3.9)
Now if the value of shear modulus G is modified, making the zone A more flexible than the
rest of the soil, the outer (zero displacement) boundary moves closer from infinity to a
distance krl from the centre of the pile. In other words displacements dampen out quickly by
making the zone softer. This modification leads to a new integral expression for the vertical
displacement of point X:
Analytical evaluation of this integral gives:
(3.11)
Theoretically, Wx from equations (3.9) and (3.11) must be equal as both the equations
represent the same definite integral having an upper limit, a radius where displacements are
zero. If that is so, the value of shear modulus for the zone A can be obtained from the
following expression:
For an assumed value of k, values of G A and VA can be obtained from the above expression.
3-27
Chapter 3 Elastic analysis of a single.pile under axwlloading
There are three finite element meshes used in the analyses presented in Section 3.5.3. Typical
values of G, G
A
, VA for the smallest mesh are 250000 kPa, 79800 kPa, 0.25, respectively,
I
with a load spread angle of 9 '" 45, k = 2, = 40 and r
1
= lOro'
ro
3.5.2 Correcting layer under the pile base
Figure 3.15 shows a pyramidal load transfer beneath the pile base. Zone B under the pile
base will deform in bearing and we know that the settlement of point D assuming soil of
infinite extent can be estimated by the following expression:
arl
Ecot9
(3.10)
If the elastic modulus for the layer B is now modified by making it softer, as we did for zone
A, the boundary at infinity moves to a distance kZI from the base of the pile and the
displacement of point D can be expressed in the integral form as given below:
hI
= f a(z)d w
D
__ Z.
E
z, B
The analytical evaluation of the integral yields:
3-28
Chapter 3 Elastic analysis of a single-pile under axio.lloading
(3.12)
kzl
rigid base
Figure 3.15: Base correcting layer
Now if we equate w
D
from equations (3.10)and (3.12J the modulus for the base correcting
layer can be found from this expression:
Again by choosing a reasonable value of k the values for the modulus can be obtained, from
the above expression, for the base correcting layer. For the three finite element meshes used
in the analyses described in Section 3.5.3, the typical values ofE, E
B
, VB for the smallest
3-29
mesh are 745000 kPa, 352400 kPa, 0.49, respectively, with a load spread angle of 8 = 45",
I
k = 2, _ = 40 and r
l
= lOro'
r!)
3.5.3 Finite element analysis results
Figure 3.16 shows an axi-symmetric finite element mesh. The pile, the soil and the correcting
layers are modelled with 182 eight-noded quadrilateral elements (CAX8) included in the
element library of ABAQUS (version 4.6 from Hibbitt, Karlsson and Sorensen, Inc.,
Providence, Rhode Island), the commercial finite element package used for this analysis. The
rigidity index G/s
u
is equal to 5*10
3
and a Poisson's ratio of 0.49 is used for the soil i.e.
an undrained or total stress analysis is performed. The mesh is quite coarse which will affect
the accuracy of the finite element results. As we said at the beginning, our main emphasis
here is to select a suitable method for the analysis of axially loaded piles - detailed discussion
of mesh refinements is omitted for that reason. The present analysis is not intended to inspect
the possibility of slippage at the interface therefore pile-soil boundary has not been modelled
using zero thickness interface elements (see e. g. Cheung et al. (1991. Although, interface
elements are helpful in obtaining stresses at the pile-soil boundary, a simple way of extracting
shear stresses at the pile-soil boundary described by Burd (1986) has been adopted. In this
method only the force at the midside node is used to determine the shear stress at the element
boundary. The total shearing force acting on the element side is3/2 of the midside nodal
force of the element and calculation of shear stress on the element boundary becomes quite
straightforward
In the present axi-symmetric finite element analysis, the shear stresses and the vertical
displacements down the pile are compared with those obtained by the integral equation
method. The results correspond reasonably well and are shown in Figures 3.17 and 3.18.
Various mesh sizes are shown in brackets at the bottom of these figures. The shear stresses
3-30
Chapter 3 Elastic analysis of a single-pUe under axial loading
1
l
or !eel fng ayer
..:::..;;...
Figure 3.16: Finite element mesh with two correcting layers
obtained using these methods agree quite well with each other except near the top and bottom
of the pile. This suggests that a more refined discretisation of soil-pile boundary in those areas
is required in both the integral equation method and the finite element method. The trend in
question was also noted by Evangelista and Viggiani (1976) who reported that the accuracy
of the results from the integral equation method improves by varying the size of the pile-soil
boundary discretisation down the pile. For the displacements, on the other hand, no similar
close agreement is observed. The discrepancy within the displacements arises from the coarse
discretisation of the meshes, the inherent approximations involved in these methods and the
approximate analysis used for obtaining the properties of the correcting layers.
Figure 3.18 shows that, for larger mesh dimensions, the agreement between the results of the
integral equation method and the finite element method is much closer than for the smaller
mesh dimensions. As we can see from Sections 3.5.1 and 3.5.2 that the properties of the base
correcting layer depend on the choice of load spread angle. The gap between the results from
the integral equation method and the finite element analysis using correcting layers will
3-31
Chapter 3
0
-1
'"
'"
~
~
-2
..:
'"
1
-3
-I ~ "0"-1:
~ 0
I-' .til
c:
~
-4
:a
,
c:
0
c:
-5
-6
0

IEM +
0.2 0.4 ~ 0.6
non-dimensidnalised depth
FEM(20r) <> FEM(30r)
I
=20
d
0.8
lJ. FEM(50r)
Figure 3.17: Shear stress distribution along the length of a flexible pile
-0.01
-0.02
-0.03
~
-0.04
'-'
~
.,
-0.05
E
.,
-0.06
loA
].
I
~ -0.07
d
'" 20
-;
(J
'e -0.08
.,
v
0.49
>
-0.09
Ep
-= 50
Es
-0.1
-0.11
o 0.2 OA L 0.6 0.8
I
non-di mensionalised depth
FEM(2Or) .. IEM () FEM(3Or) lJ. FEM(5Or)
Figure 3.18: Vertical displacement profile of a flexible pile
diminish, giving an excellent match between the two curves if a load spread angle of 47
3-32
degrees is chosen. However, some discrepancy between these two sets of results will remain.
This is not unexpected as the integral equation method and the development of correcting
layers are based on different underlying assumptions. With the revised properties of the base
correcting layer i.e. using a load spread angle of 4 7 ~ even the smallest mesh is sufficient to
obtain results of reasonable accuracy - the difference between the results from these methods
(shown in Figures 3.17 and 3.18) is almost the same as we can see in Figures 3.12 and 3.13.
This proves that the correcting layer approach significantly reduces the problem size in a
finite element analysis. In the particular problem here described, the mesh dimensions are less
than half of those previously adopted (e.g., Randolph (1977); Baguelin and Frank (1980)). A
trend opposite to that shown in Figure 3.17 can be observed in Figure 3.18. This is due in
part to the reasons given above, and in part to improper modelling of the pile-soil boundary.
3.6 Conclusions
The validity of the integral equation method (or the quasi boundary element method) for the
analysis of pile foundations is conditional on the assumption that the piles are embedded in
a linear elastic soil medium, and spaced widely enough for each pile to be outside the radius
of influence of its neighbour. The pile-soil interactions are assumed to be elastic and
considered only at the pile cap level. However, in closely spaced pile groups or even for a
single-pile, the condition of elastic interactions no longer apply, and therefore linear elasticity
is not suitable to model the soil behaviour. It follows that the semi-analytical approach given
by Randolph and Wroth (1978), (1979a) for the behaviour of a single-pile and pile groups is
no longer satisfactory. A more refined constitutive soil model is required for the real life
situation. One possible extension of the boundary element method to incorporate plasticity is
tested in this chapter. However, as we have seen, such an extension makes the method
inefficient: it loses the advantage of requiring only the boundary discretisation of the domain
(whereas FEM always requires discretisation of both surface and interior of the domain), and
it therefore loses the merit of a smaller system of algebraic equations than is needed in the
3-33
Chapter 3 Elastic analysis of a single-pile under axia//oading
finite element method. The algorithm to incorporate plasticity and the overall solution
procedure becomes similar to that of the finite element method which is already well
established. Moreover, the architecture of a boundary element computer program is not
conceptually as easy as the one for the finite element method. In addition, strong numerical
singularities arise in the formulation of the boundary element method which restrict its further
development. For all these reasons, the choice has to be the finite element method. The
correcting layer approach for the modelling of infinite boundaries in the finite element mesh
for the analysis of a single vertical pile under axial loading is quite successful. It not only
models the infinite boundaries quite well but also helps to reduce the overall size of the
problem significantly (for comparison see Figures 3.2 and 3.16). The analysis of pile groups
can therefore be computationally much less expensive.
Finally, it is worth emphasising that, in an offshore environment, piles are subjected to not
only vertical but also horizontal and moment loading. Proper interface analysis for vertical
loading was deliberately omitted in this chapter as irrelevant for the simplified problem
considered in the axi-symmetric two-dimensional finite element model.
3-34
Chapter 3 Elastic analysis 0/ a single-pile under axial loading
APPENDIX 3A
Determination of influence coefficients
free surface
element 1
"
.-..
-
jJ
,
c
v
h
e
z
~
t-
-'"
1
I
l-cj - ~
II df D
-'-- --'-
I
-
-'--
j
,
element n
-'-
Figure 3Al: Geometry of a typical pile
The geometry of a typical pile element is shown in Figure 3A-1. For the point i at the mid
height of the i th element on the periphery of the pile of radius r 0' the influence factor for
vertical displacement due to the uniform shear traction on the j th element is:
jl, 1t!2
Ii} = 4r 0 f JIm de de
(i-I)I, 0
3A-l
(3A.l)
where
it == length of the pile element, and
1 influence factor for vertical displacement due to a vertical load, which is taken from
m
Mindlin's equation, as follows:
1
1m = -----
16rrG(1-v)
Rl
(3 -4v)
+---
Rl
where
3 -4v)z 2-2cz+2c 2)
+ - - - - ~ - - - -
R3
2
z = h+c (h and c are defined in Figure 3A-1.)
z = h-c
I
(3A.2)
The integral with respect to c in equation (3A.1) is evaluated analytically and the integral
with respect to 8 is evaluated numerically using Simpson's rule.
The geometry of the pile base is shown in Figure 3A-2. For the element i on the shaft, the
influence factor for vertical displacement due to uniform loading at the base is given as:
Zit r.
fb = J J I
m
rdrd8
o 0
3A-2
(3A.3)
Chapter 3
... f .
2
1
Figure 3A-2: Geometry for integration over circular area
where I In is given by equation (3A.2), but now
c = I
Again the integration with respect to r is carried out analytically and integration with respect
to e is evaluated by Simpson's rule.
The integral for the flexibility coefficient for the displacement at the centre of the base
resulting from shear traction on a shaft element j is given by:
3A-3
Chapler 3 Elastic analysis of a single-pile under axial loading
Jf
f ~ = 2nro f 1m de.
V-I )1,
This integral can be evaluated analytically.
For the vertical displacement of the base associated with the stresses on the base itself, the
integral for the flexibilty coefficient for the base is:
where
c = I
R=r
I
The analytical evaluation gives:
where
It
4
3A-4
Chapter 4 Interface Analysis
CHAPTER 4
INTERF ACE ANALYSIS
4.1 Introduction
In order to proceed with the detailed analysis of interactions between the pile and the soil,
it is necessary to consider the interface between them. The main reason for carrying out an
interface analysis is to model the possibility of loss of contact between them near the ground
surface. Soil close to the pile under such a situation is highly constrained in front of the pile
and stretched at the back of the pile where it has the tendency to break away from the pile.
An additional advantage of having an interface element is to extract the stresses at the pile-
soil interface, which are the key to the understanding how the problem is behaving. In
Chapter 3, the objective was to select a suitable numerical method for the analysis of a single
vertical pile under axial loading considering the simplest possible case. No attention was
given to performing a complex interface analysis to extract the stresses at the pile-soil
interface.
In this chapter, a purely cohesive interface element to model small relative slip and loss of
contact at the back of the pile is developed and validated by verifying various test problems.
It is also quite useful in extracting stresses at the interface. It is important to have suitable test
problems where we have loss of contact between bodies of two dissimilar materials. One such
problem considered in this chapter is a plane strain surface footing subjected to combined
vertical and moment loading. The footing experiences an uplift force which causes a partial
loss of footing contact with the ground surface. An approximate solution has been given by
4-1
Chapter 4 Inletface Analysis
Hansen (1970) which is modified to take account of the uplift force acting on the footing
which arises due to the applied eccentric loading. The interface element has been successful
in modelling the uplift force and the results match with the modified Hansen's (1970)
solution. In the next chapter, the response of a single-pile under lateral loading in a purely
elastic, cohesive material will be studied.
4.2 Background
In a finite element analysis, the interface between bodies of two dissimilar materials needs
special attention. If, a perfect bonding is maintained between the bodies, an interface analysis
is useful to extract stresses at the interface. But if, at some point in the loading history, the
bond between the two contacting bodies breaks down then an interface analysis is necessary
to model either slip between the two bonding surfaces, such as the rigid pile and the soil
continuum shown in Figure 3.16 or in addition to slip, loss of contact as in the case of a
laterally loaded pile. Conventional continuum elements have often proved inadequate in
modelling relative slip and loss of contact between two bodies of dissimilar materials (see e.g.
Gens et al. (1989) and Bell (1991. Interface elements have been applied comparatively
successfully. They have been applied, typically, to the problem of modelling jointed rock
masses (Goodman et al. (1968, or to the simulation of particular soil-structure interface
behaviour in, e.g., piles (Zaman (1985, footings (Bell (1991 and two-layer soil systems
separated by soil reinforcement (Brocklehurst (1993).
Basically there are two types of analysis procedure to deal with the two bodies in contact with
each other. One works on the constraint approach in which the two contacting bodies are
treated independently and allowed to meet at the mutual boundary satisfying some contact
constraints. It is quite a powerful method which does not allow interpenetration of the two
contacting bodies. But it has the disadvantage that it does not automatically provide interface
stresses i.e. there is no obvious way of storing the stresses at the interface. This approach
4-2
follows an iterative solution procedure e.g. as described by Katona (1983), who models the
friction interface between the two contacting bodies as sequence of node pairs. The different
states of the interface i.e. stick, slip and debonding are imposed through a constraint matrix.
(Earlier, a similar model for the interface behaviour 'linkage elements' was presented by
Frank et al. (1982) in which only the connections between opposite nodes are considered.)
Moore and Booker (1989) and Kodikara and Moore (1991) adopted a similar approach to
model circular tube-soil interface. They suggest using a non-linear frictional-adhesive zero-
tension interface model at the boundary between the two separate bodies, rather than
employing a discrete interface element. The technique they presented is capable of capturing
large rotations and slipping as well as the occurrence or loss of contact at the interface.
(However, large relative slip is not an imponant characteristic of the problem under
consideration here, therefore, this model will not be pursued in this chapter.) The second type
of procedure is based upon the stiffness approach, as described by Goodman et al. (1968),
which falls more naturally into the finite element formulation. In Goodman et al. (1968), the
interface element is said to have zero thickness, as also in the work of Gens et al. (1989) and
Burd and Brocklehurst (1991). The relative displacements between opposite nodes become
primary deformation variables instead of movements at the nodal points, and observe elastic
or elasto-plastic constitutive behaviour. The interface element can also have a finite thickness,
as proposed by Desai et al. (1984), who called it a 'thin layer element'. The element works
more like a continuum element of small finite thickness. It represents a thin layer of material
between two bodies and seems suitable for jointed rock masses. The element exhibits a non-
linear elastic plastic behaviour. The values for the normal stiffness and the thickness of the
interface are suggested, following parametric studies. The shear stiffness is determined from
a laboratory direct shear test but no suggestions have been made for choosing cross coupling
terms. The idea of a 'thin layer element' is derived from the early work presented by
Zienkiewicz et al. (1970). They introduced quasi-continuum elements possessing a weakness
plane in the direction of the interface plane. However, it has been reported by Ghaboussi et
al. (1973) that numerical difficulties may arise from ill-conditioning of the stiffness matrix
4-3
due to very large off-diagonal terms or very small diagonal terms which are generated by
these elements in certain cases. Within the stiffness approach, the most obvious way of
modelling interfaces is to use standard continuum finite elements of small thickness. Griffiths
(1985) demonstrated the suitability of a thin eight-noded quadratic element with an aspect
ratio of up to 100, by analysing interface behaviour in various geotechnical problems.
Griffiths reported that as the aspect ratio is increased to 1000, a much tighter tolerance for
convergence is required than would be required in normal usage of elements other than
interface elements. Finally, Herrmann (1978) presented an analysis procedure which combines
the constraint and the stiffness approach.
For the problem being tackled here, the stiffness approach is obviously preferable to the
constraint approach and to the combination of them proposed by Herrmann (1978). In the
analysis which follows, a purely cohesive interface element is described. It is adapted from
Goodman et al. (1968) for undrained clays, just as Burd and Brocklehurst (1991) adapted this
approach for sands to study the interface stresses in a two-layer soil system. It will be used
to study the pile-soil interface behaviour in clays under lateral loading in the next chapter.
4.3 Six-noded isoparametric interface element formulation
In this section, a zero thickness, six-noded isoparametric interface element which is
kinematically compatible with a six-noded linear strain triangular element (used for the plane
strain finite element analyses in Chapter 5) is described. The geometry of the parent element
and the mapped element is shown as Figure 4A-l in Appendix 4A. The element geometry
consists of three pairs of nodes sharing the same coordinates. However, within each pair, the
nodes have independent degrees-of-freedom. The corner nodes are numbered 1-4, 2-5 and the
middle nodes carry node numbers 3-6 such that the difference between the bottom node
number and the corresponding top node number is three. The bottom nodes of the interface
elements are connected to the corresponding continuum nodes in the finite element meshes
4-4
shown in this chapter. As all the finite element analyses are displacement controlled, the
imposed displacements are applied to the top interface nodes.
4.3.1 Kinematic equations
Although details about the kinematics of the element can be found elsewhere (Burd and
Brocklehurst (1991, the basic equations are repeated here for convenience. Figure 4.1(a)
shows two coincident arbitrary points, a and b, lying on the element sides A and B. The
vector of the global co-ordinates of the two points in the element is x
T
= [ xa Y
a
Xb Y
b
], and
is related to the vector of global nodal co-ordinates X by:
x = [N]X (4.1)
where [N], the shape function matrix, is a function of a. The elements of this matrix are
given in Appendix 4A. (Note that the co-ordinates ( xa Y
a
) are identical to ( Xb Y
b
).)
Similarly, the vector of relative velocities [U
r
v,] is given by the matrix equation:
(4.2)
where [B] relates the vector of nodal velocities, U, to the vector of relative velocities at a
point in the interface, u" and is given in Appendix 4A.
4-5
f( a), defined by:
f( a) = 11: I + c = 0
(4.5)
where c is the cohesion at the interface. The vector of stress is denoted by a and ~ I is the
absolute value of shear stress. The dotted line in Figure 4.2 marks the tension cut-off, the
stress at which we get a sudden loss of contact. This line is not part of the yield surface Le.
we can use the equations of plasticity everywhere on the yield surface and not on the dotted
line. The plastic potential in this analysis is assumed to be identical to the yield function (i.e.
plastic flow is fully associated):
g(a) = It I +c (4.6)
where c is the cohesion of the interface.
An associated flow rule requires that the plastic potential is the same as the yield function
(Le. g(a) = f(a). Therefore, the plastic relative velocities will lie normal to the yield surface.
The vector of stress, a, of an arbitrary point in the interface can be differentiated with respect
to time (denoted by a superior dot), to give the stress rate vector, OT = [t oJ. This is then
related to the relative velocities, u
r
' by the expression:
4-7
(J
'\
~
yield surface
~
relaxed
-0- (J
elastic
plastic
-1
Figure 4.2: Yield surface with an allowable tension cutoff
or substituting equation (4.2)
cr = [D]EP [8] V (4.7)
where [D]EP is the elasto-plastic material stiffness matrix and is expressed in tenns of the
elastic component of the material stiffness matrix [D]E and the plastic component [DY thus:
(4.8)
Note that the interface element relative velocities are analogous to strain In continuum
mechanics.
The elastic material stiffness matrix is essentially the same as that proposed in the original
work proposed by Goodman et al. (1968), which is given by:
4-8
[D]' = [ ~ ~ J
(4.9)
The parameters ks and k" represent, respectively, shear and nonnal stiffness of the interface.
The plastic material interface stiffness matrix is derived from the approach conventionally
adopted for continuum elements (e.g. Burd (1986)) and is given by:
-[D]' [ # ] [ ~ ]'cDl'
[ : ~ J'rDl' [ ~ ; ]
(4.10)
This equation (4.10) involves the assumptions that the relative velocities can be decomposed
into elastic and plastic parts and that the stress rate is related to elastic relative velocities. The
yield function, f(a), is given in equation (4.5), and the plastic potential g(a) in equation
(4.6) is used to define the plastic relative velocities, using equations of the fonn:
f' = "I ag
u, fI.
at
(4.11)
(4.12)
4-9
where A is an undetermined mUltiplier.
Although it is possible to calculate rDY
P
using the standard approach but following the
experience of Brocklehurst (1993) it is kept equal to [DJE to make the calculations better
conditioned i.e. a localised initial stress approach is adopted for this part of the finite element
calculations.
4.3.3 Finite element equations
If a set of virtual velocities U are applied to the nodes of the element, then by the principle
of virtual work:
U
T
p := L (4.13)
where u: are the relative velocities compatible with the virtual nodal velocities, U' and
P, is the vector of nodal forces which are in equilibrium with the stresses aT = [t a]. The
element surface area, dA, is given by:
dA := J da (4.14)
where J is the Jacobian or transformation matrix (see Appendix 4A).
It follows, therefore, that the virtual work equation (4.13) can be fe-written, substituting the
relative velocities vector as given in equation (4.2), as:
4-10
Chapler 4 Interface Analysis
u
T
P = i UT[Bfa J do..
E _ _
(4.15)
Since this expression must hold for arbitrary nodal velocities U ., it follows that:
p = i[Bfa J do..
_ E _
( 4.16)
Equation (4.16) is differentiated with respect to time to give an expression for the nodal force
rate, giving:
(4.17)
The term in the integrand may be used to derive conventional infinitesimal displacement
stiffness matrix (Zienkiewicz (1977 by substituting equation (4.7) and writing cr in terms
of nodal velocities, which gives:
~ = L[Bf ~ J da = L[Bf [DJ [BJ ~ J da
(4.18)
where [D] will be either [Dt or [DY
P
depending upon the state of stresses (elastic or
plastic) during the analysis.
The nodal force rate (equation (4.17 can therefore be simplified to:
4-11
P = [K] U ( 4.19)
where [K] is the incremental stiffness matrix given by:
[K] = L(rBY [D1 [BJ J) du. (4.20)
4.3.4 Stress update calculations
The finite element calculations are typically conducted as a series of steps for a particular
prescribed displacement, or prescrihed force, or a mixture of the two. At the end of each load
step, corresponding to a finite time increment 6t, the updated normal and shear stresses
associated with the calculated relative velocities are detennined by a 'stress update procedure'.
This gives the out-of-balance nodal forces that are needed for incremental solution procedure
and provides a check on the general equilibrium of the system. The computed final stresses
are characterised into different states according to their position in stress space relative to the
failure envelope shown in Figure 4.3. The three possible states are denoted as Elastic, Plastic
and Relaxed as indicated in Figure 4.3.
The length of the increment of time for each load step, 6t is nonnalised into a time interval
of unity. At the start of a load step ilt=Q at the end 6t=1 In order to calculate Gauss point
stresses at the end of a load step, the response of each Gauss poim is assumed, initially, to
be elastic, regardless of initial state. The relevant final stresses 't and (J', are calculated using
the elastic material stiffness matrix [D]E in equation (4.9), giving:
1:=1: +kuilr
o s r
(4.21)
4-12
yield surface
1:
I
.. '*
..
-
possible stresses ~
'"
.. , "*
impossible stresses I
-0- 0-
elastic iii
'" ~ ..
relaxed
flo'"
plastic
-1:
Figure 4.3: Yield surface and possible cases of constitutive behaviour
(4.22)
where 'to and 0"0 are the initial Gauss point stresses, and !J.t =: 1. This assumption of a purely
elastic stress update is then checked to ensure that the final stresses, 't and 0", do not violate
the failure conditions:
1) If j't/c I > 1 then the updated stress point has yielded and must be re-analysed
plastically beyond the point at which the stress path intersects the yield surface.
2) If 0" ~ 0"1 (0"1 is the allowable tension cut-off) then separation of the interface has
occurred during the load step and the stress point is relaxed by setting the
stresses to zero; i.e. 0" = 1: = ().
Given that the initial stress state is either elastic and lies inside the failure envelope, or plastic
and lies on the failure envelope, or relaxed (see Figure 4.3), there are five possible cases of
4-13
constitutive behaviour:
Case i
The calculated stresses at the end of the assumed elastic step are correct.
Case ii
The Gauss point stress is initially on the failure envelope, but the final stresses are correct
since the behaviour involves a purely elastic increment away from the yield surface.
Case iii
As in Case ii the initial Gauss point stress lies on the failure envelope, but in Case iii the
assumed elastic increment produces stresses olltside of the yield surface. The entire step,
therefore, has to be fe-analysed plastically lIsing the elasto-plastic stiffness matrix [D]EP in
equation (4.8).
a = a + k v !::'t .
o n r
But shear stress will be:
though only when a < at' If a ~ at then 't and a are set equal to zero representing the
separation of two contact bodies. If no separation is allowed, then a = at and 1: = Ie I.
4-14
Case iv
The initial point is elastic but the elastic updated stresses lie outside the yield surface. The
proportion of the time increments, PI' )12' at which the stresses cross the yield surface, must
be determined from the equations (4.21), (4.22) and the yield function equation (4.5):
PI
= (It I - lei)
when a < 0'/' Otherwise, if a ~ 0'1 then:
k u
s r
(a - a)
kv
" r
(4.23)
(4.24)
The elastic and plastic parts of the stress increment can then be calculated separately for the
time increments )1D.t and (l-p)D.t, respectively. Where )1 can be )11 or 112 depending upon
the location of the point of intersection with the yield surface either left or right side of the
origin in Figure 4.3. The yield stresses 't and a , are calculated from the elastic equations:
y y
(4.25)
a = a + k v uD.r
yon r
l
(4.26)
and the final plastic stresses, given by the elasto-plastic material stiffness matrix [D]EP, are:
4-15
Chapter 4
where t cannot be greater than Ie I,
t = t
y
and where 0' has always to be less than 0'/'
Interface Analysis
(4.27)
(4.28)
The proportion of time increment ).12 is only considered when separation at the interface is
not allowed, i.e. when the interface stresses are not set to zero. This arbitrary condition has
been imposed in the test problem of two rigid blocks in Section 4.4 below to verify that the
present formulation of the element works also when only tensile stress is applied.
Case v
The initial point is elastic but the elastic updated stresses lie outside the yield surface to the
right side of the origin in Figure 4.3. The stresses in this case are set to zero to mark the
separation of the interface.
A careful choice of a numerical integration scheme is important to evaluate the incremental
element stiffness matrix of equation (4.20). It has been suggested that a Newton-Cotes
integration scheme provides a better solution then a three point Gaussian quadrature scheme.
Gens et al. (1989) investigated the application of one- and two-dimensional interface elements
to model a pull-out test in which a reinforcement element is extracted from a surrounding
material, they observed spurious oscillations of the shear stress when a Gaussian integration
scheme was used but which did not appear when a Newton-Cotes scheme, where all of the
4-16
integration points lie at the element nodes, was adopted. However, Gens et al. also pointed
out that the adoption of a NewlOn-Cotes integration scheme resulted in decreased accuracy
in integrating polynomial expressions. Langen and Vermeer (1990) adopted this scheme in
the pile penetration analysis in which they introduced internal interface elements near the
corners of the pile base. They reported that although this technique proves to yield smoother
stress results in the pile penetration analysis, bad results are still obtained near corners in the
finite element mesh, as encountered in many problems of soil-structure interaction. Burd and
Brocklehurst (1991) successfully adopted a three point Gaussian integration scheme in their
study of the problem of soil reinforcement: no instability in the interface stresses was
reported. Such stresses do not occur neither in the case of a reinforcement problem nor in the
case of the problem being tackled here. Therefore, following the experience of B urd and
Brocklehurst (1991), a three point Gaussian integration scheme is adopted to evaluate the
incremental stiffness matrix. The numerical instability experienced by Gens et al. (1989)
might be due to the generation of large shear stresses during a pull out test which is quite a
severe test for an interface element.
4.4 Verification of the plane strain interface element
4.4.1 Two jointed solid blocks
The plane strain interface element formulation described in Section 4.3 is tested in this simple
finite element problem which consists of two solid blocks, one sitting on the other and
separated by the new interface element, as shown in Figure 4.4.
The lower block is fixed at the base and the upper block is free. Two loading conditions were
tried: upper block sliding horizontally over the lower block and upper block being pulled
away vertically from the lower block. These two loading situations are depicted in Figure 4.4.
The two blocks were given very high strength compared to the interface strength so that
4-17
Chapter 4
b
h
interface element
I
Figure 4.4: Two jointed solid blocks
Interface Analysis
(zero thickness)
15
1.0
1
yielding could only occur at the interface. Goodman et al. (1968) showed that the joint
strength, c, and the interface stiffness values are independent of each other. Following his
experience, in this problem, therefore, arbitrary values were set for the joint strength, c, the
normal stiffness, k , and the shear stiffness, k , of the interface element. The results here are
n s
obtained through displacement controlled finite element analyses. The load displacement
curves are shown in Figures 4.5 and 4.6.
The figures show that the elastic deformations in the element took place at the prescribed
initial tangent stiffness values until the yielding occurred at the interface and the plastic
deformations continued. The computed shear stresses are equal to the joint strength, C, and
the normal stresses are equal to allowable tension, at' at the interface. This shows that the
computed stresses at the interface are consistent with the applied normal and shearing stresses.
4-18
25
20
5
k. = 50kN/m
k lOOkN/m
s n =20kN/m2
u
-------.
...
Ilum! .. I
".m." ;) / / II 7
o ______ ______
o 0.1 0.2 0.3 0.4 0.5 0.6
displacement (m)
Figure 4.5: Interface stresses against relative displacement (shearing mode)
12
10
8
4
2
o
\ \ .... interface element
7777777777
k =50kN/m
s
k
n
lOOkN/m
(J allo
w
lO
kN/m2
\
0.1 0.2 0.3 0.4 0.5 0.6
displacement (m)
Figure 4.6: Interface stresses against relative displacement (tensile mode)
4.4.2 Plane strain footing subjected to moment and vertical loading
The ultimate vertical load carrying (bearing) capacity of a surface strip footing on a cohesive
soil, in the absence of a moment loading is given by this plasticity solution based on Tresca's
4-19
Chapter 4
yield criterion (Calladine (1985)):
P wi
sB
u
Interface Analysis
1t+2
where P wi is the ultimate bearing capacity, Su is the undrained shear strength of the cohesive
soil and B is width of the footing.
As the footing receives a moment loading in addition to the vertical loading, the ultimate
bearing capacity of the footing is significantly reduced. There are no exact solutions available
for this problem. Various models have been proposed e.g. Meyerhof (1953), Hansen (1970)
and Vesic (1975) based on the concept of effective contact area of the footing, acknowledging
the fact that tensile stresses develop between the footing and the soil and footing experiences
a partial loss of contact when subjected to a moment loading. The effect of moment loading
is considered by assuming the vertical load to be applied at some eccentricity from the
centroid of the footing and the point of application of load coinciding with the centroid of the
effective (reduced) contact area of the footing with the soil. It was thought to be a good test
problem to validate the interface element by comparing the finite element analysis results with
one of these existing analytical solutions. Hansen's (1970) solution, widely used by the
offshore industry, is chosen for this comparison. This solution has been modified to take
account of allowable uplift force acting on a strip footing subjected to moment and vertical
loadings. The modified solution is used to validate the results obtained from the interface
analysis described in Section 4.4.2.2.
4.4.2.1 Analytical solution (loss of footing contact)
This section presents a simple analytical model which is a development of the earlier model
proposed by Hansen (1970) for a footing subjected to moment and vertical loading. We are
modifying Hansen's (1970) equation to model the uplift force acting on a plane strain footing
4-20
Chapter 4
i
i
Bi
I
I
i
i
i
i
p
I_ B _I
-j
0' t I Ir-----..... : --,
i
Figure 4.7: A strip footing resting on the surface of a clay layer
Interface Analysis
subjected to moment and vertical loadings. Consider a strip footing of width B subjected to
eccentric loading, P, as shown in Figure 4.7. The shaded area represents the reduced effective
contact area in compression. A small amount of tension cut-off stress, 0"1' is specified.
The vertical equilibrium of the footing can be expressed by the equation:
P + cr (B - Ii) = (n + 2)s Ii
I u
(4.29)
The moment, M, applied to the footing due to the eccentricity, e, will be:
, (B - B J ' (B - B J
M = (n + 2)su
B
2 + 0"1 B 2
(4.30)
From equation (4.29), the effective footing width under compression B can be written as:
4-21
Chapter 4 Inter/ace Analysis
P + Bcr
B = -:--_-::-:-__ 1_
(n + 2)s .. + crt
(4.31)
substituting the value of B in equation (4.30), results in the following new equation:
2M = (P + PI )(1 - P + P
l
J
B P wi + PI
(4.32)
where:
This is a modified form of the Hansen's equation for eccentrically loaded footings which
represents a parabola in the moment-vertical load space with the origin shifted to the left by
an amount equal to the actual amount of tensile cut-off force (shown in Figure 4.13). If P
l
is zero, equation (4.31) becomes Hansen's equation representing a parabola.
For any given values of tension cut-off crt and su' a combination of the moment and the
vertical load can be determined by equation (4.32) and the width of loss of contact (B - B)
can be determined using equation (4.31).
4-22
4.4.2.2 Finite element analysis
In this section an interface analysis of a smooth strip footing on a clay layer is performed
using the new interface element.
Figure 4.8: The finite element mesh with interface elements
Figure 4.8 shows the mesh with interface elements representing a strip footing which was first
pushed vertically until the collapse load was achieved, it was then rotated to predict the
reduction in vertical bearing capacity of the footing with the application of moment loading.
The displacements are applied at the upper layer of nodes of the interface element. Generally,
mesh dimensions for the prediction of vertical collapse loads accurately, are of the order of
five times the width of the footing or, in the case of a circular footing, ten times its radius
(Sloan (1981); Sloan and Randolph (1982); De Borst and Vermeer (1984. In the present test
problem the mesh dimensions are 30% smaller than usual as the aim is validation of the
interface element rather than accuracy in the predicted collapse loads. The material properties
used in the analysis are: rigidity index CIS .. = 50 and Poisson's ratio v == 0.49.
Figure 4.9 shows the predicted vertical collapse load (without any moment loading) using the
4-23
250
200
178
analytical solution
,...-...
Z
150
..::.:::
---
"0
~
100
0
-
______ :-iwnericaJ solution
50
Analytical solution
0.05 0.1 0.15 0.2 0.25
displacement (m)
Figure 4.9: Load displacement curve with interface elements
mesh with interface elements shown in Figure 4.8. The mesh results in an over-prediction of
the collapse load by just over 10%. The normal and shear stiffnesses of the interface were set
to very high values to make sure that the interface should not yield before the continuum
elements did. The values of the interface stiffness need to be chosen very carefully, taking
into account the continuum material properties. The values of ks and k" cannot be too high
or too low: too high results in ill-conditioning, too low leads to errors (Brocklehurst (1993.
The ratio of shear stiffness, G, to ks has dimensions of length and represents the equivalent
continuum thickness of the interface element when subjected to pure shear. For a rough
footing base, the problem at hand, Glk
s
should ideally be set to a very small dimension to
represent a footing base that has no movement relative to the soil until the adhesive strength
(joint strength) of the interface, c, taken to be the undrained shear strength of the continuum
s , is exceeded. Bell (1991) found that if the values of Glk were too small, it caused
u s
numerical instability when non-linear behaviour is being modelled. To represent a zero
thickness element in the direction normal to the interface, the normal interface stiffness kll
should be given a value such that the ratio Elk" is also a very small value. The ratios Glk
s
4-24
and Elkll in this analysis are chosen to be O.61l
i
where Ii is the length of the interface
element, in the light of the recommendations made by Brocklehurst (1993). The footing is
given rigid body displacements vertically followed by rigid body rotations. As we do the
rotation part we track the values of moment and vertical force. The first set of calculations
involve pushing the footing down until we get complete failure and then rotating. We find
that a failure envelope is traced out and this particular calculation must stop when the moment
reaches the apex of the failure envelope. It is possible to explore other part of the curve by
applying initial vertical loads less than the collapse load. Five load cases were tried and are
shown in Figure 4.10. The results, for zero interface tension, are plotted in Figure 4.11.
Hansen's (1970) solution is also shown for comparison purposes. Some over-prediction is
expected because the discretisation of the mesh but generally the shape agrees quite well.
0.9 ,..-----------------,.--------,
0.8
0.7
0 0.6
~
I-< 0.5
'-'
C':l
(iJ 0.4
..c:
-- 0.3
0.2
0.1
o ~ - ~ - - ~ - ~ - - ~ - - ~ - ~ ~ ~ ~ ~ ~
-0.6 "().5
....... Ioad case 1
_load case 4
-0.4 -0.3 -0.2
OIB
v
........ Ioad case 2
-e-1oad case 5
Figure 4.10: Applied footing displacements
-OJ o 0.1 0.2
.......tr.-load case 3
Figure 4.12 shows interaction curves for moment versus vertical load. Both analytical and
numerical results have been plotted. The figure shows that the footing-soil interface tension
has negligible effect. There is a very small increase in the allowable moment carrying
capacity of the footing.
Figure 4.13 shows that when a small tension is allowed at the footing-soil interface, the origin
4-25
Chapter 4 Interface Ana/)'sis
30 r--------------------------------------------------,
25
-.
9 20
Z
~
'-"
..... 15
c::
S
o 10
S
5
50 100 150 200
vertical load (kN)
__ analytical solution .....- fe sol. (load case 1) -6-- fe soL (load case 2)
_ fe sol. (load case 3) --e- fe sol. (load case 4) -A- fe sol. (load case 5)
Figure 4.11: Moment vertical load interaction curves with zero interface
tension
250
of the parabola is shifted to the left by the actual amount of allowable uplift force specified
in the analysis. The area of the footing under tension, eventually representing loss of contact,
can be detennined analytically for a strip footing (see Section 4.4.2.1). However, Bell (1991)
produced contour plots of nonnal stresses showing areas where loss of contact took place in
circular footings subjected to eccentric loading. The forces and moments are obtained either
directly by the summation of elemental nodal forces or worked out by integrating the Gauss
point stresses over the element.
The load cases shown in Figure 4.13 are the same as in Figure 4.11 except that load case 4
is not shown and in load case 5, an interface tension of lOkPa is assumed. For load cases 1-3,
we use method 1 in which elemental nodal forces are summed up to get external forces and
moments. For load case 5, two separate calculations are shown. These calculations correspond
to method 1 and 2 respectively. In method 2, Gauss point stresses are integrated over the
element to obtain external forces and moment. These two methods are carried out to compare
two different ways of getting external forces and moments and to check on stress updating
4-26
Chapter 4
30 r - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ~
25
~ 20
~
'::: 15
~
o 10
S
5
o ~ ____ ~ ____ ~ ____ ~ ____ ~ ____ _ + ~ - - A - ~ - - - - ~
80 100 120 140 160 180 200 220
vertical load (kN)
__ fe sol. (ten.= 1O.kPa) ___ analytical solution
-.to-fe sol. (ten.=2.5kPa) _fe sol. (ten.=O)
Inter/ace Analysis
Figure 4.12: Moment vertical load interaction curves showing the effect of
interface tension
30
25
-"""'
9 20
Z
..!<l
'-'
..... 15
c::
<1)
E
10
0
S
5
0
-50 0 50 100 150 200 250
vertical load (kN)
__ analytical solution ___ fe sol. (load case 1) -.to- fe sol. (load case 2)
_ fe sol. (load case 3) -e- fe sol. (load case 5: melh. I) -A- fe sol. (load case 5: melh. 2)
Figure 4.13: Moment - vertical load interaction curves (shifted)
procedures. The comparisons do not show any significant difference (0.4% in the calculation
of vertical load and 0.6% for moment calculations) in either set of results i.e. from load case
5. By contrast, Bell (1991) reported a quite significant difference in the computed stresses in
his two-dimensional cohesive interface element. The finite element analysis results in the
4-27
in the figures here show that the Hansen (1970) results are quite conservative, as also found
by Bell (1991). (The three-dimensional finite element analysis of a partially penetrated
spudcan (conical) footing presented in Chapter 8 also confinns that Hansen's (1970) solution
gives conservative estimates of bearing capacity of surface footings subjected to vertical and
lateral loading.) Although the finite element analysis results cannot be truly accurate due to
the close proximity of the zero displacement boundaries in the finite element mesh and to the
coarse discretisation of the soil continuum, the over estimation of results is nevertheless
within acceptable limits. The interface element model has produced satisfactory results in
modelling the rigid footing-soil contact problem. In Chapter 5 it will be used for interface
analysis of a rigid circular disc in cohesive soils.
4-28
APPENDIX 4A
y
5
4
4 6 5
X
1
3
2
ex
( -1,0) (G,O) (I,D)
(a)
(b)
Figure 4A-l: Mapping for the six-noded plane strain interface element
Figure 4A-l shows the isoparametric mapping between the parent and the reference elements
for the six-noded interface element
The vectors x and X, that are related by equation (4.1), are defined as:
X T ::: [x y x y]
a a b b
(4A.1)
(4A.2)
Where Xi ' Y
j
are the global co-ordinates of the i Ih node. The shape function matrix, which
4A-l
is required to keep the co-ordinates (xa ' Y) and (x
b
' Y
b
) congruent, is of the form:
1/2 /3
0

h
12
13 0
0 0 0 0 0
[N] =
o 0 o 0

o h
12 13
0

0
00/
1
/
2
/
3
where:
a(a - 1)
11 = 2
a(a + 1)
13 = ---:2::---
The matrix [8] of equation (4.2) is given by:
1 [-/IX' -h< -1;/ -/2y' - / ~ - / ~ ' Ilx 11y' IzX' 12Y' li/ I ~ ' ]
[8] = J 11y' -/1x f
2
y' -fTt' f ~ ' -Iii -hy' f
1
X ' -lzY' fzX' - f ~ ' fii'
4A-2
(4A.3)
(4AA)
(4A.5)
(4A.6)
(4A.7)
Chapter 4
where
x =
d.x
da
= dy
Y
da
Interface Analysis
4A-3
Chapter 5 Analysis of Single Piles under Lateral Loading
CHAPTER 5
ANALYSIS OF SINGLE PILES UNDER LATERAL LOADING
5.1 Introduction
The detennination of the ultimate lateral soil resistance is essential in the design of laterally
loaded piles. In practice, such design is often done using p-y curve analysis with p-y curves
derived empirically (see above, Section 2.3). An alternative, the use of numerical methods to
evaluate a p-y curve, is explored in this chapter. It contains a description of plane strain finite
element analyses that were carried out to predict the ultimate soil resistance when a circular
pile is loaded laterally in a cohesive soil. It is shown that with the help of these analyses, not
only is the ultimate soil resistance predicted accurately but also a theoretical p-y curve can
be derived which match closely with the p-y curve proposed by Matlock (1970). (Frictional
soils are dealt with in Chapter 7.) In these analyses, the assumption is made that the depth
is sufficiently large for behaviour around the pile to be modelled adequately using a plane
strain finite element fonnulation i.e. when the soil at failure moves round the pile in one
plane (Fleming, Weltman, Randolph and Elson (1985)). Numerical analysis of pile groups
follows in Chapter 6.
The ultimate soil resistance which develops along the laterally loaded pile in a cohesive soil
can be calculated by means of the theory of plasticity, assuming failure of soil takes place
along plane or spiral slip surfaces and according to Broms (l964a) hypothesis, at a depth of
approximately 3.0 times pile diameter below the surface where plane strain conditions hold.
5-1
Randolph and Houlsby (1984) presented an exact solution using the method of characteristics
in which they match both upper and lower bound solutions for a rigid circular pile loaded
laterally in a purely cohesive soiL This analysis, excellent as it is, does not provide any
infonnation at all about the size of the actual displacements: it gives only the ultimate lateral
load (Pu =1 L94S
u
d for a rough pile-soil interface and Pu =9. 14Sud for a smooth interface). But
when designing a pile foundation, it is necessary to have as precise an idea as possible about
how far the foundation may move before the ultimate load point is reached. If such
displacements can be satisfactorily calculated it is possible to design the foundation
accordingly. A plane strain finite element analysis does provide infonnation on both the
ultimate lateral load and the required displacements i.e. a p-y curve. Although, it provides a
theoretical p-y curve which is more reliable than the empirically constructed curve, there is
a problem regarding the size of the mesh required to capture displacements exactly. The
analytical solution for a laterally loaded rigid circular disc in an ideal elastic soil medium
given by Baguelin et al. (1977) (reported in Chapter 2) suggests that for an infinite elastic
medium the displacements are infinite i.e. with the increasing size of the finite element
meshes the displacements will not converge to a single value but will be continuously
increasing with the mesh size. If we find the solution to be correct then the problem is how
to model the infinite domain or far-field behaviour of the soil. The important question is
whether we really need to do it when the main interest lies within the near-field. This problem
is addressed in detail, in the following Section 5.4.
In what follows, two out of the three pile loading phases, namely the pile installation phase
and the phase after the soil recovery period Le. after soil equilibration or consolidation phase
when external load is applied, are investigated in respect of how the pile and the soil interact
and affect the ultimate soil resistance. First: it is assumed that there are no initial stresses
present in the soil; then the effect of initial stresses which are due to the self weight of the
soil and after that the effect of initial stresses caused by pile installation are taken into
account. The ultimate soil resistance in both the cases is compared. Secondly: after the end
5-2
of the consolidation phase, when excess pore pressures developed during installation phase
dissipate and the soil recovers, its strength increases as well as the in situ stresses in the soil;
the effect of both the in situ stresses and the increase in shear strength on the ultimate lateral
soil resistance is studied. In all the calculations, the soil is assumed to be elastic perfectly
cohesive material. However, the possibility of using a non-linear elastic perfectly cohesive soil
model is also explored (see Section 5.2). Finally: pile-soil interaction is studied using the
cohesive interface element developed in Chapter 4. The possibility of loss of contact of soil
at the back of the pile (near the ground surface) is modelled using the interface element. (The
ultimate lateral reaction in frictional soils is discussed in Chapter 7.)
The most common constitutive law for modelling undrained clay is the linear elastic perfectly
plastic Tresca model (Calladine (1985)). However, under axi-symrnetric conditions,
singularities (i.e. the yield function is not differentiable) are encountered at the corners of the
yield function (Sloan (1981 )). Even though corrections can be made for these special cases,
the altemative von Mises plasticity yield criterion (Calladine (1985)), shown in Figure 5.1 is
used to model the cohesive soiL
The figure shows the yield criterion in the 1! plane - a section through the yield locus in the
(q-p) stress space. Where q is the difference of maximum and minimum principal stresses and
p is the volumetric stress. The plastic strain rates are derived from a fully associated flow rule
which gives zero plastic dilation rate and is therefore suitable for the modelling of undrained
clay behaviour as discussed by Wroth and Houlsby (1985). A fully associated flow rule
requires that the plastic potential and the von Mises yield functions are exactly identical when
plotted together using consistent axes, and the plastic strains are then normal to the yield
surface.
In elasto-plastic soil models, it is assumed that yield surfaces exist which bound the region
of stress space which can be reached by a soil without incurring irrecoverable deformations.
5-3
von Mises
1t plane
Figure 5.1: Yield criterion in the 1t plane (after Burd (1986))
Within the yield surface, deformations are entirely recoverable and elastic. It has been
common practice to evaluate soil moduli in the fOffi1 of an equivalent linear elastic response
over a wide range of strains. With the recent advancements in laboratory strain measuring
techniques/methods, it is possible to measure soil stiffness more accurately. Evidence from
laboratory tests demonstrates that the pre-yield soil behaviour is inelastic throughout loading
except at very small strains (Atkinson (1993)). Real soils show more gradual yielding process.
Their response is highly non-linear and influenced by many variables such as mineralogy,
fabric, structure, stress state, loading and drainage conditions. It is revealed that the
characteristic soil stiffness-strain curve resembles S-shape (see Figure 5.2). Three
characteristic regions of soil stiffness can be identified. There is a region where the behaviour
is elastic but this corresponds to wry small strains (usually less than about 0.001 %). There
is a region of large strains, greater than about 1 %, where the behaviour is inelastic and there
are large irrecoverable plastic strains. Between these there is a region where the stiffness
decreases with strain and the behaviour is highly non-linear. Figure 5.2 shows that soil is only
linear and elastic within the range of very small strains corresponding to changes of stress of
a few kPa. At large strains, soil is yielding and high proportions of the total strains is
5-4
Chapter 5
G
elastic
0.001
very smaJl small
Analysis of Single Piles under Lateral Loading
1
Cam clay
(yielding)
In e %
large
Figure 5.2: The 'S-shaped curve' relating stiffness to magnitude of strain
(after Atkinson (1993))
inelastic. The question then arises as to the nature of the straining within the intermediate,
small strain region. One possibility is to regard this as essentially elastic, but non-linear and
to use curve fitting techniques to obtain an empirical expression relating shear modulus to
strain. This is the approach followed by Duncan and Chang (1970) and Jardine et al. (1991);
the method requires complex laboratory tests in which stress paths mimic the in situ paths and
numerical analyses which should start and restart at each change in the direction of a stress
path. An alternative approach is to regard soil behaviour in the small strain region as inelastic
with yielding and hardening with kinematic surfaces inside the yield surface. The simplest
approach seems to be to adapt the Cam clay models by including additional surfaces (e.g.
Mroz et al. (1979); AI-Tabbaa and Wood (1989); Atkinson and Stallebrass (1991). A similar
model is proposed by Houlsby and Chow (1994) and is explored in this chapter in addition
to the linear elastic plastic model. Although Simpson (1993) presented a novel idea of
expressing non-linearity with an analogue of bricks on strings pulled by a man. The man was
taken to represent the point in the strain space of a small element. Each brick represented a
proportion of the element. Movement of a brick represented plastic strain and the elastic strain
is given by the difference between the movement of the man and the sum of the movement
5-5
of the bricks, each weighted by the proportion of the soil it represented. Thus the S-shaped
curve can be modelled in a stepwise fashion.
5.3 Choice of the finite element
The difficulty of using the finite element method to do deformation analysis at constant
volume - collapse load calculations (under vertical loading) incompressible or nearly
incompressible materials) - is widely acknowledged. The difficulty arises from the imposition
of additional kinematic constraints on the nodal velocities by the incompressibility condition.
The condition reduces the number of degrees-of-freedom in the finite element mesh. As a
result, the finite element method produces an over-stiff response and predicts the exact
collapse load with an error margin which may be unacceptable. The same difficulty arises
when the finite element method is used to determine the response of laterally loaded single-
piles in soft clay - because the incompressibility condition again applies, producing the same
over-stiffness and unacceptable error margin.
The incompressibility condition also affects calculation of stress field in a solid body. The
result, when the finite element method is used, is spurious variation of stresses across each
element in the finite element mesh. These inaccuracies in the predicted stresses are confined
to the mean normal stresses, the deviatoric stresses being relatively unaffected (see Naylor
(1974)).
There are several methods of dealing with the incompressibility condition in the finite element
calculations. An excellent review of these methods by Burd and Houlsby (1990) can be
summarised as follows:
1) Quadrilateral element with reduced integration. This is adopted to evaluate the element
stiffness matrices. The eight-noded quadrilateral element with four-Gauss point integration
5-6
rule is nonnally recommended (Kooijman (1989.
2) A variant is quadrilateral element with full integration for deviatoric strain tenns. Reduced
integration scheme is used only for the volumetric strain tenns (only they are affected by the
incompressibility condition).
3) 'Mixed methods'. In this approach, the primary variables in the fonnulation of the finite
element are not all of the same type. Hernnann (1965) described a variational fonnulation in
which mean pressure is taken as an additional variable. This was discussed further by
Nagtegaal et ai. (1974).
4) The hybrid method is an alternative fonn of 'mixed method' in which stresses as well as
displacements are used as primary variables.
It has been shown that the variants of the finite element method based on reduced integration
or mixed method approaches are essentially equivalent (see Zienkiewicz et al. (1983) which
has a good review of the application of mixed methods to incompressible analysis).
A quite different way of approaching the application of the finite element method to problems
involving material incompressibility is described by Sloan and Randolph (1982). They showed
that, as the order of polynomials used to describe the strains within an element is increased,
the rate of increase of incompressibility is less than the increase in degrees-of-freedom
introduced by the additional nodes. If the total number of degrees-of-freedom in a finite
element mesh is held constant, therefore, the total number of constraints may be reduced by
using higher order elements. A satisfactory finite element solution should, therefore, be
possible without having to resort to reduced integration. This hypothesis had been confinned
by Sloan (1981) who performed a series of small strain collapse load calculations using 15-
noded triangular elements.
5-7
A similar study was carried out by Bell (1991) in which he defined a new parameter, namely
'free degree s-of-freedom , , to explore the suitability of three-dimensional finite elements for
collapse load calculations.
According to Sloan and Randolph (1982) and, later, Bell (1991), a linear strain triangular
element is quite suitable for performing collapse load calculations. The appropriateness of the
element for the prediction of the ultimate lateral resistance is tested in this and the following
chapter. As we shall see in Section 5.4, a very fine discretisation within one pile diameter
around the pile is required. To predict the lateral resistance accurately, the minimum number
of free degrees-of-freedom are found to be 1000 (see Table 5.1). If this number is held
constant, fewer higher order elements are required i.e. mesh can be less fine.
The reduced-integration approach has been tried elsewhere (see Kooijman (1989), and the
mixed variable approach, as Zienkiewicz et aZ. (1983) showed, is equivalent in practice. But
there are other reasons than novelty in favour of the choice of a triangular element:
1) The mathematical functions defining the displacements within the element are complete
polynomials unlike other finite elements.
2) Aspect ratio of the element is not as important as in quadrilateral elements where high
aspect ratios cause problems in the numerical solution.
3) The reduced-integration approach is of course inadequate in large displacement problems
where geometric non-linear effects have to be considered (Burd (1986), Brocklehurst (1993)).
While the problem we are considering is not a large displacement problem, as soon as a gap
opens up at the back of the pile between it and the soil, large displacements occur but without
causing large strains i.e. strains remain infinitesimal.
5-8
5.4 Finite element analysis of a single-pile
A rigid disc approximation of a single-pile is adopted in the finite element analyses for the
prediction of limit loads for elastic, perfectly plastic soil behaviour. The finite element meshes
used in the analyses are shown in Figures 5.3-8.
~ __________ ~ ~ ~ ~ ' - ~ L L ~ I 5
- - ~ ~ ~ P
Figure 5.3: Coarse finite element mesh (CM!) for a rigid circular disc/ 77
Figure 5.3 and Figure 5.4 show coarse finite element meshes whereas Figures 5.5 to 5.8 show
some refined meshes. Since the problem is symmetric, it is only necessary to analyse one half
of the disc and indeed the consideration of anti-symmetry allows the analysis to be reduced
to a quarter of the disc.
Six-noded linear strain triangular elements are used. (The suitability of these elements for the
prediction of collapse loads for incompressible materials was demonstrated by Sloan et al.
(1982) and later by Bell (1991).) When performing limit analysis or collapse load
calculations, the incompressibility condition imposes additional kinematic constraints on the
nodal velocities which reduce the number of degrees-of-freedom in the finite element mesh.
Choosing the mesh dimensions is an important consideration in a finite element analysis,
especially when modelling an incompressible material behaviour. The close proximity of the
5-9
Chapter 5 Analysis oj Single Piles under Lateral Loading
L<;/
Figure 5.4: Coarse finite element mesh (CM2) for a rigid circular disc
/P7/ 0 ... /P7/
Figure 5.5: Refined finite element mesh (RM1) for a rigid disc
mesh boundaries results in fewer degrees-of-freedom compared to total number of constraints
in a finite element mesh. This severely hampers the quality of the finite element solution.
Error in the finite element solution arises not only because of the closeness of the mesh
5-10
boundaries to the loaded area but also because of incompressibility constraints (discussed
earlier, Section 5.3), and because of coarse discretisation of the mesh. In such circumstances
the finite element solution often overestimates the exact collapse load by an unacceptably
large amount (Nagtegaal et al. (1974); Sloan (1981); Sloan et ai. (1982)). (For a detailed
review of the methods of dealing with the condition of incompressibility in finite element
analyses for the prediction of collapse loads, see Burd and Houlsby (1990).)
Figure 5.6: Refined finite element mesh (RM2) for a rigid circular disc
The most obvious way to get degrees-of-freedom in a finite element mesh sufficient to predict
collapse loads with reasonable accuracy, is to set the boundaries far apart from the loaded
area and thus increasing the number of elements.
In order to find the optimum distance to set the outer boundaries, the finite element analysis
needs to be repeated by varying the distance of the outer fixed boundary each time.
Otherwise, the finite element mesh becomes very large for the computer to process and a
large computer memory space is required, even for a simple problem. The optimum
dimension of the outer boundary is decided by setting an acceptable percentage error in the
finite element solution when compared to the closed form solution (if one is available), and
5-11
e element mesh (RM3) for a rigid
Figure 5.8: Refined finite element mesh (RM4) for a rigid circular disc
bearing in mind that the available closed form solutions (lower bound) are for infinite soil
domain. Based on such numerical experiments, Yegian et af. (1973) suggested that the radius
of the zone of the pile influence needs to be approximately eight times the pile diameter. The
5-12
Chapter 5 A nalysis of Single Piles under Lateral Loading
zone of influence depends very much on the type of calculation. For an elastic calculation we
need quite a large area however, there is a question mark how large? (The elasticity solution
given by Baguelin et al. (1977) (reported in Chapter 2) suggests that displacements increase
monotonically with the increase in zone of influence.) and for a plasticity calculation we need
much smaller area (Randolph and Houlsby (1984. Another way of avoiding these boundary
effects is to couple the finite element mesh to external boundary elements. The main
drawback of this approach however, is that, the boundary element matrices are fully populated
representing coupling of all the nodes on the boundary as a result the coefficients of the
stiffness matrix are altered: the whole matrix is disturbed and its banded nature is lost which
loss reduces the efficiency of the analysis. Neither this possibility nor the possibility of
combining infinite elements around the finite element mesh (e.g. Chen and Poulos (1993))
was explored for the current problem. As mentioned earlier, we have to answer the question -
do we really need to model infinite boundaries for this problem? The problem under
consideration is not a conventional finite element problem in which it is necessary to discount
the boundary effects from the calculated solution. This problem is peculiar in the sense that
the elastic solution given by Baguelin et al. (1977) suggests that before focusing our attention
on modelling infinite boundaries we need to think about the choice of material constitutive
behaviour. In this situation, initially it was thought that in order to predict reasonably accurate
limiting loads and to prevent the finite elements from being locked-up we need to have a
softer layer around the finite element mesh. (Since the use of a softer layer in the analysis of
a large displacement cavity expansion problem has been quite successful (Burd and Houlsby
(1990, there will be an added advantage in setting up correct initial stresses in the finite
element mesh that would be generated due to pile installation.) Therefore, an outer soft
(correcting) layer around the soil was introduced as was done in the analysis of a single-pile
under axial loading (see Section 3.5). For the current problem, it was difficult to calculate the
parameters of the outer soft layer. The parameters were chosen arbitrarily to be half those
used for the inner soil layer. (The analytical solution for a rigid circular disc loaded laterally
in an ideal elastic soil medium given by Baguelin et al. (1977) was considered in order to
5-13
Chapfer 5 A nalysis of Single Piles under Lateral Loading
detennine the soil parameters for the outer soft layer in a manner similar to the one used by
B urd and Houlsby (1990) to determine the soil parameters of the soft layer in the analysis of
a large displacement cavity expansion problem. It was found that the elastic solution given
by Baguelin et al. (1977) is unreliable for the reason that, it gives infinite displacements for
a rigid disc loaded laterally in an infinite medium, which is physically not right. In order to
check the reliability of the solution, a series of plane strain finite element analysis for purely
elastic material constitutive behaviour with increasing mesh sizes (with position of outer
boundary varying from 20r 0 -120r 0' where r 0 is the radius of the circular disc.) were carried
out. The results are shown below in Figure 5.9 which seems to confirm the solution given by
1.2 ,----------------------------.
c.o
0.8
*
linear elastic
0
-
'"'0
* 0.6
'"'0
C':$
..s
0.4
non-linear clastic
0.2
0
0 om 0.02 0.03 0.04 0.05 0.06
d/2R
Figure 5.9: Effect of mesh boundaries on the initial elastic stiffness
Baguelin et al. (1977). The solution appears to suggest that we can choose an infinite number
of outer boundaries where displacements are zero and there is no unique zero-displacement
boundary. The choice of each boundary will give a different value of initial stiffness i.e.
variable initial response of the problem. The initial stiffness will decay gradually to zero as
the boundary will move outward to infinity. The introduction of an outer soft layer, therefore,
can in no way discount the boundary effects leaving the initial response less certain. Its ability
to prevent locking up of the finite element meshes and help in accurately predicting the
ultimate load in an incompressible material for the particular problem under consideration was
5-14
also investigated. Figure 5.10 below shows two curves one with an outer soft layer and other
1.4 ,...-----------------------------,
"0
* ::::l 1.2
U')
* -..:t
0;
--'
,......
::0 0.8
t';\
o
ii 0.6
...
9 0
U .4
~
"0 0.2
t';\
.9
0.1
old
__ without correcting layer
0.2 0.3
with correcting layer
--
Figure 5.10: Effect of correcting layer on the pre-yield response
OA
without a soft layer. It was found that the outer soft layer did not help in achieving a better
accuracy in the predicted ultimate soil reaction but it had a detrimental effect on the elastic
stiffness (i.e. stiffness was not the same as given by the elastic analytical solution given by
Baguelin et al. (1977)). (We shall see later in this chapter that the density of elements around
the disc and mesh grading are rather more important for the accurate prediction of the
ultimate load.) This simple investigation suggested perhaps it is not possible to decide on the
limit of the mesh size using elastic perfectly plastic model and it was decided to repeat the
above analyses with varying mesh sizes using non-linear elastic plastic model developed by
Houlsby and Chow (1994). These results are also shown in Figure 5.9 and we can see this
model is also unable to help solving this particular problem. Although the material model is
non-linear elastic, the problem arises due to the tiny part of elasticity which it has got right
at the start of the analysis. Figure 5.9 also shows that it is certainly better to use a non-linear
elastic model with an initial high stiffness and then decaying stiffness with increased strain
levels but there is no analytical solution available to-date to validate the finite element results.
It was therefore decided to keep the soil model simple (as suggested by Wroth (1984)) and
use elasto-plastic model in the following analyses in such a way that it gives us realistic
solutions. Mesh boundaries can be chosen such that the initial stiffness is representative of
5-15
a particular soil deposit and at the same time match with Baguelin et al. (1977) solution.
Perhaps, it is not always possible to apply complex elasto-plastic soil models with features
like kinematic hardening to every geotechnical problem. We shall see in the later Section 5.8
on lateral loading of rigid disc at the end of the consolidation phase that our choice of mesh
boundaries and soil model turned out to be quite successful. The p-y curve obtained from this
analysis (a loading case close to real field conditions) shows a surprising close match with
the p-y curve suggested by Matlock (1970).
Although the present analysis is of a different nature to vertical collapse load calculations:
collapse load calculations for lateral loading. The dimensions of the mesh on either side are
taken as ten times the radius of the disc, following Sloan (1981) and B urd (1986) for two-
dimensional finite element meshes used in collapse load calculations for vertical loading. A
correcting layer of the same dimensions was then added round the mesh, thus setting the fixed
boundaries at a distance of twenty times the radius from the centre of the disc. Later, we
found that there is no advantage in having a correcting layer round the finite element mesh
even though keeping the outer fixed boundaries at the same distance (20r
o
on either side)
from the centre of the disc. A complete rough contact between pile and the soil was assumed,
with full interface tension behind the pile (no loss of contact of soil with the pile). The soil
(soft clay) had a rigidity index, I, ( where IT = G/S
u
) of 50, Su = 30kPa and Poisson's ratio
of 0.499, representing a nearly incompressible material. A series of displacement controlled
finite element analyses for the different meshes shown earlier (Figures 5.3-8), were done
using the finite element program 'OXFEM', developed by the Civil Engineering Group,
Department of Engineering Science, the University of Oxford. A full description of the
program can be found in Burd (1986), Yu (1990) and Chow (1994). A displacement finite
element approach was adopted. Although, OXFEM was specifically designed to conduct large
displacement, large strain calculations for both the plane strain and axi-symmetric analyses
of geotechnical problems, it is an added feature of this program that it can also perform small
5-16
1.5 ..... ---------------------------,
0.1 02 0.3 0.4 0.5
old
__ C ~ l 1 __ CM2 -4- RM 1 _ RMl -e- RM3 -A- RM4
Figure 5.11: Load factor displacement curves
displacement, small strain analyses. The results are shown in Figure 5.11. The predicted
limiting load is non-dimensionalised by the analytical solution given by Randolph and
Houlsby (1984) and is termed as load factor. Each curve in the figure shows an over-
estimation in the predicted limiting load. The percent over-estimation is given below in Table
5.1.
Table 5.1: Summary of finite element analyses
Mesh No. Number of Total number Geometric Inrornpressibi- Total number Poisson's ratio %.ge over-
Elements of degrees of constrainLtr; lily constraints of f ''''' degrees estimation
freedom of freedom
eMl 174 762 98 522 142 .49 22.56
CM2 348 1494 ]50 1044 300 .49 ]7.5
R:-'lI 432 ]&42 178 1296 368 .49 12.23
RM2 860 3590 230 2580 780 .49 9.03
RM3 696 2906 186 2088 632 .49 6.13
RM4 1376 5642 202 4128 1312 .499 4.5
The solution obtained by using mesh RM3 is superior than mesh RM2 despite the fact that
5-17
Chapter 5 Analysis of Single Piles under lAteral Loading
it has less number of free degrees of freedom than mesh RM2. This superiority is due to the
mesh density as greater number of elements are packed near the pile. A further (40%)
increase in the mesh dimensions of RM4 did not show any improvement in the results. A
small change in the value of Poisson's ratio v from .49 to .499 greatly improved the
numerical solution. The result is shown as the curve corresponding to mesh RM4 in Figure
5.11 (also see Table 5.1).
~ I G - ..
: 800.0000
I 600.0000
I 400.0000
I 200.0000
0.0000
= -200.0000
I -400.0000
I -600.0000
I -800.0000
Figure 5.12: Horizontal stress distribution around a laterally loaded pile
Figure 5.12 shows the horizontal stress (kNlm 2) distribution in x-direction around a laterally
loaded pile. The horizontal stresses are high close to the pile and decay very rapidly within
a radial distance of one pile diameter from the face of the pile. Figure 5.13 shows the soil
yielded zone around a laterally loaded pile. (The figure shows various zones between 0-1
representing elastic and yielded zones respectively.) The extent of this yielded zone do match
very well with the characteristics mesh adopted by Randolph and Houlsby (1984) in the
determination of upper bound solution for a laterally loaded pile in a cohesive soil using
plasticity theory. The prediction of small scattered yielded pockets in addition to the main
yielded zone around the pile is a result of coarse finite element discretisation of the problem
5-18
I
I
I
I
Figure 5.13: Yielded zone of soil around a laterally loaded pile
and inability to model accurately a completely incompressible material behaviour in a finite
element analysis. A vector plot of displacements around a laterally loaded pile is shown in
Figure 5.14.
5.5 Effect of mesh density on solution quality
It has been observed that the problem of a rigid disc moving laterally in a cohesive material
is highly constrained near the disc. The quality of the results is dependent on the mesh
grading and element density within the area adjacent to the disc. A significant over-estimation
of the analytical limit load (Randolph and Houlsby (1984) was obtained when the meshes
were coarse. Subsequent analyses, using refined finite element meshes, revealed that a very
high density of elements is required within a distance of one diameter on either side of the
disc in order to obtain a good solution with an acceptable level of error. In the mesh
refinement process, each mesh fulfilled the criteria for the two-dimensional collapse load
5-19
... ...
-

-
Figure 5.14: Vector plot of displacements (at the end of the analysis) around a
laterally loaded pile (RM4)
calculations put forward by Sloan et al. (1982) and later modified by Bell (1991). namely that
the number of degrees-of-freedom should be greater than the number of constraints in a finite
element mesh consisting of a particular type of element. During the mesh refinement process,
it was found that a stage is reached when further mesh refinements do not improve the
solution. This is due to the fact that the shape functions used in the finite element do not
capture accurately the strain fields around the disc (Burd and Houlsby (1992. (In his thesis,
Sloan (1981) had also made a similar emphasis that, for the accurate prediction of collapse
load in a fairly refined mesh, it is necessary not only that the ratio of the number of degrees-
of-freedom to the number of constraints (including the incompressibility ones) should be
greater than or at least equal to unity, but also other criteria (such as ability to match strain
gradients correctly) should be taken into consideration when choosing a finite element mesh.)
Mesh refinement beyond the stage when any further refinement did not improve the solution
quality. was pointless because this would have made the problem unnecessarily large.
occupying too much computer memory space.
5-20
5.6 Effect of varying soil parameters
The soil parameters involved in the finite element plane strain analyses are undrained shear
strength of the soil Su' shear modulus G and Poisson's ratio v. The analyses carried out by
varying the rigidity index, lr show that the higher the rigidity index, the stiffer the response
and vice versa. But the non-dimensionalised load factor remains unaffected. This means that
if we have soil from two different locations with the same rigidity index but one is stiffer
than the other, the load factor-displacement response curves for the two locations will overlap
on the same plot but the limiting loads p,. will be different. The soil with a higher shear
strength will have a higher limiting load. This result confirms the Randolph and Houlsby
(1984) solution for p u which is independent of the value of G.
'"C
*
1.4
rJ5 1.2
*
'<:I"
~
.....
::0 0.8
~
.
II 0.6
...
.8 0.4
u
~
'"C 0.2
C'::I
.s
__ 0/Su=30
_0/Su=45
0.1
D/d
__ 0/Su=35
--e--- 0/Su=50
0.2 0.3
_____ 0/Su=40
--e,- G/Su=60
Figure 5.15: Load factor displacement curves (variable rigidity index)
0.4
Thus, curves with variable rigidity index shown in Figure 5.15 provide initial non-linear
response of the soil prior to reaching the ultimate load. If we non-dimensionalise the
horizontal scale in Figure 5.15 by lr' the curves in the figure overlap (as shown in Figure
5.16) which is another indication that the ultimate lateral load is independent of shear
5-21
stiffness. A similar observation was made by Chen and Poulos (1993), namely that shear
stiffness of the soil has little effect on the limiting load: in their analyses, there was a
tendency in the stiffer soils to achieve slightly greater value for p u;'
1.4 .-----------_____________ ........
-0
*
::l 1.2
U')
*
~ I
.......
......
:0 0.8
<':l
o
II 0.6
;....
o
] 0.4
-0 0.2
<':l
.9
5
__ G/Su=30
__ G/Su=45
10 15 20
8*G/Su*d
__ G/Su=35 -.to- G/Su=40
-B- G/Su=50 -A-- G/Su=60
Figure 5.16: Non-dimensionalised load factor displacement curves
25
5.7 Analysis of a single-pile assuming initial cavity expansion to model
installation
A complete static analysis of piles comprises three phases. The first one is the pile installation
phase, when the pile is driven into the ground. The process of pile installation is difficult to
consider theoretically as the problem is indetern1inate to an extreme degree. As a preliminary
to this study, arbitrary stresses were generated in the mesh by prescribing radial displacements
prior to lateral ones to the inner boundary nodes of the mesh (representing the circular disc).
The effect of initial stresses in the ground on the ultimate lateral soil reaction is shown in
Figure 5.17. There are three load displacement curves. The curves with initial stresses
(representing different Ko values) have a flexible initial response compared to the third curve
without initial stresses. and they reach the same limiting load as for zero initial stresses but
require higher displacements. The negligible difference between the peaks of these three
5-22
1.4 r---------------------------,
"0
* ;:::I 12
U')
* -::t
0;1
,.......
........
:0 0.
8
cd
o
ii 0.6
~
o
t) 0.4
~
~ 0.2
o
0.1 0.2 0.3
old
____ without initial stresses __ initial stresses due to radial displacement (19% of radius)
-.a.- initial stresses due to radial displacmcnl (33.25% of radius)
Figure 5.17: Effect of initial stresses on the load factor
0.4
curves is due to finite element discretisation and the approximations inherent in the finite
element method. A similar observation was made by Brown and Shie (1991) in their finite
element analyses by varying Ko values in their three dimensional analyses. Although, the
presence of initial stresses do not effect the ultimate lateral resistance but they have an effect
on the initial response of the pile i.e. at loads less than ultimate loads. The stresses in the soil
will depend on the past loading history of the ground. As the pile is driven, these stresses are
altered. Theory derived from an analogy with the monotonic expansion of a cylindrical cavity
(to model stresses in the ground due to pile installation) was proposed by Carter et al. (1979).
More recently two-dimensional strain path solution have been developed by Baligh (1986).
We will not discuss the strain path method in this section and will use the former theory as
its suits computationally to the plane strain finite element analyses described in this Chapter.
(However, some recent studies reported in the literature seems to support the strain path
method e.g. Jardine and Potts (1988) and Bond and Jardine (1991).) Carter et ai. (1979)
showed that these stresses can be taken as stresses caused by expanding a cylindrical cavity
from a zero initial radius. (In contrast to this real situation, axi-symmetric numerical
calculations need to begin with a finite cavity radius to avoid infinite circumferential strain
5-23
that would occur for an initial cavity radius of zero.) Hill (1950) has presented closed fonn
solution for the creation of a cavity from zero initial cavity radius in a rigid plastic material
which follows Tresca's yield criterion (Calladine (1985). Hill (1950) has shown that limit
cavity pressure is reached immediately, before any displacement occurs. Gibson and Anderson
(1961) has presented the closed fonn solution in an elastic-perfectly plastic Tresca material.
Carter et al. (1979) has shown that Gibson and Anderson's (1961) solution approaches Hill's
(1950) solution at larger radii. The authors have found that the doubling of cavity radius is
adequate (to generate stresses in the soil which are created by expanding a cavity from zero
radius) and expanding beyond twice the cavity radius is simply to cause a further growth of
annular region of yielded soil. Carter et al. (1979) have also shown that the final stresses are
independent of the initial overconsolidation ratio OCR. The effect of these initial stresses,
after the expansion of the cavity, on the limiting lateral resistance of the soil subjected to the
applied lateral load is studied in this section. The problem of cylindrical cavity expansion in
incompressible material (material under consideration in this chapter) has also been studied
by Sagaseta (1984) for an elastic-perfectly plastic von Mises material using large strain
analysis and presented close fonn solutions for the stresses. In order to model the effect of
installation stresses on the lateral response of a rigid disc, the stresses (see Figure 5.18) given
by Sagaseta (1984) are set up in the finite element mesh by expanding a cavity by a ratio of
4 before loading the disc laterally. Although, there is a yielded soil around the pile due to
cavity expansion which seems to suggest that there will be zero incremental stiffness upon
loading the pile laterally soon after installation. Figure 5.19 shows a finite value of stiffness
Le. there is still an initial elastic response before reaching the value of ultimate lateral
reaction. This elastic stiffness perhaps arises from the elastic region bounding the plastic zone
and due to change in stress path (i.e. once the cavity is expanded) which soil panicles have
to follow in order to mobilise a collapse mechanism to achieve the ultimate lateral reaction.
The analytical solution of Randolph and Houlsby (1984) for a rigid disc is an exact solution
whereby both lower bound and upper bound solutions match. In the derivation of the
analytical solution, a complete tension at the back of the disc and a complete adhesion all
5-24
7
6
5
4
r - ~ ~ ~ - - ~ ~ - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ~
-I
0, I
__ G
r
0,2 0,3 0.4 0,5 0,6
rolr
0,7 0,8
--+- G
e
Figure 5.18: Distribution of stresses after pile installation
0,9
-.t.-
G
z
1.4 ,-------------___________ ---..
:;p1.2
::l
:irI
a,
,..... 0.8
-
II
l3 0,6
U
~ 0.4
-0
C':I
.9 0,2
0,]
__ zero initial stresses
0.2 0,3 0,4 0,5
old
--+- iniLial stresses due to cavity expansion
Figure 5.19: Effect of installation stresses on the load factor
0.6
I.l
around, are assumed. The presence of initial stresses due to installation, which are assumed
to arise from the expansion of a cylindrical cavity in an elasto-plastic material, is supposed
to relieve some of the tension at the back of the disc. It seems that the stresses arising from
the cylindrical cavity expansion are not sufficient for complete tension relief at the back of
the pile but do result in a partial tension relief, which is compensated by the increased friction
5-25
Chapter 5 Analysis of Single Piles under Lateral LQading
at the sides of the disc. Thus, no effect of the initial stresses on the limiting load is observed
in Figure 5.17. Figure 5.19 shows that when the disc is displaced laterally after expansion of
the cavity, the cavity expansion stresses have negligible effect on the soil limiting resistance.
However, pre-yield behaviour is affected by the presence of initial stress and plastic flow
around the pile. The numerical results of this section confinn the Randolph and Houlsby
(1984) solution which is derived as independent of the current (total) stress level. It is
important to note here that (according to Baguelin et al. (1977) elasticity solution) these
displacements are also affected by the choice of zero-displacement boundary in the finite
element mesh. As shown earlier, the strong dependency of displacements on zero-
displacement boundary can be reduced by adopting a non-linear soil model as described
earlier in Section 5.4 (see Figure 5.9).
5.8 Analysis of a single-pile assuming radial variation of soil strength
around the pile
The subsequent phase following installation is soil consolidation around the pile. The excess
pore pressures generated in the previous phase dissipate over time increasing the effective
stresses around the pile. Along with this rise in stresses, the undrained shear strength of the
soil within the plastic zone also increases over time compared to the rest of the surrounding
soil (see Figure 5.20). This constitutes a radial inhomogeneity of the soil. The shear strength
increase with radial distance from the face of the pile, after the end of the consolidation phase
can be given by a logarithmic distribution, which was presented by Randolph, Carter and
Wroth (1979).
The third and final phase, after consolidation, is the actual loading stage when the soil has
completely recovered from all the initial disturbances and regained strength to carry the
applied loads. This situation is modelled by setting up the stresses and the increase in
strength, after the end of the consolidation phase obtained from Randolph, Carter and Wroth
(1979), in the mesh. In the present study, the logarithmic strength increase was, as suggested
5-26
6
5
4
o ~ ~ . - ~ - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ~
-1
o
0.1 0.2 0.3 0.4 0.5
rolr
0.6 0.7 0.8 0.9 l.l
--O"r -e-O"e
-.01.- 0";,;
_Su
Figure 5.20: Distribution of stresses after the end of consolidation
by Randolph, Carter and Wroth (1979), incorporated to begin with, in the finite element mesh
by dividing it into various sub-regions and specifying the new shear strength values in a
stepwise fashion (obtained from the logarithmic curve at the end of consolidation, as
mentioned above - see Figure 5.22) of these regions. However, a more appropriate way to
achieve this variation in shear strength with radial distance was to do consolidation analysis
using a work-hardening soil model (the Modified Cam-Clay model) before running any lateral
load analysis. Studying the effect of the shear strength increase on the limiting lateral
resistance of the soil in this way, there was found to be a 35% increase in the latter. The
finite element mesh is shown as Figure 5.21. There are 396 linear strain triangles with 348
free degrees of freedom. Variable shear strength used in the analysis and the results are shown
in Figures 5.22-23 respectively. (The shear strength at the end of consolidation is denoted by
cr'). The pile radius is r() and r is the radial distance from the centre of the disc.)
The results in Figure 5.23 show that if the pore water in the soil around the pile is allowed
to drain, a higher load is required for failure than otherwise. This implies that pile groups
would also require higher loads at which the plastic zones around the individual piles will
5-27
Chapter 5 Analysis of Single Piles under LAteral Loading
/<.//
/977 --.... / P77 ._
Figure 5.21: Finite element mesh with variable shear strength around the pile
1.7 ..---------------------------,
1.6
1.5
1.4
1.2
1.1
0.9 '--__ -'-__ --' ___ -'-__ --' ___ -'-__ --'-__ ---.J
o 2 4 6 8 10 12 14
rlro
Figure 5.22: Variation of shear strength around the pile
start to overlap causing greater pile-sail-pile interaction. Later on a simplified analysis was
performed in which the logarithmic variation of shear strength and the initial stresses after
consolidation were incorporated in the finite element program. This modelling was achieved
by fitting mathematical expressions to the stresses and the increase in shear strength predicted
by Randolph, Carter and Wroth (1979) and then implementing in the finite element solution.
Thus radial inhomogeneity was duly accounted for in the analysis. (Another approach to
5-28
Chapter 5
Analysis oj Single Piles under Lateral Loading
-0
* ;::l
(/)
2
;'1.5
0;
......
......
---
-0
C'::l
o
-
II
I-
o
U
..;:: 0.5
-0
C'::l
o
o ~ - - - - - - - - ~ - - - - - - - - ~ - - - - - - - - ~ - - - - - - - - - - ~ - - - - - - ~
o 0.1 0.2 OJ 0.4 0.5
Did
__ constant Su
--+- variable Su
Figure 5.23: Load factor displacement curves (variable shear strength)
model radial variation of the soil shear strength is to relate it with the change in the
hydrostatic pressure in the soil as the finite element proceeds. Randolph, Carter and Wroth
-0
*
L2
< J5 1
*
'<t
0\
~ 0 . 8
......
--- ""g
06
-
II
gO.4
~
-00.2
C';$
.2
o ~ ______ ~ ________ ~ ______ ~ ______ -J ________ ~ ______ -J
o 0.2 0.4 0.6
Did
__ predicted curve after consolidation
-.A.- predicted curve after installation
Figure 5.24: Comparison of p-y curves
0.8 1.2
___ empirical curve
(1979) suggested an approximate relationship between the change in the hydrostatic pressure
and the soil shear strength at the end of the consolidation phase. This relationship when
implemented in the finite element analysis provides a more realistic response of laterally
loaded single-piles as well as pile groups.) The mesh dimensions were kept at 20r
a
from the
5-29
Chapter 5 Analysis of Single Piles under lAteral Loading
centre of the pile (see Figure 5.21). The results are shown in Figure 5.24. The figure shows
that the p-y curve obtained numerically match closely with the p-y curve proposed by Matlock
(1970). Thus this pilot study shows that better p-y curves can be derived using plane strain
finite element analysis. Although a non-linear model seems more appropriate in this case, as
we can see that after the end of the consolidation, soil has undergone a loading-unloading
cycle, and its behaviour will tend to be strongly non-linear rather than linear upon re-Ioading.
(The strain response on unloading and reloading is hysteretic. There is no longer one to one
relationship between stress and strain in the region. The soil stiffness is influenced by the
recent stress history. Therefore, nested yield loci model (Houlsby and Chow (1994) seems
to be more suitable to model soil behaviour at the end of the consolidation phase.) Finally,
Figure 5.24 compares the load displacement responses both immediately after installing the
pile and after the end of the consolidation phase. The figure shows that the pile capacity
immediately after installation is almost half of that at the end of consolidation phase.
In the analysis of Randolph, Carter and Wroth (1979). the choice of shear modulus G in the
soil model adopted, influences the generation of excess pore pressures. They pointed out that
if the shear modulus was tied to current effective stress level (rather than to the past
maximum stress level) then lower pore pressures were generated. They did not recommend
that shear modulus be linked with current in situ stresses which results in low values of I
r
at high values of OCR. This seems to suggest that the initial effective stresses in the soil may
affect the resulting shear strength increase after the end of consolidation around the pile and
in tum the ultimate limiting soil resistance of the soil. Another point made by Randolph,
Carter and Wroth (1979) about their analysis concerns taking account of the assumption of
plane strain conditions and ignoring the shear stresses in the vertical plane (r-z planes). Shear
stresses will be present not only due to residual driving stresses but also to balance the high
total vertical stresses which occur during pile driving. Field measurements have shown that
excess pore pressures, greater than the existing effective overburden stresses, are generated
during pile installation. The resulting changes of total vertical stress must be balanced by
shear stresses in the soil around the pile in order to maintain vertical equilibrium. The effect
5-30
on the resulting increase in soil shear strength around the pile of the presence of these shear
stresses, and the resulting out of plane (vertical) soil movement, has yet to be investigated.
5.9 Interface analysis
A purely cohesive interface, zero thickness element with some allowable tension cut-off value,
developed in Chapter 4, is used here to model pile-soil interaction. The need to use this type
of element arises for the following reasons. Firstly, it is capable of capturing the deformations
close to the pile where there are very high stress concentrations and strain gradients i.e. it
smooths stresses at the discontinuous boundary. Secondly, it is capable of studying tension
at the back of the pile. Thirdly, it is possible to vary the adhesion between the soil and the
pile. Although this possibility was not tried out. It was found that even a coarse finite element
mesh with a layer of interface elements next to the pile predicted the exact limiting load given
by Randolph and Houlsby (1984). Earlier, Kooijman (1989) also showed that the predicted
limiting lateral load using a refined finite element mesh was the same as from a coarse mesh
with the addition of a thin layer of continuum elements or interface elements to model the
pile-soil interface. In the analysis of a rigid disc displaced laterally in a cohesive soil, the soil
is highly constrained next to the rigid disc so does the mesh close to the rigid disc, which
explains the over-stiff response in the absence of an interface element. In the absence of an
interface element, the number of free degrees-of-freedom are not enough in the coarse mesh
to predict an exact collapse load, which results in a gross over-estimation of the limit load
(see Section 5.4). The presence of the cohesive interface elements provides additional degrees-
of-freedom in the coarse mesh without any increase in the incompressibility constraints, and
so corrects the over-prediction by eliminating the over-stiff response. The finite element mesh
used in the analysis is shown in Figure 5.25. The use of interface elements need a careful
choice of interface stiffnesses. Choice of lower stiffnesses makes the interface flexible and
the predicted limiting load remains lower than the exact limiting load given by Randolph and
Houlsby (1984). A predicted lower value for the limiting load should not mean that the finite
5-31
element method under-predicts the limiting lateral load (which can not be right). The under-
prediction can be improved (though not completely) by using a higher order integration
scheme for the interface elements i.e. by over integrating the stiffness matrix or with the
choice of appropriate interface stiffnesses as is the case for the analysis described in this
section. Interface shear and nonnal stiffnesses used in the present analysis are 71500 kNlm
and 143000 kNlm respectively. Sometimes the choice of interface stiffnesses makes the
interface stiffer than required which brings ill-conditioning in the calculations and instabilities
in the final solution.
/<}/
Figure 5.25: Finite element mesh with a layer of interface elements around the rigid
pile,
The load displacement curve with full tension at the back of the pile is shown in Figure 5.26.
5.10 Effect of no-tension behaviour of pile-soil interface
As the lateral load on the pile increases, the strength of the soil at the back of the pile is
5-32
1.1
'"0
1 -
... ..
... ... ... ... ... ...
* ::l
0.9 (/.:l
-
* ~ 0.8 -
0\
--
0.7 :-
.......
-
0.6 '"0
ro
0
0.5
-
-
,.
II
I-.
0.4
0
.....
0.3
u
~
0.2
'"0
ro
0.1
0
......
0
0 2 3 4 5 6 7 8 9 10 11
old
__ zero tension ___ full tension
Figure 5.26: Interface analysis for a single rigid pile
Figure 5.27: Vector plot of displacements when a gap opens up at the back of a
laterally loaded pile
exhausted and unable to take any further tension. At this stage of tension cut-off, a gap
between the pile and the soil opens up with a sudden decrease in the lateral load capacity of
the pile. (A vector plot of displacements when a gap between the pile and the soil opens up
5-33
is shown in Figure 5.27.) The interface element predicts the limiting soil reaction which
matches with Randolph and Houlsby's (1984) solution (see Figure 5.26) and if there is a
tension cut-off, the value of limiting soil reaction reduces to half of its value with full tension
i.e. as in the exact solution given by Randolph and Houlsby (1984). The tension cut-off value
specified in the analysis is 250 kPa. The adhesion between the soil and pile can easily be
varied by changing the value of the shear strength of the interface.
5-34
Chapter 6 Plane Strain Analysis oj Pile Groups
CHAPTER 6
PLANE STRAIN ANALYSIS OF PILE GROUPS
6.1 Introduction
Groups of piles are generally arranged in square/rectangular grids, or in concentric circles,
with the spacing of the piles determining the size of the foundation. A pile group may contain
battered piles and may be subject to simultaneous axial, lateral, moment and possibly,
torsional loads. In this chapter we will be dealing with laterally loaded pile groups only. In
the absence of analytically derived solutions for the ultimate soil resistance in a laterally
loaded pile group, the focus will be mainly on the numerical methods currently in use. These
methods can be categorised into three main types as described below.
In Type 'one' methods, the piles are treated as beams and the soil is represented by non-linear
springs. The characteristics of the springs may be determined by empirical rules e.g. the p-y
curves proposed by Matlock (1970) for saturated soft clays or by plane strain calculations as
described in this chapter.
Type 'two' methods treat both the pile and the soil as a 'continuum'. The behaviour of the
piles grouped in this case may be detemlined by any suitable method e.g. the integral
equation!boundary element method or the finite element method. (Generally speaking the
boundary element method is most appropriate for elastic soil behaviour and the finite element
method is most appropriate for cases where non-linear soil models are assumed. Detailed
6-1
Chapter 6 Plane Strain Analysis of Pile Groups
discussion of previous work on these methods is given in Chapter 3. The obvious problem
with these methods is that the amount of computation required is immense in any multi-pile
group. Therefore. it is usual to perfonn the calculations for two piles and to work out the
interaction factor ex. between them. The analysis of a two-pile group is then extended to any
number of piles in a group as proposed by Poulos and Davis (1980) and is known as elastic
interaction approach for the prediction of pile head displacements.) To model elasto-plastic
soil behaviour, the finite element method is the best among those available, although it too
has limitations. Its limitations derive from the problem of knowing where to fix a 'boundary'
to the soil continuum, and from the problem of knowing which of the soil models will give
the best results in any given situation. It follows from these limitations that care need to be
taken in applying this method. But does it matter whether the finite element analysis is done
for three-dimensions or only for two-dimensions? Three-dimensional analysis is prohibitively
expensive in tenns of volume of calculations and corresponding computer memory involved.
However, the technology is developing all the time and it may be possible to run a three-
dimensional analysis more economically in the not-too-distant future. Until a new generation
of fast high memory computers becomes available, a two-dimensional analysis is clearly
preferable.
Type 'three' methods use interaction factors to describe the interactions between the piles in
a group. It is assumed that an interaction factor is valid for any two neighbouring piles in the
group regardless of their position in the group and regardless of the total number of piles.
Interaction factors may be computed using a type 'two' method. The interaction factor
approach is certainly applicable at the pile cap level (pile head). That is good enough insofar
as it provides the vertical settlement and lateral deflection of the pile head. However, it
ignores the variation in deflection and bending moment along the length of the piles. Piles
have to be considered as having a constant cross-section, and the pile head restraints have to
be either fully fixed (i.e. no rotation) or fully free (i.e. no bending moment). The basic
6-2
assumption is that all piles in the group behave identically i.e. the piles are spaced equally
around a circumference of a circle, each pile displaces equally and carries equal load (a
symmetrical group) and the principle of superposition is valid. The limitations of the elastic
approach to pile interaction under lateral loading were exposed by the model pile tests carried
out on a centrifuge by Barton (1982). These tests showed that the leading piles of groups
embedded in sand carry a higher proportion of overall applied load than the trailing piles even
at load levels lower than the failure load. For pairs of piles at close spacing (s/d'=2) , Barton
found that at a lateral displacement of 0.2d, the applied load was shared 60% to the front
piles and 40% to the rear piles. By contrast, Matlock et al. (1980) had reported results of tests
on pile groups in clay, in which applied load was found to be shared quite uniformly. Clearly,
the adequacy of the interaction factor approach depends on a number of factors such as load
level, soil type, and most importantly pile spacing. The assumption that the interaction factor
between the piles remains elastic, and the principle of superposition applies, is not valid at
all for closely spaced pile groups: the elastic interaction approach works well at small
displacement levels. (We shall see in the following Sections 6.2-3 (Figure 6.5 and Figure
6.11) that at close spacings, the load shared by individual piles in a 2-pile group and a 2x2
pile group is equal only at small displacements.) Such simplification of the problem, together
with the other assumptions and approximations entailed, means that the resulting solutions
based on elastic interaction approach are not reliably accurate.
Chapter 5 analysed pile-soil interaction for a single-pile and demonstrated that a two-
dimensional idealisation corresponds very well with real behaviour beyond a depth of 3.0
times the pile diameter below the surface (Broms (1964a, and when the mesh boundaries
are more than 10 times the pile diameter from the centre of the pile. The same conditions
should apply in the finite element analyses of pile-soil interaction for pile groups. It is the aim
of this chapter to test whether that is the case.
The chapter is divided into five further sections as follows: Section 6.2 describes the finite
element analyses of 2-pile groups under two loading conditions - when the load acts in a
6-3
direction parallel and perpendicular to the line joining the pile centres of the 2-pile group. A
comparison is made between soil resistance around a single-pile and a 2-pile group in order
to determine the most effective spacing to meet the specified loading conditions. Analytical
expressions are also suggested as a way of comparing and checking results from the finite
element method. Section 6.3 presents plane strain analysis for a 2x2 pile group, and it
confirms the conclusions drawn from Section 6.2. Section 6.4 describes interaction factors
between piles where the effective pile spacing (Figure 6.7, Sections 6.2 and 6.3) is not
possible. It is found that the inter-relationships between the piles spaced otherwise is not
linear at all, rather it is highly non-linear. An iterative analysis procedure is proposed to
predict load-displacement response of a symmetrical pile group. The analysis procedure makes
use of non-linear interaction factors determined in Section 6.4. Section 6.5 describes plane
strain analyses of pile groups to obtain theoretical p-y curves which take account of non-linear
pile-soil-pile interactions. Section 6.6 tries to explain theoretically why certain results obtained
from field tests carried out on a 5-pile group and a lO-pile group drift from what is
theoretically predicted by others. The methods of back analysis that have been tried to date
were discussed in Chapter 2 and the difficulties with them pointed out. The numerical
approach used here enables us to quantify the interaction attributed to the overlapping plastic
zones around the piles. Section 6.6 also suggests an alternative method in the absence of the
ideal three-dimensional finite element analysis. A two-dimensional finite element analysis
combined with the p-y curve analysis enables the engineer, particularly one designing pile
foundations for offshore structures, to design with a greater assurance that the foundation will
be sound and economical.
6.2 Plane strain analysis of 2-pile groups
The limiting load for a pile in a rigid 2-pile group where both piles undergo equal
displacements could be higher or lower than that for a single-pile depending upon the
6-4
direction of loading. In the discussion following, a 2-pile group is analysed for two loading
conditions - when the applied load is parallel to the centre line of the pile group and when
it is perpendicular to the centre line. The finite element meshes used in the analysis are shown
in Figures 6.1 and 6.2. The symmetrical nature of the problem means that only half of the
geometry needs to be discretised.
The mesh dimensions are chosen in a way similar to that adopted for the single-pile analysed
in the preceding chapter. The element type, the boundary conditions and the material
properties are also the same.
/P77 III- /P77
Figure 6.1: Finite element mesh for a 2-pile group (loading parallel to the centre
line of the pile group)
The curve for the single-pile in Figure 6.3 shows some considerable discrepancy (%error =
14) with respect to the analytical solution given by Randolph and Houlsby (1984). This is due
to the coarse discretisation of the mesh used in the analysis. A mesh of similar density is used
in the analysis of the two-pile group. Therefore, a comparable error does exist in the curve
for the two-pile group. The errors are also due in part to the particular type of finite element
chosen for the analyses. Six-noded triangular elements are quite suitable for the collapse load
calculations in nearly incompressible materials (i.e. undrained clays with which we have been
6-5
/<}'/
/977 ~ /P77
Figure 6.2: Finite element mesh for a 2-pile group (loading at right angle to the
1.1
~
......
1
0.
u
0.9
Ob
c::
'-0
0.8 CIl
C<:I
....
.B
0.7
"d
0.6
C<:I
0
0.5

C<:I
0.4
.9
II 0.3
.9
u
0.2
~
0.1
"d
C<:I
0
.9
0
centre line of the pile grou p)
0.1
__ single pile
d
I' 1
~
... e. ~ = O
- - ~ -- ~ ' . - . - . ------- : ~ ' - . - . ' - - ~ - - 0.49
s=
0_2
o Id
0.3
______ avo load/pile in a 2-pile group
Figure 6.3: Load displacement curves (2-pile group; B=O)
0.4
dealing), and have the additional advantage (as compared to the quadrilateral family of
elements) that the polynomial functions defining the displacements within the element are
complete. However, they fail to capture the strain fields around the piles and need finer
6-6
discretisation to do so. At a centre-to-centre pile spacing of two pile diameters (see Figure
6.3), if the loading is parallel to the centre line of the pile group, the limiting load for a pile
in a group is nearly half that for a single-pile as described in Chapter 5. The reduction is due
to the overlapping plastic zones between the two piles which modifies the stiffness of the soil.
Randolph and Houlsby (1984) showed that the expected plastic zone in front of a laterally
loaded pile extends to a distance of one diameter from the face of the pile. This means that
there will be a complete overlap of plastic zones at a centre-to-centre pile spacing of two pile
diameters (see Figure 6.4). At a spacing of three pile diameters the boundaries of the plastic
zones will be just touching each other; they move further apart until, at a spacing of about
six pile diameters, the lateral reaction approaches that of a single-pile. (The pile spacing at
7d is derived analytically: see below, equation (6.1) and Appendix 6A.)
pile
'-'-'-'- CL
spacing = 2*d
Figure 6.4: Overlapping plastic zones in 2-piJe group (loading parallel to the centre
-line)
The load shared by the individual piles in a 2-pile group (spacing = 2d; loading parallel to
the centre line of the pile group) is shown in Figure 6.5. Load shared by each pile is uniform
at small displacement levels.
6-7
. 1.4
'0.
o 1.2
"Eb
C
.-
</)
C\l
....
.E
~ 0.8
o
~ 0.6
o
-
II
I-< 0.4
o
....
u
~ 0.2
~
o
- 0
o 0.5 1.5 2
old
__ group
-.- front pile -.t.- rear pile
Figure 6.5: Load distribution among the piles in 2-pile group
(spacing = 2d)
2.5
The limiting load for a pile in a group may be higher than that for a single-pile, if the loading
is perpendicular to the centre line of the pile group. The increase is due to the mobilisation
of a larger plastic zone around the pile group which behaves like a large diameter single-pile
instead of two piles in a group. This is shown in Figure 6.6 where angle ~ is the direction
of loading. Figure 6.7 shows the results for variable pile spacing up to 10 diameters. Besides
the finite element analyses, two simple analytical expressions can be used (with reservations
to be explained below) to work out the spacing at which the pile group will attain the limiting
load for a single-pile. The two expressions, when the load is in line with (parallel to) the pile
group centre line and when it is perpendicular to that line are, respectively:
and
p = ~ + 11.94
d {3
6-8
(6.1)
(6.2)
In both cases s is the pile spacing and, d is the diameter of the pile. (For details, see
Appendix 6A.)
o1.2r---------------------------, ,.......
's.
o
~
.....
"J
c<:i
..8 0.8
"0
c<:i
o 0.6
~
..9
II 0.4
.....
o
. ~
- 0.2
-g
o
.--.
o - - - - ~ - - ~ - - ~ - - - ~ - - ~ - - ~ - - ~ - - ~
o 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
__ single pile B /d -..... av.load/pile in a 2-pile group
Figure 6.6: Load displacement curves (2-pile group; B =90)
Whereas expression 6.1 is quite successful in that it suggests that at a pile spacing of
approximately 7 times the pile diameter, the pile group attains the limiting load for single-
pile, expression 6.2 is not so successful. Possibly, as the expression is derived on the basis
of size of a plastic zone for a single-pile whereas in reality we have two piles and their
interaction will effect the size of the plastic zone. It predicts higher limiting loads; and it does
not fully trace the variation of the limiting load with the pile spacing. It works only as far as
sid = 3, at which point the limiting load drops abnlptly to the value for a single-pile.
The results show that when the load is applied in line with the two piles, the limiting load
decreases with closer spacing; as the spacing is increased, the value of the limiting load
gradually reaches that for a single-pile. The reason for this is that plastic zones around each
pile move further apart as spacing increases. The converse is the case when the load is
6-9
30
o
o 2
_ fe sollllion for single pile
o
2x2 pile group @ = 0)
4 6
sid
8
o 0
10
__ 2-pile group (13=0) -.tr.-- 2-pile group (i3=90)
o 0
-e- theoretical boundary (13=0) -A-- theoretical boundary (13=90)
Figure 6.7: Variation of ultimate lateral reaction with variable pile spacing
12
perpendicular to the centre line of the pile group: i.e. the limiting load increases at closer
spacings - the piles tend to act as a large diameter single-pile - and decreases gradually to the
value of a limiting load for a single-pile at wider spacings. Note that the 2-pile group is stiffer
for the loading case when than the one when
The dotted line in Figure 6.7 showing a possible curve for the 2-pile group when P = 90.
The spurious data point falling below this curve is due to numerical scatter: as explained
earlier, the discretisation of the meshes used was too coarse to give ideal results. It is quite
certain from the analytical solution (6.2) shown in Figure 6.7 that at pile spacing of three pile
diameters, there can be no overlapping of plastic zones (see Figure 6.8) and therefore the
limiting load must be the same as for a single-pile. In other words, the numerical results
(using (6.2) are sound on this poim.
There is also non-convergence between the results for a single-pile and those for the pile
group loaded parallel to the centre line. As shown in Figure 6.7, this non-convergence is very
small.
The results summarized in Figure 6.7 clearly imply that the pile group is much stronger when
6-10
plastic zone plastic zone
plastic zone plastic zone
spacing = 2*d
I'
Figure 6.8: Plastic zones in a 2-pile group (loading perpendicular to the centre-line)
the applied load is at right angles to the centre line. It follows that for pile groups in
rectangular patterns, where piles cannot be spaced far apart, there will be weakness on one
side of the pile system. In this situation, and assuming spacing at more than two pile
diameters, a circular arrangement of piles avoids this potential weakness. Returning to the
results in Figure 6.7: ideally, piles should be arranged in a circle of diameter 6 times the pile
diameter with a centre-to-centre spacing of 3 times the pile diameter.
6.3 Plane strain analysis of a closely spaced 2x2 pile group
Plane strain analysis of a 2x2 pile group confirms the danger, alluded to above, of block
failure when piles are closely spaced. The finite element mesh for loading parallel to the
centre line of 2x2 group is shown in Figure 6.9. As before, the mesh dimensions are fixed as
in Chapter 5: the finite element type, boundary conditions and the material properties are the
same as calculation in previous section. The results are summarised in the three curves in
6-11
Chapter 6 Plalle Strain Analysis of Pile Groups
Figure 6.9: finite element mesh for a 2x2 pile group
u 1.1
:-;:::i
0..
1
u
e.o
0.9
c::
.-
U'.l
0.8
<'!:S
....
.8
0.7
"'0
0.6
<'!:S
0
~
0.5
"'0
<'!:S
0
0.4
.-
II
0.3
....
0
...-
0.2 u
~
"'0
0.1
<'!:S
.9 0
0 0.05 0.1 0.15 0.2 0.25
o /d
__ single pile -+- avo load/pile in a 2-pile group
-.t.- avo load/pile in a 2x2 pile group
Figure 6.10: Load displacement curves (2x2 pile group)
Figure 6.10. As mentioned earlier, the curve for the single-pile shows a significant error of
about 14% in the predicted limiting load; so a comparable error does exist as well for the
6-12
curves for the 2-pile and 2x2 pile group. In both groups the piles are spaced at 1.5 times pile
diameter, centre-to-centre. Each pile develops a plastic zone which, at this spacing, overlaps
with the plastic zones of the others. The results clearly show that the 2-pile and 2x2 pile
group reach the same limiting load per pile, which is less than half the load for a single-pile.
If we accept the error in the finite element results and assume that it could be reduced to
acceptable limits (see above, Section 5.2) by using a finer mesh, then the results tell us that
in the case of a rigid 2x2 pile group, the interaction between any two piles in line with the
applied load is more significant than the interaction between any two piles whose centre line
is at right angles to the applied load. The results also show that, initially, i.e. when the
displacements are smaller, the response of each pile in the 2x2 pile group is more flexible
than the one in the 2-pile group. Taken together, the results imply that in large
square/rectangular pile groups where piles are arranged in longitudinal and transverse rows,
the interaction between the piles in individual pile rows parallel to the direction of loading
is going to be critical in the design of the foundation. Therefore centre-to-centre pile spacing
should not be so small as to provoke block failure. The load shared by individual piles in a
2x2 pile group (spacing = 2d) is shown in Figure 6.11. The load distribution among the piles
is uniform at small displacement levels.
6.4 Determination of pile-soil-pile interaction factors
In this section, non-linear pile-soil-pile interaction factors are discussed. A single plane strain
finite element analysis of a pile represented as a rigid circular disc proved remarkably helpful
in determining the interaction factors for variable pile spacings and for lateral loading at
various departure angles (13).
Figure 6.12 shows the finite element mesh used. As before, half of the geometrical
arrangement is enough for the analysis because of symmetry. A very fine discretisation of the
mesh is necessary close to the pile due to high concentration of stresses. The numerical over-
6-13
2
"'"' 0.5
~
..s
--g
o
-
O ~ - - - - - - - - ~ - - - - - - - - ~ - - - - - - - - ~ - - - - - - - - ~ - - - - - - ~
o 3 4 5
__ front pile
-.-group
-.to- rear pile
Figure 6.11: Load distribution among the piles in a 2x2 pile group (spacing = 2d)
/ P77 I' 01
Figure 6.12: Finite element mesh used for the determination of interaction factors
estimation of the limiting load as detennined by analytical method is below 5% i.e., the
adopted mesh is better than the one used in the previous chapter. (The load displacement
curve obtained from the finite element analysis using this mesh is shown in Figure 6.13.)
6-14
1.2 .--------------------------------------------------,
:;p analytical solution
::l
C/.l
* ~
0; 0.8
,.......,
,.......,
--
'"0
B 0.6
-
II
;.....
00.4
......
C)
~
~ 0.2
o
-
0.1 0.2 0.3
8/d
Figure 6.13: Load displacement curve
0.4 0.5 0.6
Variable centre to centre pile spacing in Figure 6.12 is represented by establishing nodal
points in radial directions at departure angles ( ~ ) of 0, 22.5, 45, 67.5, 90 degrees. These
nodal points represent pile spacings of up to 10 pile diameters.
When the 'disc' is moved laterally, it induces lateral movements in the same direction in the
established nodal points. These lateral movements are in addition to movements which occur
as a result of other loads (if any) acting on those points. An interaction factor can therefore
be defined as the ratio of the additional displacement caused by an adjacent pile to the
displacement caused by its own loading. This type of interaction has been studied by
Mandolini and Viggiani (1992) for the case of a pair of axially loaded vertical piles.
Figure 6.14 shows interaction curves for loading at various departure angles, at a spacing of
1.5 times the diameter of the pile (similar curves can be drawn for other spacings). As we can
see, the interaction factors are linear for small displacements (Le. elastic interactions for small
\
load levels) and after that increase non-linearly with the applied displacement/load. In other
words, the figure shows that the interaction factors are not unique for all load levels. Curves
for departure angles ( ~ ) = 0, 22.5 converge toward the same level of ultimate interaction.
6-15
0.4
-0.2
-0.3
o 0.1 0.2 0.3 0.4 0.5 0.6
l/d
P == 22.5
0
-,t.- P ::: 45
0
--
P = 90
Figure 6.14: Interaction curves
Interactions with <P) :: 45, 67.SO, 90 reach a maximum value which is lower than the
interaction values for (P)= 0, 22.5 which diminish as the load increases. The diminishing
values for departure angles (P)= 67.5, 90 become negative and reach maximum negative
value. This trend is as might be expected. When the lateral load approaches a limiting value,
the soil tends to flow round the piles which causes them to move in a direction opposite to
the applied load. This trend explains why side-by-side pile groups are stiffer. The movement
of each of the piles impedes the other.
The results suggest that, in order to keep the effect of pile-soil-pile interaction to a minimum,
the optimum angular pile spacing is 45.
The non-linear variations of interaction factors with variable pile spacing for lateral loading
for various departure angles are shown in Figure 6.15. These factors correspond to o/d=20%,
the value at which the soil starts yielding. The elastic interaction factors given by Poulos and
Davis (1980) are also shown for comparison.
The diminishing trend is striking and it can be seen from the figure that, for departure angles
6-16
Chapter 6 Plane Strain Analysis of PUe Groups
0.8 r------------------------------.
0.7
0.6
......
c.o 0.5
-N
c.o 0.4
II
i:::5 0.3
0.2
0.1
o
o
__ =0
o
o

Figure 6.15: Interaction factors
2 3 4
sid
___ =: 22.5
o
o
-e- = 90
5 6 7 8
o
=45
= 0, 22.5, the response of a pile is still about 10% more flexible than a single-pile,
even at a spacing as large as seven pile diameters. This flexibility occurs when only one pile
is loaded in the pile group; if all the piles are loaded at the same time, the response will be
even more flexible.
The results taken together show the complex nature of the interactions among closely spaced
piles in a pile group: the assumptions of the interaction factor approach (type 3 methods, see
Section 6.1 above) are not applicable in this situation. However, for a fixed headed pile group
in which all the piles displace equally, the method proposed by Poulos and Davis (1980) for
the analysis of a general pile group can be used to determine the load displacement curve of
a circular pile group adopting an iterative procedure (Section 6.4.1), using the non-linear
interaction factors given Figures 6.14-15.
6.4.1 Approximate prediction of load-displacement curve for circular pile groups
In this sub-section we will describe how non-linear interaction factors derived above can be
6-17
used to construct an approximate load displacement curve at pile cap level of a rigid circular
pile group in which piles are located only at the periphery of a circle and displace equally.
We will make use of the following expression proposed by Poulos and Davis (1980):
(6.3)
where
8 is the displacement of a single-pile under unit horizontal load
8,1: is the displacement of pile k
H
j
is the load on pile j
U
kj
is the interaction factor between piles k and j
and the total load on the group is H G given by:
II
(6.4)
In the case of equal displacements, the n equations from the above expression together with
the equilibrium equation may be solved for the unknown loads and the group displacement.
If we use interaction factors at failure given in Figure 6.15 we can find out the ultimate lateral
load taken by the group. The ratio of average load per pile in the group to the ultimate load
of a single-pile will give us the capacity reduction factor for that particular group. Applying
this capacity reduction factor to the load-displacement curve of the single-pile will give us
a new load displacement curve (modified p-y curve). (The derivation of p-y curves using the
finite element method is described in Chapter 5), Now, a load displacement curve can be
constructed for the pile group following an incremental procedure in which for each value of
6-18
Plllne Strain Analysis of Pile Groups
incremental displacement, an average load per pile in the group is obtained from the above
expressions 6.3-4, using the appropriate interaction factors for that particular displacement
from Figure 6.14. (Mathematical expressions to the interaction curves shown in Figure 6.14
can be fitted and appropriate values of interaction factors will be calculated directly from
these expressions during the analysis.) But in order to have the correct stiffness each
calculated load-displacement couplet will be matched iteratively with the modified p-y curve
obtained above. Validation of this method is left to future work. In the following Section 6.6,
we back analyse (using the same concept) field tests carried out by Matlock el al. (1980) on
two circular pile groups using the p-y curve analysis. The p-y curves used in this analysis take
account of elasto-plastic pile-soil-pile interactions within the pile groups. Section 6.5 describes
how these curves are obtained using the finite element method.
6.5 Plane strain analysis of multiple pile groups
The results from the preceding Section 6.4 suggest that any analysis which does not consider
the non-linear pile-soil-pile interactions among closely spaced pile groups is liable to predict
pile group response inaccurately. Many methods currently in practice fail in this respect,
except the finite element method. In this section, plane strain finite element analyses of
multiple pile groups are described. Plane strain idealisation is made on the basis that in the
event of horizontal displacement of the pile, any vertical displacements will be negligible by
comparison. However, plane strain idealisation is not valid for the whole length of the pile
but only after a certain depth (approx. 3d). Theoretical p-y curves are obtained from these
analyses which take account of non-linear pile-soil-pile interactions within the pile groups.
The effect of pile-soil-pile interaction can be incorporated by introducing a strength
modification factor to the undrained shear strength of the soil to make up the reduction in the
soil resistance in a pile group. Strength modification factors are obtained through plane strain
finite element analyses described in Section 6.5.1-2. Although, we have discovered in the
6-19
Plane Strain Analysis of Pile Groups
previous chapter that we do not require a soft layer for the analysis of a single laterally
loaded pile, we have used a correcting layer in these analyses used for the determination of
strength modification factors. A correcting layer is not going to affect these factors.
Field lateral load tests on two pile groups of 5 and 10 piles, arranged in a circle in both cases,
carried out by Matlock et al. (1980) are chosen here for analysis. In these tests, all piles were
of steel, 13.4m long, with an external diameter of 168mm, and a wall thickness of 7.1mm.
They were driven to a penetration depth of 11.6m below the bottom of an excavated pit. The
pile centre line spacing were 3.4 and 1.8 pile diameters for the 5-pile and 10-pile groups
respectively. In order to simulate pile head restraints typical of an offshore platform, pile head
deflections were enforced at two elevations, one just above the groundline, and the other near
the pile top to simulate head restraint of a pile cap. During the loading tests, pile group
deflections at two supports above the groundline and the total group load, as well as
individual pile head shear and bending moments along the individual pile shaft, were
measured at various load levels.
In the plane strain finite element analysis which follows, static monotonic loading is
considered and the soil properties (Su = 21kPa; G = 2070kPa) are taken from the paper by
Matlock et al. (1980). A nearly incompressible material behaviour is assumed for the soil
(undrained soft clay) by taking Poisson's ratio of the soil to be 0.49.
6.5.1 Plane strain analysis of S-pile group
In this section, a 5-pile group is analysed as a two-dimensional problem to obtain a strength
modification factor for the group. The finite element mesh used in the analysis is shown in
Figure 6.16. The lateral soil reaction in the 5-pile group is shown in Figure 6.17.
Exploiting, as before, the symmetry of the problem, half of the geometry is discretised.
Boundary conditions and direction of loading are also shown in the figure. The average
response of a pile in a group is flexible compared to the isolated single-pile. This is due to
6-20
Figure 6.16: Finite element mesh of the 5-pile group
"'d
* ::l
CI)
* -<;t
0'-
,......
,......
-
"'d
ro
.s
II
I-<
0
.....
u
~
"'d
ro
.2
1.2
0.8
0.6
0.4
0.2
O.! 0.2 0.3
bId
__ av. load/pile in the 5-pile group
0.4 0.5 0.6
__ isolated single pile
Figure 6.17: Average lateral soil reaction in a pile group
the interaction of plastic zones around individual piles which tend to overlap and the ultimate
resistance offered by the soil is reduced. Figure 6.17 shows the flexible response of a pile in
the group, whereas Figure 6.18 shows that, if the soil strength is modified by 5%, the average
soil reaction in the group is the same as for a single-pile. In other words if we apply a
6-21
P/Q.ne Strain Analysis of Pile Groups
-0
~ 1.2 ,---------------------------,
-0
o
E 1
'-'
::J
CI)
;. 0.8
~
::: 0.6
:0-
C':l
..s 04
l!.
o
U 0.2
~
~ 0 L-___ ~ ___ ~ ___ ~ ___ ~ ___ ~ ___ ~
.9 0 0.1 0.2 0.3 0.4 0.5 0.6
b/d
__ av. load/pile in the 5-pile group __ single isolated pile
Figure 6.18: Load displacement curve (modified shear strength)
capacity reduction factor of 0.95 to the ultimate resistance of a single isolated pile we can
obtain an average soil resistance of a pile group.
Figure 6.19: Finite element mesh for the IO-pile group
6-22
6.5.2 Plane strain analysis of to-pile group
The finite element mesh of a lO-pile group is shown in Figure 6.19. Again, half of the
discretised geometry; and boundary conditions and direction of loading are shown. A
behaviour similar to the 5-pile group is observed. The ultimate soil resistance developed in
this group is shown in Figure 6.20. The response 'as expected' is much more flexible than
the 5-pile group since the pile spacing is 1.8 times pile diameter causing overlap of soil flow
around the piles. Figure 6.21 shows that if the soil strength is modified by 38%, the average
soil reaction in the group is the same as for a single isolated pile.
1.2.-----------------------,
"'0
* :::I
(/)\
*
7
0"-
.......: O,R
,......
-
"'0
15 0.6
-
II
8
0
.
4
]
-0 0.2
C'::l
.9
0.1 0.2 0.3
Old
__ avo load/pile in the lO-pile group
0.4 0.5 0.6
___ isolated single pile
Figure 6.20: Average lateral soil reaction in a pile group
6.6 Back analysis of field tests
Back analysis of the field load tests done by Matlock et al. (1980) on a single-pile and two
pile groups (5-pile and lO-pile) in soft clays, is presented in this section. Various attempts
have been made to back analyse these tests using methods with different degrees of rigour
e.g. Leung and Chow (1987) and Kooijman (1988) described in Chapter 2.
In this section a simplified approach to back analyse field tests using plane strain finite
6-23
-;::)
~ 1 . 2 r - - - - - - - - - - - - - - - - - - - - - - - - - . . . . ,
-g
E I
'-"
::l
CI')
::;. 0.8
0;
..-
...... 0.6
-
-;::)
<":I
.2 0.4
!!.
o
t) 0.2
~
~ 0 L-____ ~ ___ ~ ___ ~ _____ ~ _____ ~ _____ ~
o 0 0.1 0.2 0.3 0.4 0.5 0.6
-
old
__ av. load/pile in the 10-pile group
__ single isolated pile
Figure 6.21: Load displacement curve (modified shear strength)
element analysis is presented. It is shown how the p-y curve analysis can be used to mimic
the response of a pile group. (The variation of in situ stresses, soil stiffness and shear strength
with depth is duly accounted for in the analysis). A single plane strain analysis is found to
be sufficient to obtain a series of p-y curves along the length of the pile. In this section, a p-y
curve analysis was performed independently to model these tests. Soil parameters are taken
from the paper by Matlock et al. (1980). The undrained shear strength Sit of the soil at the
groundline was taken as 13.1kPa, increasing to 23.3lkPa at a depth of DAm, and thereafter
increasing with depth at a rate of 2.258kPa per metre. Based on the above information and
using a rigidity index, Ir ( where IT = GIS
u
) equal to 100 for soft clay, an approximate
linearly increasing shear modulus profile was used. The self weight of soil, y, was taken
20kPa per metre depth and earth pressure coefficient at rest, K
o
' was unity.
It was not necessary to carry out a series of plane strain finite element analyses at every metre
depth along the length of the pile group to obtain theoretical p-y curves. We have shown in
the previous chapter (Section 5.6) that the variation of IT can be represented by a single curve
6-24
if the horizontal axis (abscissa) of the load-displacement curve is non-dimensionalised with
I . It is further revealed in this section that load-displacement curves with variable K values
r 0
can also be non-dimensionalised in a similar way. Figure 6.22 shows such a curve and a
series of p-y curves are obtained from this single curve. Appropriate strength modification
factor (strength modification factors as detennined in sections 6.5.1-2) is then applied to these
p-y curves. Thus single plane strain analysis yields p-y curves which become iaput
-0
90
80
70
60

!?2.
"'g40
..9
30
20
10
00
5 10
Figure 6.22: Non-dimensionalised p-y curve
15 20 25 30
c*G/Su*d
for a standard non-linear beam-column equation which is solved numerically using iterative
solution procedure. This was done and the results are shown in Figures 6.23, 6.24 and 6.25.
(This method is similar to the one proposed by Kooijman (1989). Although, the method
suggested by Kooijman (1989) is more integrated than the one described in this section, the
p-y curves used in the former method are not as reliable as being used in the present method.
The theoretical response curves obtained by p-y curve analysis were then compared with the
field results. The p-y curve response is seen to be more flexible than the field results. The
obvious reason for this behaviour is the invalidity of plane strain assumption in the upper
regions of the pile.
In the absence of a complete three-dimensional finite element analysis of multiple pile groups,
6-25
35
30
_ ...........
,.-..,
25
Z
~
.........-
20
"0
(';j
...9 15
"'@
2
10
5
0.1 0.2 0.3 0.4
Old
-II- p-y curve analysis
__ experimental results
Figure 6.23: p-y curve analysis of a single-pile group
150 ,---------------------------,
50
0.1 0.2 0.3 0.4
Old
__ p-y curve analysis __ experimental results
Figure 6.24: p-y curve analysis of 5-pile group
the analysis presented here is a viable and useful way of working out lateral load capacity of
multiple-pile groups. Plane strain analysis enables a designer, from the predicted ultimate
lateral load capacity of the pile group, to optimise his design rather than designing over-
conservatively. Similarly, the modified p-y curve analysis proposed here enables a designer
6-26
Plilne Strain Analysis of PUe Groups
300 r------------------------------------------------,
250
;--., 200
Z
~
'-' 150
"0
cd
o
-100
50
0.1 0.2 0.3 0.4
old
__ p-y curve analysis __ experimental results
Figure 6.25: p-y curve analysis of lO-pile group
0.5 0.6
to work out flexural stiffness of the individual piles in a group. Three-dimensional finite
element analysis is inherently more economical than full scale field tests on pile groups, and
it allows for soil parameters to be varied which, plainly, is hardly practicable in the field.
Thus three-dimensional finite element analysis has much to recommend it. The drawback, as
we have stated, is that it is simply too bulky and too expensive in computing time to do at
present: the necessary improvements in computer technology are still some way off. Two-
dimensional plane strain analysis, with the reservations we have expressed above, is a present,
practicable alternative to both three-dimensional finite element analysis and field tests. It is
worth remembering that even when (if) the technology allows for three-dimensional analysis
of the 5 to 10 pile pile-groups usual in offshore applications, the two-dimensional alternative
will remain relevant for pile groups with more than 10 piles, usual in on-shore/inland pile
foundations.
6-27
Chapter 6 Plane strain analysis of pile groups
APPENDIX 6A
Analytical solutions for a 2-pile group
Pile-soil-pile interaction in a multi-pile group is a very complex phenomenon. It is difficult
to quantify analytically the influence of this interaction on the response of a multi-pile group.
A simplified analytical approach is here presented to study the behaviour of a 2-pile group.
We consider two loading cases: loading parallel to the centre line of the group and loading
perpendicular to the centre line.
Loading parallel to the centre line
S S
r r
2 2
Figure 6A.I: 2-pile group (loading parallel to the centre line)
Consider a group of two rigidly connnected circular discs of diameter d shown in Figure
6A-l
6A.1. Both the discs are pushed laterally by an equal amount parallel to the centreline of the
group. The two discs are located at a distance s apart.
It is assumed that at failure, the discs and the soil between them will move together as a
whole block. The forces acting on this whole block are shown in Figure 6A.l.
Considering the horizontal equilibrium of the forces, the net force F acting on the group will
be:
F = 2u + 5
r
where
't = 5
u
and
The horizontal equilibrium equation can be rearranged as:
F s 11.94
=
+
25 d d
13
u
or
p
s 11.94
= +--.
d
13
In the above expression P is the ultimate load taken by the group. It is clear from the above
expression that P depends upon the spacing between the piles. The group response will be
more flexible than the response of a single-pile at close spacings and will become equal to
that of a single-pile at a spacing of 6.89d.
6A-2
Loading normal to the centre line
S
r
2
I
\
\
\
"
I
d
,
(
s+d
\
'-i--
S
2
r
Figure 6A.2: 2-piJe group (loading perpendicular to the centre line)
The arrangement of the 2-pile group is shown in Figure 6A.2. In this situation, when the two
piles are too close to each other, a block failure is unlikely and the piles tend to behave as
a single larger diameter pile. Considering the equilibrium of horizontal forces acting on the
group (see Figure 6A.2), the net force F acting on the group will be:
F=S
r
where
The horizontal equilibrium equation can be rearranged as:
F = 11.94(1+":')
2Sud f3 d
6A-3
Chapter 6
or
p = 11.94 (1 < ~ : _ ) .
!3 d
Plane strain analysis of pile groups
The above expression shows that the response of a 2-pile group will be more stiff than the
response of a single pile at close spacings in contrast to a 2-pile group with loading parallel
to the centre line of the group_
6A-4
Chapter 7 Analysis of laterally loaded piles in frictional soils
ANALYSIS OF LATERALLY LOADED PILES
IN FRICTIONAL SOILS
7.1 Introduction
In order to proceed with the analysis of laterally loaded piles in frictional soils. we need to
determine the ultimate lateral soil reaction as discussed in Chapters 5 and 6. It is difficult to
determine the ultimate lateral reaction in a frictional soil due to dilatancy effects. When the
pile is pushed laterally, the soil in front of the pile tends to dilate and as a result, both the
normal stress and the shear stress in the soil continues to increase indefinitely. Likewise, the
lateral reaction of the soil also continues to increase due to increase in volume. In this
situation, the ultimate reaction cannot be defined and, in the absence of any analytical solution
available for frictional soils, there is no absolute certainty about the ultimate lateral reaction.
(It is not right, however, to assume that the volumetric strains will continue to increase
indefinitely, since dilation gradually diminishes with the increase in mean effective stress and
decrease in initial density.) This problem has not received much attention in the past.
However, Lane and Griffiths (1988) carried out a series of finite element analysis of laterally
loaded single piles in frictional soils using the Mohr-Coulomb model to represent the soil
constitutive behaviour. They adopted a non-associated flow rule (i.e. plastic strain increments
are not normal to the yield surface) with zero dilation as they found that dilation caused
'locking' of the meshes and an over-stiff, over-strong response at high pressures. Although
they did carry out parametric studies taking into consideration both rough and smooth pile-soil
7-1
contact, these results cannot be validated in the absence of an analytical solution or
experimental data. Lane and Griffiths (1988) also considered no-tension criteria (by inclusion
to the Mohr-Coulomb model) and showed the interesting result that a pile with zero tension
at the back of the pile had stiffer response compared to the one with full tension, and that
there is no well-defined ultimate lateral reaction. They explained that no-tension criterion
causes stress redistribution away from the back of the pile and this increases the apparent
stiffness of the soil in front of the pile, leading to its increased load response. As briefly
mentioned earlier in Chapter 5, the result seems to confirm our view that once the gap
between pile and soil opens up at the back of the pile, the problem of a laterally loaded pile
is no longer a problem of small displacement but one of a large displacement. (Large and
small displacement problems generally refer to geometric effects. This, however, is not strictly
a large displacement effect.)
In this chapter the ultimate lateral soil reaction is calculated using the Matsuoka (1976) model
to represent material constitutive behaviour of frictional soils, having considered both the
associated (plastic strain increment is normal to the yield surface) and the non-associated flow
rules. The finite element calculations are carried out for both constant as well as variable
friction and dilation angles. For the case of variable friction and dilation angles, the dilatancy
of frictional soils is considered to be a function of relative density ID and mean effective
stress as proposed by Bolton (1986). The dilation of the soil gradually decreases with
increasing stress levels and reducing ID during the analysis. The model captures successfully,
the post peak behaviour of the soil and tackles the problem of 'mesh locking'.
7.2 Background
In order to understand the subtle stress-strain behaviour of frictional soils, it is necessary to
define the angles of friction <1>- and dilation'll. The angle of friction <1>- describes the strength
of the soil and is expressed as the ratio of a shear stress to a normal stress:
7-2
sin4>' ==
where (j 1 ,(j3 are the maximum and minimum principal stresses. The superscript n is used
to show that these are 'effective' and not 'total' stresses. However, in fully drained analyses
(Le. pore pressures are zero) as is the case here in this chapter, 'total' stresses are equal to
'effective' stresses.
The angle of dilation 'V is expressed as the ratio of a volumetric strain rate to a shear strain
rate:
sin'V
-(1 +3)
El -3
in the case of plane strain (2=0), where 1'3 are the maximum and minimum principal
strains. The superior (.) represents strain rates. The term rate is used in the same sense as in
plasticity theory, in which the time increment is artificial. The (-) sign arises simply from the
convention that compressive stresses and strains are taken as positive in soil mechanics.
The angle of dilation is related to the density of the soil as dense soils have large dilation
angles. If we perform a simple shear test on a dense sand, the sand usually dilates i.e.
expands in volume as the test proceeds. The dilation usually takes place after a small initial
compression and denser samples have the tendency to expand more rapidly. A peak is usually
observed in the stress-strain curve followed by a reduction in shear stress at large strain. The
corresponding peak and large strain friction angles are denoted by Qp' and <1>cv' respectively.
The subscript' cv' indicates that the shearing takes place at constant volume i.e. no dilation
occurs.
If the same test is carried out at a number of normal stress levels, the angle of peak friction
decreases with increasing stress level. The peak friction angle approaches the large strain
friction angle <1>' cv at very high stress levels. Barden et al. (1969) showed from drained, plane
strain compression tests on sand at the same density that the angle of peak friction and the
7-3
angle of dilation were significantly affected by the confining stress. Similarly, Stroud (1971)
carried out simple shear tests on sands. These tests, covering a range of initial stresses and
densities show that denser soils exhibit higher friction angles but for a given density, the
angle of friction reduces slightly with increasing stress level. These test results suggest that
the angles of friction and dilation, density and stress levels are inter-related.
Various sophisticated theories explaining the relationship between the friction and dilation
angles include those given by Taylor (1948), Schofield and Wroth (1968), Rowe (1962) and
de Josselin de long (1971, 1988). The broad basis of these theories is the sawtooth model in
that the friction angle is the sum of the angle of friction at constant volume and the angle of
dilation. Bolton (1986) presented a comprehensive review of the experimental data on friction
and dilation angles, and suggested a very simple empirical fit to the data. This empirical
expression is given by:
(f p <l> -ct; + O. 8'V max
where 'V
max
is the peak angle of dilation.
Bolton showed that the above expression agrees well with the stress dilatancy relationship of
Rowe (1962). (A detailed discussion on these theories is given by Houlsby (1991).)
Similarly, the link between the angle of dilation, density and pressure has also been explored.
Based on the Critical State concept, Wroth and Bassett (1965) suggested that the rate of
dilation is a function of specific volume and pressure. Wroth (1990), as reported by Houlsby
(1991), replotted shear tests data obtained by Stroud (1971) and suggested simple expression
for the variation of the angle of dilation with distance from the critical state line:
sin'V = a(r- V,) (7.1)
where
V = specific volume,
7-4
A. = slope of the Critical State Line whose position is fixed by r,
r = value of specific volume at some reference pressure P
a

The value of a was found to be typically near unity.
Recently, Wroth's concept is further elaborated in terms of relative density In rather than
specific volume by Houlsby (1991):
Sill'" = a(r -V...,) + a ( V m ~ - V""nlID - Mln(:J
(7.2)
where
e -e
max
----,
emax-e
min
e = current void ratio,
emax.e
min
= maximum and minimum void ratios.
Expression (7.2) can be simplified as:
where A, Band C are simple functions of well-defined properties of the soil.
Earlier, Bolton (1986) presented the following empirical expression in which he combined the
relative density and pressure into a relative dilatancy index JR'
"'- - "'- = 31
'+'p 'f'cv R
where
The correlation is found to be applicable within the range O<J
R
<4. Bolton's (1986) empirical
7-5
expressions can be combined to derive the following expression for the angle of dilation
(plane strain) in radians:
\jI = -0.11 +0.59[D
(7.3)
Both Bolton's (1986) expression (7.3) and the theoretical expression (7.2) based on Wroth's
concept, given by Houlsby (1991) fit well with the range of experimental data available on
sands. Although, theoretical expression (7.2) is equally convenient and axial capacity of piles
in frictional soils has been determined by Houlsby (1991). The finite element calculations,
presented in this chapter using Matsuoka model for variable friction and dilation angles, make
use of Bolton's empirical expressions to relate dilation angle as a function of pressure and
density. (Similar finite element calculations have been carried by Nutt (1993) to show the
influence of the variable friction and dilation angles on the cone pressuremeter tests in
carbonate sands.)
In order to carry out the finite element calculations based on the variable friction and dilation
angles defining the Matsuoka yield surface, it is necessary first to calculate appropriate
material stiffness matrix and update (Gauss point) stresses during the course of analysis, as
the analysis is usually performed by dividing the applied loads or displacements or both into
a number of steps. New subroutines are added to the finite element code 'OXFEM' to carry
out these tasks without altering the structure of the program.
7.3 Material stiffness matrix
In this section, we will describe how the material stiffness matrices are derived for the
Matsuoka (1976) model. In the following Section 7.4, we can see some important parts of the
stress updating procedures which are adopted in the finite element analyses presented in this
chapter.
7-6
7.3.1 The Matsuoka yield surface constant friction and dilation angle (model 1)
As explained by Burd (1986), the Mohr-Coulomb yield surface contains singularities at the
comers of the yield surface i.e. yield function is not differentiable at the comer points. It is
therefore preferable to choose the Matsuoka surface in which the yield function is
.cr
3
-cr
I
a - Mohr-Coulomb
b Matsuoka
Figure 7.1: Yield surfaces in the' 7t' plane
b
differentiable everywhere. Both these surfaces are shown in Figure 7.1
I
cr
2
The Matsuoka yield function, in terms of stress invariants, can be defined as:
where
77
Chapter 7 Analysis of lilterally loaded piles in frictional soils
<l>' Ie is the triaxial compression friction angle.
A material is said to have yielded when the stresses satisfy the yield criterion:
From the nonnality of the plasticity theory, originally developed for metals, the plastic strain
increments are in the direction normal to the yield surface defined by j so that:
=
dj
A component of the normal to this yieJd surface is given by:
and similarly other components can also be derived.
For granular soil behaviour, however, it is more appropriate to use a separate functiong(O")
tenned as plastic potential is used to describe incremental plastic strain:
dP
= dg
del'
When a dilation angle \jf is specified, the plastic potential 9 can be calculated from the yield
function j. This can be done using one of the two following approaches. In Figure 7.2, the
yield function j is shown in "C : (J space where 't is the shear stress with a constant friction
angle <jJ. Full association occurs when the plastic strain increment dg/dO" coincides with the
nonnal to the yield surface j for a particular stress level. Zero association, which results in
zero volume change, occurs when the plastic strain increment is vertical. Following the
7-8
~ ~ b f
00" -- =0
o f b a'
00"
~ I
,
x
a
Figure 7.2: Yield function and plastic potential in ' 1:-0' ' stress space
approach proposed by Burd (1986), the plastic strain increment can be calculated from linear
interpolation between the full and zero association cases. In an alternative approach, proposed
by Yu (1990), an assumption is made that the plastic potential function 9 is similar to the
yield function f but (ftc is replaced by ViC' Since f and 9 must intersect at the current
stress state, the apex of surface 9 is shifted from the origin by a distance -X so that:
where
and I is the unit vector. This condition is used to derive a cubic equation in terms of the
principal stresses from which X can be solved. (The approach proposed by Yu (1990) is
adopted in the finite element calculations (Section 7.5) involving the Matsuoka yield criterion
with variable friction and dilation angles (Section 7.3.2).)
Once the plastic strain increments are derived, the plastic stiffness matrix [DjP defined in
7-9
equation (4.10) can be calculated. The stress changes can then be related to the strain changes
by an elastic-plastic stiffness matrix defined as:
[DjEP = [DjE
IDlE [ ~ ] [ # r IDlE
[ddH IDlE [ ; ~ - ]
where [DJE is a symmetric matrix of elastic constants: bulk modulus K and shear modulus
G.
K+4G
K_2G K_2G
01
3 3 3
K_2G
K+4G
K_2G
0
[DJE = 3 3 3
K_2G K_2G
K+4G
0
3 3 3
0 0 0 G
The friction angle Q- Ie and dilation angle \Vte remain constant throughout the analysis and
hence the yield surface continues to extend indefinitely with increasing mean effective stress.
7.3.1.1 Finite element analysis of a single-pile
In the finite element analyses described in this section, a rigid disc approximation of a single-
pile is adopted as it has been done earlier in Chapters 5 and 6. The lateral response of a rigid
disc in a frictional soil is calculated using displacement controlled finite element analyses. A
7-10
constant friction angle (I>' tc of 32 with dilation angle 'V of 0.5
0
10 , 15
0
is used in the
analyses. The soil has a Poisson's ratio v of 0.3 and GI Po" ratio of 300 (where Po" is an
initial compressive stress).
Since the problem is symmetric, it is only necessary to analyse one half of the disc. Six-noded
linear strain triangular elements are used. A complete rough contact between the pile-soil is
assumed. The finite element mesh is shown in Figure 7.3. The outer boundaries of the mesh
are fixed at a distance of 20r
o
from the centre of the disc; where ro is the radius of the disc.
The fixed boundaries are thought to be far enough to have no detrimental effect on the load
displacement response.
Using the material model described in Section 7.3.1, the lateral response of a single disc in
frictional soils is calculated and the results are shown in Figures 7.4 and 7.5.
L<jL
Figure 7.3: The finite element mesh for a rigid circular disc
The results in Figures 7.4 and 7.5 involve the calculation of plastic potential as proposed by
Burd (1986) and Yu (1990) respectively. Burd's (1986) approach considers dilation as a
degree of association and is denoted by ''Y'. (The material properties considered are the same
7-11
150 r - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ~
,.-....
Z
C
"'0
co:j
o
100
- so
0.005
-.-y=o
_ y= 0.4
0.01 O.Ql5
displacement (m)
0.02
y= 0.2
y= 0.6
0.025
Figure 7.4: Lateral response of a single-pile in sand
140 .---------------------------------------,
120
100
~ 60
o
.-
40
20
0.005 0,01 0.015 0.Q2
llf= 0
0
displacement (m) 0
__ 't' -.- \j! - 5
o 0
-6- \j!= 10 - \j!= 15
Figure 7.S: Lateral response of a single-pile in sand
0.Q25
as mentioned above except that instead of using dilation angles, an equivalent degree of
association y of 0, 0.2 0.4 0.6 is used in the analyses.) It can be seen that using a non-
associated flow rule, the ultimate lateral load is obtained only when the dilation angle is zero.
However, Lane and Griffiths (1988) have reported 'mesh locking' for this case using the
Mohr-Coulomb model. In case of dilation, plastic volumetric strains continue to increase and
7-12
Chapter 7 Analysis of io.leraUy loaded piles in frictional soils
it appears that the ultimate lateral load cannot be achieved.
7.3.2 The Matsuoka model - variable friction and dilation angles (model 2)
It is unrealistic to assume that the peak friction angle remains constant in frictional soils at
increasingly high normal stresses, or that volumetric strains continue indefinitely. As was
mentioned earlier (see Section 7.2) dense frictional materials have higher rates of dilation but
for a given density, the angle of dilation reduces with increasing stress levels. The
relationships between (j>' p - <1>' cu' IR and 'Vrnax as defined earlier in Section 7.2, are used here
with the Matsuoka yield function (Section 7.3.1) to model a more realistic yield surface as
shown in Figure 7.6. The failure surface remains unbounded, and any increase beyond a
critical stress level, the friction angle remains constant at <1>' cu and the dilation is zero.
Variations in <1>' and 'V are accounted for in the model by specifying <1>' cu' V
max
' V min and
initial I
D
, denoted by I
D
o
where ID is continually updated during the calculation. An
Figure 7.6: Curved yield surface in ' "C-O' ' stress space
7-13
Chapter 7 A nalysis of laterally loaded piles in frictional soils
additional tenn is also necessary in detennining the nonnal to the failure surface aJ/acr'
given above because of dependency of <1>' on cr'. This tenn is given by:
=
where
48 IiDtan<Jl ,
I} cos
2
<Jl'
I D = _V_max_-...:;..[(_V m_aJC_-=-=I D_O(_V-::-:max;---_V_ml_)_) (l_+_tJ-..c.. J
V -V.
max rrun
tv are the current volumetric strains.
7.4 Stress update calculations
An important feature of the finite element program 'OXFEM' is the numerical solution
scheme (Modified Euler scheme (Burd (1986. It requires an equilibrium check at the end
of each load step. The equilibrium is maintained by making sure that the work done by the
internal stresses is equal to the applied load which is done by evaluating an integral of the
fonn:
where [Bl is obtained by appropriate differentiation of the shape functions of the element.
This feature of the program 'OXFEM' requires accurate stress updating at the end of each
load step. As described by Burd (1986), three separate operations, in general are needed to
7-14
Chapter 7 Analysis of laterally loaded piles in frictronal soils
calculate the updated Gauss point stresses. Firstly, it is necessary to calculate the strain ratest
which when integrated over the time step, are compatible with the incremental displacements.
The second stage is to obtain the incremental Gauss point stresses by integrating the
constitutive equations. In some cases, updated stresses may be obtained from a closed form
solution (e.g. elastic behaviour, and von Mises plasticity excluding the Jaumann terms) but
for the Matsuoka model it is necessary to use a numerical integration scheme. An 'error
control' algorithm is used, based on a Runge-Kutta scheme in which the updated stresses are
evaluated by dividing the calculation step into a number of 'sub-increments'. In order to
improve the efficiency of this calculation, a procedure is used in which the size of the 'sub-
increments' are varied automatically in order to perform the integration at each of the Gauss
points to a specified accuracy. In this algorithm, which is described by Sloan (1984), an
estimate is made of the errors occurring at each' sub-increment'. The magnitude of the error
is then used to determine the size of the following 'sub-increment', or, if the error is
unacceptable, the size of the sub-increment is reduced, and the calculation repeated. The third
and the final stage of the 'stress update' calculation is to correct the stresses back to the yield
surface, as a numerical scheme has been used to integrate the constitutive law.
7.4.1 Calculation of yield surface intersection
As we have seen in Chapter 4 that if the material becomes plastic during the calculation step,
it is necessary to find the point in the increment at which the yield surface intersection occurs
(see Figure 7.7). The increment is then split into two parts and the 'stress update' calculation
for the elastic and plastic parts performed independently.
For the case of Matsuoka yield surface defined by a constant friction angle, a cubic equation
of the form:

is derived (see Burd (1986)) and solved for the appropriate root. The coefficients a
l
, and
7-15
Chapter 7
Figure 7.7: Yield surface intersection
Analysis of Interally loaded piles in frictional soils
1
(j
i+l
(j
involve the calculation of the constant ~ defined earlier. For the case of Matsuoka model with
variable friction angle, the value of ~ varies with the angle of friction. The angle of friction
varies with time and is a function of logarithm of hydrostatic pressure which is also a
function of time.
where
The coefficients c
1
' c
2
and c
3
are defined in Burd (1986). As we can see that the cubic
equation initially derived by Burd (1986) to calculate yield surface intersection, now takes a
much more complex form due to the dependence of the angle of friction on pressure and
relative density. Therefore, we need to resort to numerical methods to calculate yield surface
intersection. An interpolation method (false position method) is used to calculate the
intersection. The rest of the stress update procedures i.e. integration of the constitutive
equations and the yield surface correction are the same as described by Burd (1986).
7-16
7.5 Verification of the Matsuoka model (2)
7.5.1 Plane strain analysis of a unit square soil element
In this test problem a square element is given initial isotropic compressive stresses and is
sheared by imposing additional tensile forces in the x-direction. The boundary conditions and
the applied loading on the element are depicted in Figure 7.8. The material properties used
(a) (b)
Figure 7.8: Unit square soil element
in the analysis are: G = 100OOkPa; Po' = 800kPa; Poisson's ratio v = 0.3; IDa = 0.9; V
max
= 2.79 and V min = 2.49. The validation of the model is carried out by plotting out variation
of the friction and dilation angles during the analysis (see Figures 7.9-10).
As we have seen in the formulation of the Matsuoka model, the yield surface is defined by
the triaxial friction angle. Therefore, in order to carry out the plane strain analyses using this
model, we need to find out the equivalent plane strain friction angle. Burd (1986) presented
such correlations for the zero and full association of the plastic potential. We have used the
correlation for the zero association (i.e. no dilation) to see that the plane strain friction angle
7-17
50 .-------------------------.,
'W
0
c::
0
~ O
u
:.s
I\)
20
en
c::
~
10
( ) ~ - - - ~ - - - - ~ - - - - - ~ - - - - ~ - - - ~
o (),I (),2 (J,} 0.4 (),S
displacement (m)
Figure 7.9: Variation of friction angle
calculated from the correlation, matches with the critical state friction angle at the end of the
analysis. A triaxial critical state friction angle (1)' cv equal to 3l.04 is used in the analysis.
The corresponding plane strain critical state friction angle is worked out to be of 34. The
friction angle at large displacement, shown in Figure 7.9, matches exactly with this value.
()12
01
0,08 -
'" c
.;:;
b
001> -
V>
u
'E
OJ
).j
E
;:l
'0
>
002 -
!l
0.02
0 0,1 0.2 0.3
strains in x-direClioll
Figure 7.10: Variation of volumetric strains
7-18
--
04 0,5
Figure 7.10 shows the variation of volumetric strain tu with the strain in the x-direction
obtained from a force controlled finite element analysis. We can see from the figure that the
square element after having an initial compression tends to dilate and as the stresses increase
and the relative density decreases during the analysis by the application of further tensile
force, dilation gradually disappears.
7.S.2 Laterally loaded single-pile (after Lane and Griffiths (1988)
In this section, the analysis for a laterally loaded rigid disc analysis is carried out using the
Matsuoka model (2). The finite element mesh used in the analysis is the same as described
in Section 7.3.1.1. A similar analysis, using Mohr-Coulomb model, was carried out by Lane
and Griffiths (1988). In their analysis, the soil had a friction angle <l>' of 30, cohesion c
of lOkPa and Poisson's ratio of 0.3. The load displacement response is shown in Figure 7.1l.
The figure shows that the ultimate load value is not well-defined. The load displacement
curve obtained, using the Matsuoka model (2), is also shown in Figure 7.11. The following
35
30
25
'"0 20
'* 0
-e:
~ 15
o
--<
10
5
o ~ - - - - ~ - - - - ~ - - - - ~ - - - - ~ - - - - - - ~ - - - - ~ - - - - ~ - - - - ~
o 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8
old
_Lane and Griffiths (1988) -e- The Matsuoka model (2)
Figure 7.11: Comparison of load displacement response
7-19
material parameters are used in this analysis: the initial relative density I
D
o
of 0.3,G/po
ratio of 1000, the initial and final void ratios of 1.49 and 1.79 respectively. The model not
only traces back exactly Lane and Griffiths's (1988) curve but also gives a value for the
ultimate lateral load.
7.6 Laterally loaded single-pile - the effect of initial relative density IDo
In this section the analyses described in Section 7.3.1.1 are repeated using the Matsuoka
model (2). The finite element mesh and the material properties are the same as given in
Section 7.3.1.1. The analyses are carried out for values of initial relative density 10 = 0.33
and 0.9 representing relatively loose and dense soils. The results are plotted in Figures 7.12-
13. It is shown in Figure 7.12 that loose sands do not dilate and the Matsuoka model (2) can
reproduce the obtained curve using the Matsuoka model (1) with zero dilation. However, if
the soil is initially too loose then it may cause the friction angle at first to drop down below <l>' cv
140
120
100
,.-....
g
80
Matsuoka model (2)
'"0
60
C':S
0
-
40
20
0.02 0.04 0.06 0.08 0.1 0.12 0.14
Matsuoka model (I) displacement (m)
__ Psi = Odeg .......- Psi = 5deg
-.t.-Psi = 10deg _Psi = 15deg
Figure 7.12: Load displacement response (l DO = 0.33)
7-20
and then recovering back to (t> cv during the course of the analysis. On the contrary, dense
sands dilate and the Matsuoka model (2) is capable of predicting the ultimate lateral load in
a dilating sand using non-associated flow rule (see Figure 7.13), However, in the absence of
an analytical solution, it is difficult to validate these results.
7.7 Conclusions
The Matsuoka model (2) is a more realistic model to predict the response of laterally loaded
piles in sands. The finite element analysis can be used to predict the load displacement
response (a theoretical p-y curve). The theoretical p-y curve captures the non-linear soil
11() r--------------------------...
\1at(\Joka modd (2)
HXI
~
Z RO
~
-
-0
,.,
t;O
..Q
~ o
20
o ~ - - ~ ~ - - - ~ - - - ~ - - - - - ~ - - - - ~ ____ L_ ___ ~
()
11.112 1I1I1 () 06 o.og !11 0.12 n l ~
\1nI'uoka Jnodel (I)
displacement (m)
___ Psi Odt'g -+- Psi = )deg -.- Psi = IOdeg -0- Psi I)deg
Figure 7.13: Load displacement response (l DO = 0.9)
behaviour and seems more reliable than the empirically constructed p-y curve proposed by
Reese et al. (1974) which is essentially a bilinear curve based on a constant friction angle
(Yan and Byrne (1992). The theoretical p-y curve can be non-dimensionalised as shown in
Figure 7.14 and a series of p-y curves along the length of the pile can be obtained by a single
finite element analysis. These curves can then be used as inputs to the standard subgrade
7-21
Chapter 7
'\:l
,* c
0..
Analysis of laterally loaded piles in frictional soils
120 r - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ~
100
80
0z:. 0 60
~
ro
...s
40
20
O ~ - - ~ - - - - ~ - - - - ~ - - - - ~ - - - - ~ - - ~ ~ - - ~ - - - - ~
o 100 200 300 400 500 600 700 800
G*6/p!d
Figure 7.14: A theoretical p-y curve
reaction analysis (see Chapters 5-6) to predict the pile head displacement of a single-pile as
well as pile groups. A further work similar to that carried out in Chapters 5-6 is accordingly,
proposed in this direction.
7-22
Chapter 8 Bearing capacity of a jackup spudcan footing
CHAPTER 8
BEARING CAPACITY OF A JACKUP SPUDCAN FOOTING
8.1 Introduction
The ultimate bearing capacity of a flat surface footing, as briefly discussed in Chapter 4, is
dependent on the footing-soil interaction under vertical, horizontal and moment loading which
governs the stability of the footing. In this Chapter we shall consider the stability of a jackup
spudcan footing and calculate its ultimate bearing capacity using the finite element method,
though without pursuing the detailed interface analysis undertaken for a flat footing.
The growing use of mobile jackup units on spudcan footings in deep offshore waters has
raised a great deal of concern about the overall stability of a jackup unit in hostile
environmental forces. The stability of a jackup unit is greatly influenced by the performance
of a typical spudcan footing under storm loading conditions. Therefore, as a result, much
research has been focused towards the study of a spudcan footing behaviour subjected to
combined vertical, horizontal and moment loading. Recently, over a period of several years,
a Joint Industry Study coordinated by Noble Denton Associates (NDA) has been considering
many aspects of jackup stability including the soil-structure interaction (see Osborne, Trickey,
Houlsby and James (1991)).
Mobile jackup units consist of a floatable drilling platform, supported on three or more legs
which can be raised or lowered (see Figure 8.1). A detailed description of the installation
procedure is given by Tan (1990). The platform legs can either be supported separately as
shown in Figure 8.1a or they can be supported on a single shared mat (Arnesen et al. (1988)).
8-1
Chapter 8
Bearing capacity of a jackup spudcan footing
WIND--........
~ - - PLATFORM
WAVE-----...
~ - - - - t - - - - - - - LEG
CURRENT --____
SPUD CAN
(a)
o
(b)
Figure 8.1: (a) Jack-up platform and environmental loads (after Tan (1990))
(b) Spudcan footing
8-2
Most modem jackup platforms are of the fonner type and have approximately conical shaped
footings with a protruding tip at the centre. These are commonly referred to as 'spudcan'
footings (see Figure 8.1 b).
The mobility of jackup units allows for relocation and re-use of the platform. Before actual
oil and gas exploration activity begins, the legs of the platform are preloaded to almost twice
the working load. During this loading operation, spudcan footings move into the seabed to
a depth of almost twice their diameter. In soft clays, a 30m penetration of a 15m diameter
spudcan is not uncommon. The bearing capacity of spudcan footings increases with the depth
of penetration below the seabed. The increase of bearing capacity with footing embedment
is an important factor, since even a small embedment can significantly increase the bearing
capacity of such offshore footings CHoulsby and Martin (1993. However, this is not usually
the case in coarse materials as penetration depths are typically very small (see e.g. Dean et
al. (1993.
In addition to the vertical load of the jackup structure, a spudcan footing is also subjected to
large horizontal forces due to severe environmental conditions reaching, in stormy weather,
up to 30% of the total vertical load and thereby producing a correspondingly large moment
(Poulos (1988. However, moment loading will not be pursued further here. Environmental
loads, as mentioned elsewhere in this thesis (see Chapter 1) are cyclic in nature and are
beyond the scope of this thesis.
In this chapter we are making a conventional assumption, that a jackup leg is pinned at the
foundation level and that therefore a spudcan footing offers no moment restraint to the leg
(see Figure 8.2a). This means that, only the interaction between horizontal and vertical
loading needs to be studied. Although, in recent studies moment restraint at the footing has
also been considered which would reduce the bending moment in the lower leg guide (see
Figure 8.2b), which is often a critical feature of the structure in extreme loading conditions:
a full interaction between vertical, horizontal and moment loading has been studied by e.g.
Houlsby and Martin (1993) for clays, and by Dean et al. (1993) for sands.
We discuss here the ultimate bearing capacity of a partially penetrated spudcan footing
8-3
Chapter 8
WAVE FORCE A
8
WAVE FORCE
10) NO FIXITY
SEABED
~ v -
H ~ " - r M
V
I 0) PART IAL FIX ITY
Figure 8.2: Bending moments in the legs of Jack-up units (after de Santa Maria
(1988))
8-4
subjected to combined vertical and horizontal loading, and describe a three-dimensional finite
element analysis of the footing. (Although, in the preceding chapters we preclude the use of
a complete three-dimensional analysis for deep foundations due to the enormous amount of
computational effort involved. The same is true here for a surface footing but the size of the
problem is not as large as a multi-pile group foundation.) The aim of the analysis was to carry
out an independent investigation of the behaviour of the spudcan footing under combined
vertical and horizontal loading, and compare the results with the existing theories. A parallel
laboratory investigation of the same geometry is described by Hambly (1992).
8.2 Background
The bearing capacity solutions of Meyerhof (1953), Hansen (1970) and Vesic (1975) are
commonly used to determine the ultimate bearing capacity of plane strain footings. These
solutions are modified by introducing shape factors to cover circular geometries (see e.g.
Hambly (1992), Dean et al. (1993. In the case of spudcan footings, their embedded circular
area in plan (i.e. plan area at ground surface) is used for the bearing capacity calculations. In
order to investigate the applicability of the bearing capacity solutions to spudcan footings,
theoretical and experimental studies of conical shaped footings (based on centrifuge model
tests) have been undertaken by Cambridge University Engineering Department (CUED) and
Andrew N. Schofield and Associates (ANS&A) over a period of several years. Work has been
reported in various publications including M.Phil and Ph.D thesis. One such extensive study,
carried out by Tan (1990) is reported in Dean et al. (1993). In this study, a series of tests of
cones of angles 0 (flat), 13. 25, 45 and 60
0
was conducted on sand at 28.3mm
diameter model footings at 56.6g, representing 1.6m prototypes in the drum centrifuge. It is
shown in the study that spudcan footings can be treated as equivalent cones enclosing the
same volume. The load displacement response of a spudcan footing and an equal volume cone
is almost identical. The ultimate bearing capacity of these footings for vertical loading only,
combined horizontal - vertical (inclined) loading, vertical - moment (eccentric) loading and
8-5
lateral - moment loading has been determined. The results are compared with existing bearing
capacity theories and analytical solutions where possible. The study shows that the vertical
bearing capacity does not vary much for cone angles between 0
0
-20
0
However, for flat
circular footings under combined vertical horizontal loading (referred as shear sideswipe
tests at constant penetration), most of the experimental data lies outside the CUIVes proposed
by Meyerhof (1953), Hansen (1970) and Dean et al. (1993).
Hambly (1992) carried out laboratory tests in sands on two model footings - flat and spudcan.
The purpose of this investigation was to compare the ultimate bearing capacity of flat and
spudcan footings. In the series of tests perfonned by Hambly (1992), the footings were
preloaded until the penetration for the flat footing was 10% of its diameter (Vesic (1975 and
for the spudcan footing, a penetration of its bearing area at the ground surface equivalent to
the bearing area of the flat footing. Hambly (1992) found that the preload force that mobilized
the ultimate bearing capacity of the sand under the flat circular footing was insufficient to
mobilize the ultimate bearing capacity of the sand under the spudcan footing. Hambly (1992)
showed that the ultimate bearing capacity of a partially penetrated spudcan footing is twice
the preload force required for the flat footing. He obtained enhanced sliding resistance of the
spudcan footing using twice the preload force for the flat footing with Hansen's (1970)
theory. He attributed the difference between the preload force and the ultimate bearing
capacity to the much greater volume of sand displaced by a flat footing mobilising ultimate
bearing resistance, as compared to the volume displaced by the partially penetrated cone of
the spudcan. He suggested 'partial penetration factors' to determine the ultimate bearing
capacity from the vertical preload force. He concluded in his study that a further laboratory
investigation is needed to validate his results. Although, in the absence of such an
investigation it is difficult to comprehend Hambly's (1992) results, the experimental results
described by Tan (1990) resembles more closely the actual preload operation of the spudcan
footing. The footing penetrates into the seabed as the load is applied and continues to
penetrate until there is no further penetration i.e. the footing is locked vertically. The final
value of the load, thus, represents the ultimate bearing capacity of the footing.
8-6
The three-dimensional finite element analysis of the spudcan footing here described was done
using the commercial finite element code' ABAQUS'. For the given penetration depth of the
spudcan footing (i.e. the depth at which the ground surface bearing area of the spudcan
footing became equal to the bearing area of the flat footing), the finite element analysis
predicts a vertical preload which is the ultimate bearing capacity of the partially penetrated
spudcan footing i.e. theoretically equal to the ultimate bearing capacity of the flat footing of
the equivalent bearing area. The numerical results obtained from the finite element analysis
(which can also be found elsewhere: B urd (1991 compares well with Hambly's (1992)
laboratory results for both the flat footing and the spudcan footing with a partial penetration
factor of 2. A comparison with the results from Hansen's (1970) solution showed the latter
(as mentioned earlier: see Section 4.4.3) to give rather conservative estimates of the bearing
capacity of such offshore footings.
Figure 8.3: Spudcan geometry
8.3 Description of the problem
i
i
i
i
i ~
L22m
:... "'1
3.35m
2.13111
The dimensions of the spudcan footing (the 'cone' and its protruding central tip) assumed in
the finite element calculations are shown in Figure 8.3. The face of the 'cone' is inclined at
8-7
an angle of 13.46
0
to the horizontal (later in the text referred as cone angle' n') with an apex
angle of 152.7. The dimensions of the central tip are also shown in Figure 8.3. The spudcan
is partially embedded (depth=2.13m) in dense sand with a bearing area diameter equal to
6.71m. The following load paths are analysed:
1) Initt.a1 vertical loading to determine the preload, for the given
penetration depth.
2) Reduction of vertical load with increase in horizontal
load to cause yielding in dense sand
3) A further reduction of vertical load followed by further increase in
horizontal and vertical loads.
In each of these load paths, as indicated above, moment loading is not taken into account.
The finite element analyses were carried out using a linear elastic, perfectly plastic soil model
in which failure is defined by the Drucker-Prager (1952) failure criterion. The failure surface
is circular when viewed in the 1t plane, as shown in Figure 8.4a. The Drucker-Prager model
is an inaccurate theoretical representation of granular material behaviour. The failure surface
can either be inner matching or outer matching to the Mohr-Coulomb irregular hexagon
shown in Figure 8.4a. In Figure 8Ab, the Mohr-Coulomb and the Drucker-Prager yield
criterion are plotted in a stress (q-p) space. The figure shows that for an inner matching circle,
the Drucker-Prager yield line matches well with the Mohr-Coulomb yield line for that part
of the stress space which represents the stress states as they are in a triaxial extension test (Le.
in the fourth quadrant) . Similarly, for an outer matching circle, the Drucker-Prager yield line
8-8
A
triaxial compressiro
o
()
A Drud:.aPrtge:r circle inner matching to MohrCoulomb he)(jgon C, eXLended [)ruckel-Purer circle
a B _ Mohr-Coulomb t.e"gon (b) D - Drucker-Pra,er <,,,,Ie outt:r mltcrun, to t.eulon
Figure 8.4: (a) Yield surfaces in '7t' plane (after ABAQUS manual)
(b) Yield surfaces in (q-p) stress space
matches with the Mohr-Coulomb yield line for the stress states as they are in a triaxial
compression test (i.e. in the first quadrant). The Drucker-Prager model does not yield a unique
match for both the Mohr-Coulomb yield lines in the stress space. It was difficult to choose
parameters for the Drucker-Prager model to depict a given Mohr-Coulomb friction angle, as
the parameter that determines the inner matching circle to a given Mohr-Coulomb hexagon
does not yield the same outer matching Mohr-Coulomb hexagon that one would expect from
an inner matching Drucker-Prager circle, i.e. it does not yield the same (chosen) Mohr-
Coulomb friction angle, and the difference between these two Mohr-Coulomb angles is quite
substantial (see Figure 8.Sa). It seems to suggest that if the failure condition is described by
triaxial extension tests, the consequence is that the angle of friction in compression tests is
smaller than the angle of friction in extension tests (see Atkinson and Bransby (1978.
Similarly, with the outer matching Drucker-Prager circle to a given Mohr-Coulomb hexagon.
the parameter leads to a non-existent inner matching hexagon for the chosen outer matching
Drucker-Prager circle, which cannot be right. (This situation is depicted in Figure 8.Sb.) We
need to have a plasticity model which provides a smooth approximation to the Mohr-Coulomb
surface. The Matsuoka (1976) model shown in Figure 8.5c. by contrast. more accurately
8-9
Chapter 8
(J
3
(l
1
A . hn.!tgnn
n.- cin::k mne, tl!8.!,-hirJ! In'!\'
A .. 'mh h .. 'tII J!11n "Iller Itll'lh hlng In '1\'
(a)
Figure 8.5: Yield surfaces in 'It' plane
Bearing capacity oj a jackup spudcan footing
u
1
A- Moht-CouJomb
B
(c)
(l
I
A - MllhrTownmh hC'J.3gnn
,\ nOB c'(iSillll!
(h)
represents granular material behaviour and resembles the shape of the Mohr-Coulomb
hexagon in the 1t plane. This model, however, is not available in 'ABAQUS'. It was
therefore necessary to determine appropriate parameters for the Drucker-Prager model to
represent an equivalent Mohr-Coulomb friction angle for the soil. In all of the finite element
calculations, the results of which follow, the parameters of the failure surface were chosen
in such a way that the surface coincides with the triaxial extension points (inner matching
circle) of the Mohr-Coulomb envelope for a Mohr-Coulomb friction angle of 37.5. That
8-10
surface corresponds to Curve (B) on Figure 8.5. An inner matching Drucker-Prager circle was
deliberately chosen since it was more likely that the stresses inside the soil underneath the
spudcan footing would reach a stress state which is better represented by a triaxial extension
test, i.e. the horizontal load is dominant to the applied vertical load on the spudcan footing.
We will see later in this chapter that this choice (inner matching Drucker-Prager circle) turned
out to be quite successful in predicting the failure load for the spudcan footing under
consideration. The predicted failure load is close to the one estimated analytically. The
plasticity model is fully associated, i.e. the angle of dilation is set as equal to the angle of
friction. In order to avoid numerical problems in the analysis procedures, it is necessary to
include a small amount of cohesion in the soil model: a value of 2.5 kPa is adopted in all of
the calculations here.
8.5 Finite element model
As we have already seen in the preceding Chapters 3,5-7, it is difficult to model rigid
boundaries in a finite element analysis. The presence of rigid boundaries in close proximity
to the loaded area has a detrimental effect on the quality of the resulting solution. Moreover,
the sheer volume of computations involved in a three dimensional analysis is too costly in
computer time. It is therefore necessary, in order to contain the bulk of the problem, to
nominate an appropriate location for the rigid boundaries in the finite element model.
Theoretical calculations of bearing capacity were carried out for this purpose by Houlsby
(1991), using an axi-symmetric method of characteristics analysis, described in Houlsby
(1982), Houlsby and Wroth (1982, 1983). This method gives a lower bound solution to the
incremental collapse problem. (Various different stress fields (characteristics) can be studied
which satisfy the conditions for a lower bound solution but the only stress field that gives a
true lower bound is the one which can be extended throughout the whole soil mass.) The
theoretical extent of the failure zone was determined and the rigid boundaries in the finite
element model were then set sufficiently far beyond this zone. In the axi-symmetric analysis
8-11
done here, the footing is modelled as a simple cone, assumed to be fully rough, and a Mohr-
Coulomb friction angle <l> of 37.5 is assumed for the soil. The mesh of characteristics
obtained from this analysis is shown in Figure 8.6. It indicates that the zone of failure, for
3.35m
rigid, rough footing
--;----,.---,.---.,.... .. ,/
Figure 8.6: Characteristics mesh for vertical loading (after Houlsby (1992))
the case of vertical loading (load path 1), extends to a depth of about 5m and to a radius of
about 16m. The rigid boundaries in the finite element mesh were therefore set at a radius of
20m from the centre of the footing and 7m below the ground surface. The linear dimensions,
at which the boundaries are set, are, respectively, 25% and 40% greater than the
corresponding dimensions of the theoretical failure zone. The elastic parameters used in the
model have a Young's modulus of 500 MPa and a Poisson's ratio of 0.2.
8.6 Finite element analyses
In this section, the finite element analyses for the load paths (mentioned in Section 8.3) are
described. Load path 1 determines the vertical bearing capacity of a spudcan footing in the
8-12
absence of any horizontal loading. Load paths 2 and 3 are complementary and detennine the
ultimate bearing capacity of a spudcan footing subjected to combined vertical and horizontal
loading.
8.6.1 Vertical bearing capacity of a spudcan footing (vertical preload to cause yield)
In this section, the vertical bearing capacity of a spudcan footing is detennined using the axi-
symmetric and the three-dimensional finite element analyses and the results are compared
with Hansen's (1970) bearing capacity solution. We describe first Hansen's solution and then
the finite element analyses.
8.6.1.1 Vertical bearing capacity solution
The bearing capacity solution for a circular footing, given by Hansen (1970) can be simplified
(for vertical loading only), and written as presented by Dean et al. (1993):
V=A .(J:...N.y,.B)
p 2 y p
(8.1)
where Bp is the penetrated diameter, Ap is the plan contact area, y , is the soil effective
(buoyant) unit weight, and Ny is an axi-symmetric self-weight bearing capacity factor. Bp
can be written in tenns of penetration depth D pas:
B =2D cotO
p p
where Q is the cone angle (see Figure 8.7a).
Substituting B in the bearing capacity solution (equation 8.1), we get the following equation:
p
8-13
Chapter 8
>
V
B
(a) Conical footing at partial penetration p

t3 100
cz::
o
t....
60 DECREE CONE (FHL)
(e:::Q.8., 8-28.3 mm. g-56.6)
FULL PENETRATION UN(--.......
(b) Initial cubic nature of curve before full penetration (after Tan (1990

\00
Z
........
60 CEGREE. CONE (FT -4L)
(1:",0.64, 8=28.3 mm. <;;:=56.6)
> t3
cz::
o
Yo.
50
Figure 8.7: (c)
2
4 6 6 10 12
(Thousonds)
DISPLACEMENT 1 (mm')
Plot of force against the cube of displacement
(after Tan (1990)
8-14
Chapter 8
/
/
20 40 EO
'"
C:ore An<;!e Q
,.
, .'
. .-
- - I .. -

o 20 40 60
(:,
Cone Q
Figure 8.8a: Theoretical axisymmetric Ny values for cones on cohesionless soil (after Tal
(1990
50
SEMI-ROUGH FOOTINGS (TEST SERIES m)
(8=28.3 mm, 9=56.6)
40
; =.;H. a '" 19.3"
30
...
l'-
; =31",
z
20

10 ... MEDIUM DENSE (e=0.80)
LOOSE (e.:0.95)
0
0 20 40 60
CONE ANGLE n
Figure RRb: Errect uf ("one on Ny of scmi-nmgh footings; thcoretical and
experimcntal compa.-ison (after Tan, 1990)
8-15
Thus for a given cone angle Q, V is directly related to the D /' if Ny is approximately
constant then V and D 3 relationship is a straight line. Such a straight line relationship is
p
shown in Figure 8.7c for a 60 cone (Tan 1990) on fine Leighton Buzzard sand tested in the
drum centrifuge at 56.6g. However, for a given depth of penetration, the value of Ny varies
with the cone angle and the variation depends on the surface roughness of the footing. The
theoretical variation of Ny with the cone angle using the method of characteristics is shown
in Figures 8.8a . Figure 8.8b shows the same variation of Ny with the cone angle calculated
from the centrifuge test results using equation (8.1) and the comparison with the theoretical
results. The results agree well in the region <1>=31 and 8=17 for loose sand (void ratio 0.95)
and <1>=34 and 8=19 for medium dense sand (void ratio 0.8). (In these figures <l> is the
friction angle and 8 represents the degree of roughness of the cone.) It can be seen from the
Figures 8.8a that for a fully rough footing (as is the case here in this chapter) the value of Ny
stays almost the same for cone angles between 0 -20 . But for cone angles greater than 25
there is a significant increase in the value of Ny. These results show that a conical footing
with cone angle Q=13.46 has a slightly higher ultimate load compared to a footing with
Q=00 (flat). But there is a need for further experimental investigation in order to prove that
the load is 100% higher as we learn from Hambly (1992).
8.6.1.2 Finite element analyses results (load path 1)
A preliminary analysis was carried out using the axi-symmetric mesh shown in Figure 8.9.
The overall dimensions of the mesh were set at 20m radius and 7m depth. The mesh consists
of a total of 21 'reduced-integration' 8-noded quadrilateral elements. The footing is taken to
be perfectly rough. The analysis consisted of two stages. In the first stage, the appropriate self
weight stresses were applied to the model, using a soil unit weight of 10.2 kPa/m (submerged)
8-16
(i. 3.351ll
- ~ - -
----,-
7m
-- ._'-
20m
Figure 8.9: Axi-symmetric finite element mesh
40 r------------------------------------------------------------,
-
..... -
3.5
)0 -
25
"'Cl
C':l
.Q
20
0.0
.5
.....
0
0
~
15 -
10
5
o
- 0.01
0.02 _ O.OJ 0.04 0.05
Vertical displacement (m)
Figure 8.10: Vertical load-displacement curve (axi-symmetric analysis)
and a Ko value of 1.0. In the second stage, a set of vertical prescribed displacements were
8-17
applied to the footing.
A plot of the footing load-displacement curve obtained from this analysis is given in Figure
8.10. The footing load at the end of the calculation is 38 MN, with the vertical displacements
reaching 47mm.
A further analysis of the vertical failure load was studied using the three-dimensional mesh
shown in Figure 8.11. This mesh has roughly the same dimensions as the axi-symmetric mesh
shown in Figure 8.9 and consists of a total of 120, 20-noded 'reduced-integration' hexahedral
elements with a total of 719 nodes.
r
r
Figure 8.11: Three-dimensional finite element mesh
In this mesh it proved necessary to flatten the apex on the axis of symmetry of the footing
in order to allow the mesh to be built entirely from the hexahedral elements. (An earlier
analysis based on the unmodified footing geometry and a mesh consisting of a mixture of
hexahedral elements and wedge elements had been unsuccessful.)
When viewed from above, the mesh is seen to be built up of a set of concentric half octagons.
The footing is therefore, implicitly, assumed to be octagonal in plan. The geometry of the
mesh is set up so that the actual vertical area of the footing is equal to the area of the
8-18
Chapter 8
50
40 -
~ 30
"d
cj
.s
b.O
. 20
o
~
10
---.-.. ----
o ~ ________ ~ ~ ________ -L __________ -L __________ ~ ________ ~
o 0.01 002 0.03
0.04 0.05
Displacement (m)
Figure 8.12: Vertical load-displacement curve (three dimensional analysis)
corresponding octagon.
The load-displacement curve obtained from this analysis is shown in Figure 8.12. It indicates
a final load of 41 MN, with the vertical displacements reaching 40mm.
The value of bearing capacity factor Ny for the current problem is worked out to be 52 using
equation 8.1. The value corresponds to a friction angle of 36 which is quite close to the
friction angle of 37.5 , used in the current finite element analyses. The discrepancy is below
5% which may be a result of the choice of soil model which is unable to model sand
behaviour accurately (see Section 8.4). The finite element analyses seem to confirm validity
of Hansen's (1970) equation to calculate the vertical bearing capacity of spudcan footings.
8-19
8.6.2 Reduction of vertical load (vertical unloading) with increase of horizontal load to
cause yield
The bearing capacity solutions and the finite element calculations described in this and the
following sub-section are regarding the bearing capacity of the spudcan footing subjected to
combined vertical - horizontal loading.
8.6.2.1 Bearing capacity solutions
Various empirical equations to calculate bearing capacity include equations proposed by
Meyerhof (1953) and Hansen (1970) for plane strain surface footings which are modified for
circular footings by the introduction of a shape factor. These equations are equally applicable
to spudcan footings since area of the footing in plan is used in the bearing capacity
calculations.
Meyerhof's equation: (8.2)
Hansen's equation:
[
(
J
'25] H = V 1- V
V
M
V
M
V
M
(8.3)
It is possible to obtain an interaction diagram in terms of (H/V M) and (V/V M)' These
interaction loci are shown in Figure 8.13. Also shown is an interaction locus based upon the
empirical equation given by Dean et at. (1993):
8-20
H =(a { V 11- V J
V
M
13 JV
M
V
M
where the constant (alP) is taken as 0.56. It can be seen that the Meyerhof and Hansen loci
are in general conservative in relation to available Ig data for strip and circular foundations
(Ingra and Baecker (1983)). Also plotted in Figure 8.13 there is data for test FT2K obtained
0.2
0.15
0.1
0.05
i) 56.6mm dia foolIng
Meyerhof (1953)
Han",,, (1961b)
Dean .et al. (1993)
FTlK zero deg. cone (semi-rough)
+
Ingra and ilaecker(1983)
i)
_,I. -+
".+
'.
o
o
,
0.2 0.4 0.6
V/V
M
0.8
Figure 8.13: Interaction loci for surface footings with inclined load (after Dean et al.
(1993))
from a drum centrifuge which stands alone from the rest of the data, This is a special test
8-21
known as sideswipe test in which the footing is allowed to penetrate vertically until it is
locked and there is no further penetration. It is then displaced laterally to bring the sand to
collapse. The bearing capacity calculations described in the next sub-section do not represent
sideswipe tests but depict the bearing capacity of a spudcan footing at working loads.
As mentioned earlier that V M is influenced by the depth of penetration D p' It implies that
Dp will also have an affect on H/V M (see Figure 8.14a where Hp is the peak lateral
resistance). Figure 8.14b shows the variation of H /V M with wedge angles. In the absence of
satisfactory method of calculation for the axi-symmetric case, this figure shows a series of
stress field calculations for the lateral loading of plane wedges. The variation of H /VM with
wedge angle and experimental results corrected for penetration are shown in Figure 8.14c.
Figure 8.14b implies that the peak lateral resistance does not vary significantly for wedge
angles less than 30 which seems to match with experimental observations for conical
footings for cone angles up to 30 (Figure 8.14c). However, for cone angles greater than 30
there is a significant increase in lateral capacity. In the light of these theoretical and
experimental results, it is difficult to confirm that a spudcan footing of cone angle 13.46
0
will
have more sliding resistance compared to an equivalent area flat footing unless some future
research will show otherwise. 0Ne shall see later in Figure 8.23 that the results of the current
finite element analyses for the spudcan footing match quite well with Hambly's (1992) results
for the flat footing and the spudcan footing with a partial penetration factor of 2. This
agreement between the two results, also seems to confirm what we learn from Figure 8.l4a
i.e. the ultimate lateral capacity of a spudcan footing with a cone angle n 3.46
0
is almost
the same as for an equivalent bearing area flat footing.)
8.6.2.2 Finite element analysis result (load path 2)
One of the limitations of the ABAQUS' package is that it is not possible to unload from the
8-22
Chapter 8
Bearing capacity of a jackup spudCQII footillg
E
>
.......
c.
.:;.
0.5.-----------------------,
0.4
BEST FIT LINEAR EQUATION
{ 0.3
I
0-
:I:
0.2
o. ,
.A FT4K
FT2K

a 0.2 0.4
o/B
0.6 0.8
Figure 8.14a: Effect of depth of overburden Oil lip IV", (after Tan 1990))
0.3 WEDGE FOOllNGS (FULLY ROUGH)
v.2
E

Q.
I
0.1
; ":30 ------;,P
; =35-
- EMPIRICAL EQUATION
COMPUTED POINTS

o 20 40 60
WEDGE ANGLE (DEGREES)
Figure 8.14b: Effect of wedge angle on peak of interaction locus (after Tan, 1990)
06.-------------------__
(9
a
28J mm. 6)
/
0.4
0 . .)
SPUO-CAN
02
I
0.1
0
0
20 0 fa
CONE ANGLE n
Figure 8.14c: Effect of cone angle on 'corrected' (H plV M) (after Tan (1990
8-23
Chapter 8
-
,....
c::
2
E1.5
'-"
.....
,....
~
E
o
u
~ I
0.
(/)
"0
-~
.S:!
.....
~ 0.5
>
o
i
r-------I
I
- ~
I
I
!
o 0.5 1 1.5 2
Horizontal displacement (mm)
Figure 8.15: Vertical-horizontal displacement curve
:
I
I
I
~
I
I
I
I
2.5 3
9 ~ - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ,
8
~ - - - - - - - - - - - - - - .. - - - - - - - - . - - - - - . ~
7 -
6
3
2
o ~ - - - - - - - - L - - - - - - - ~ - - - - - - - - ~ - - - - - - - - ~ - - - - - - - - ~ - - - - - - - - ~
o 0.5 U 2 2.5 3
Figure 8.16: Vertical loadhorizontal displacement curve
preload described in load path 1. In order to proceed to the analysis of load path 2, therefore,
it was necessary to apply a vertical load lower than the preload and then apply suitable
horizontal loads. In fact a vertical displacement of 1.59mm was followed by a horizontal
8-24
Chapter 8
5 r - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - ~
4
0.5 1.5 2 2.S 3
Figure 8.17: Horizontal load-horizontal displacement curve
9
8
7 -
6
-
--
Z
5 ~ -
:B
........-
>
4 ~
3 ~
2 -
I
.-
0
0 2 3
4
H (MN)
Figure 8.18: Vertical-horizontal load curve
displacement of 2.5mm (see Figure 8.15).
The venical load response during this loading path is plotted against horizontal displacement
8-25
Chapter 8
(Jl
....
"
0 0
C\J C\J
I-
Z
W
ILl ::I:
ll) W
a:
ill u
lD z:
H
(TJ
C\J
",
n..
.... uJ
(Jl l-
I U1
Cl-
W W
U1 ::I:
I H
l-
ll)
I
co
I
""
UI
> a:
::J
:J: Ul 0
W ::J LJJ
Figure 8.19: Deformed mesh
"l
C> U
4: 0
0 [IJ a:
(TJ 4: n.
8-26
in Figure 8.16. At zero lateral displacement, the vertical load is 8.1 MN (i.e. 19.8% of the
maximum vertical load indicated in Figure 8.12) and this is reduced to 7.0 MN as the
horizontal displacements are increased. The variation of horizontal load with horizontal
displacement is shown in Figure 8.17. The horizontal load increases monotonically until a
value of 4.2 MN is reached. A plot of vertical versus horizontal load is given in Figure 8.18.
(The final values of vertical and horizontal load from Figure 8.18, which represent the failure
loads for the spudcan footing, are non-dimensionalised by the preload force obtained in
Section 8.6.1.2, and are plotted in Figure 8.23. It can be seen that the results lie quite close
to Hansen's (1970) curve for a flat circular footing.)
A plot of the deformed mesh at the end of the analysis, viewed from above, is shown in
Figure 8.19.
8.6.2.3 A further reduction of vertical load followed by further increase in horizontal
and vertical loads (load path 3)
A further calculation was carried out in which a vertical displacement of 0.66mm was
prescribed to the footing followed by further vertical and horizontal displacements as shown
in Figure 8.20. The calculated horizontal and vertical loads are shown in Figures 8.21 and
8.22. (The final values of vertical and horizontal load from Figure 8.22. which represent the
failure loads for the spudcan footing. are non-dimension ali sed by the preload force obtained
in Section 8.6.1.2, and are plotted in Figure 8.23. It can be seen that, again, the results lie
quite close to Hansen's (1970) curve for a flat circular footing.)
8.7 Comparison of results
Figure 8.23 shows Hansen's (1970) solution for flat circular footings, in the form of a
parabolic curve in a vertical load versus horizontal load space (interaction diagram). Both the
finite element and Hambly's (1992) laboratory test results are also plotted in this figure. As
8-27
-------.
".
L----.---
.--------
---------
V--
o
o 0.5 I U 2 25 3
Figure 8.20: Vertical-horizontal displacement curve
.5
4
o ~ - - - - - - - - ~ - - - - - - - - ~ - - - - - - - - ~ - - - - - - - - ~ - - - - - - - - ~ - - - - - - ~
o o.S I U 2
Figure 8.21: Horizontal load-horizontal displacement curve
we can see, the results match reasonably well, which confirms what the finite element
analysis predicts: a vertical preload which is the precise ultimate bearing capacity of the
8-28
Chapter 8 Bearing capacuy of a jackup spudcan footing
7
r-
Z
~
'""d
4
("j
..9
c::I
3 -
u
'e
:u
>
2
.-
o ~ - - - - - - - - - - - - ~ - - - - - - - - - - ~ - - - - - - - - - - - - ~ - - - - - - - - - - ~ - - - - - - - - - - - ~
o 2 3 4
Horizontal loud (MN)
Figure 8.22: Vertical-horizontal load curve
0.25 .----------------------------------------,
0.2
O.tS
H
V
M
0.1

0.05
t;,
0
__ Brinch Hansen's n.t footings
-6- En ite element load path 2
-e- fmite element load path 3
- Hambly's n't footing.
__ Hambly's spudc," footing with partial penelnotion faclDr
-.t.- Hambly'. spudcan footing without partial penetration faclOr
o ~ - - - - - - ~ - - - - - - ~ - - - - - - - - ~ - - - - - - ~ - - - - - - ~ - - - - - - ~
o 0.2 0.4 0.6 0.8 1.2
Figure 8.23: Horizontal-vertical load interaction curve
5
spudcan footing (Le. equal to the bearing capacity of a flat footing of equivalent bearing area).
The figure also shows that Hansen's solution gives conservative estimates of the bearing
8-29
Chapter 8
() '/
.. -
:l
~
::c
0.15
0.1
0.05
a
0
*
56.orr ..:n dia ~ o o l i n ~
Hansen (19611:
Dean ~ l al. (1993)
Ff2K zero deg. cone I.SC:":11-roU gh)
+ IngrJ and ilaecker (1983)
-&- fmite elem::nt load path 2
-e- fUlite dem::nt load path 3
-....... lIamoly's $pudcan fooung WIthout pa.rt:a! ,:>enc!nuon facmr
__ HamolY' ~ spuacan :OOtl.'l& WIth parwu pe!'lc:'''l\uon (acUJr
~
++-
.. '
-+-
+. *+
+
.* *
.4-
0.2 OA-
vh
M
0.6 0.8 1
Figure 8.24: Interaction loci for surface footings (after Dean et at. (1993
capacity of such offshore footings. The finite element analysis results do not match with
Hambly's (1992) laboratory result for spudcan footing without a partial penetration factor.
Both the current finite element and Hambly's (1992) laboratory results are replotted on
another interaction diagram containing data from other sources (see Figure 8.24 which is
actually the same Figure 8.13 with additional data points). It can be seen that Hambly's test
result for a spudcan footing lies far away from the rest of the data points. This inconsistency
is resolved by Hambly (1992) with the introduction of partial penetration factors. However,
8-30
the validity of these partial penetration factors need further experimental evidence.
8.8 Conclusions
The finite element method predicts a vertical preload force which is the ultimate bearing
capacity of a spudcan footing and match with Hambly's (1992) laboratory result for model
flat circular footing and the result for the spudcan footing, using a partial penetration factor
of 2. (Although, further investigations to validate partial penetration factors suggested by
Hambly' (1992), are proposed.)
The finite element results fits well within the centrifuge test data reported by Dean et al.
(1993). However, Hansen's (1970) theory gives conservative estimates of the sliding/bearing
resistance of a spudcan footing on sand.
8-31
Chapter 9 Conclusions
CHAPTER 9
CONCLUSIONS
9.1 Concluding remarks
Most offshore oil extraction platfoffi1s are pile supported fixed structures. Groups of closely
spaced piles around each of the four columns of the structure pin it to the sea-bed. Careful
study of pile-soil-pile interactions in these pile groups is necessary in order to predict their
behaviour under severe gravity and environmental (lateral) loads. The analysis of a pile
foundation involves the solution of a three-dimensional problem in which two different
materials act together to withstand the external applied loads. Each of these materials exhibits
different stress-strain characteristics. It is difficult to put together an analysis procedure able
to capture closely the combined behaviour of the two materials and so determine the true
response of a pile foundation. Various idealisations and possible approximations need to be
made in order to represent the material constitutive behaviour and the geometry of the
problem. The methods that have been developed since around 1967 have different degrees of
rigour depending upon the approximations involved in the analysis. In addition to these
approximations, numerous other factors exist which can significantly affect the pile foundation
response. These include: the sequence of loading, involving initial disturbance of the soil in
close proximity to the pile due to installation, and subsequent re-consolidation of the soil in
the intermediate time of recovery prior to actual loading; non-linearity and inhomogeneity of
the soil; plastic deformations in the soil around the pile; relative slip between the pile and the
soil at their interface under axial loading; tension and gap formation at the back of the pile
9-1
under lateral loading. The only sufficiently versatile analysis method, capable of considering
all these factors satisfactorily, is the finite element method. Of course a full, three-dimensional
finite element analysis requires a huge amount of computational effort; nevertheless, it is quite
economical compared to full scale field experiments on pile groups. Given the present
impracticability of either a full three-dimensional analysis or comprehensive field experiments,
an alternative was presented in Chapter 6.
The objective of the thesis was to investigate a suitable method for the analysis and design
of pile foundations. To this end the integral equation methOd/quasi boundary element method
was formulated in order to compare its perfonnance with the finite element method. Although,
when considering elastic soil behaviour, it is clearly superior, in tern1S of computational effort
required, over the finite element method, it does not give complete stress distribution inside
the soil medium. If adapted to meet this demand, it starts to function in a way similar to the
finite element method. It also becomes less efficient than the finite element method when a
more refined (i.e. elasto-plastic) soil behaviour needs to be considered. It was not deemed
worthwhile to try to develop the quasi boundary element method into a general two-
dimensional boundary element method combined with elasto-plastic soil modelling - the finite
element method was already available to do this. However, it was necessary to find a new,
computationally efficient approach for performing the finite element calculations. The new
approach, presented in Chapter 3, also models the rigid or infinite boundaries in the finite
element mesh.
It was found, further, that a two-dimensional idealisation of single-piles and pile groups
subject to lateral loading is quite successful in predicting their response, just as it is for the
case of axial loading.
Three mains areas of single-pile and pile group behaviour were studied in this thesis. These
are: axial load capacity and axial response of piles (9.2); interface analysis (9.3); and lateral
9-2
load capacity and response of single-piles and pile groups (9.4).
9.2 Axial load capacity and axial response of piles
Chapter 3 of this thesis demonstrated that, for the preliminary estimation of vertical settlement
of piles, the analytical approach presented by Randolph and Wroth (1979) is very practicable
and useful. The axial capacity cannot be predicted accurately unless the distribution of shear
stress (skin friction) along the pile-soil interface is precisely known. Neither the subgrade
reaction approach/load transfer method (t-z curve analysis), nor the elastic boundary element,
can determine this with sufficient accuracy. Only the finite element method is competent to
deal with the complex pile-soil interface behaviour and predict the axial load capacity with
a high degree of accuracy. A new approach to deal with the size of the finite element mesh
and the infinite boundaries was also presented in this chapter. With this new approach it is
possible to study the behaviour of closely spaced pile groups in which, in reality, the pile-soil-
pile interactions are non-linear.
9.3 Interface analysis
After highlighting the importance of pile-soil interface behaviour in Chapter 3, a cohesive
interface element was developed in Chapter 4 to detemline the stress distribution at the pile-
soil interface. This element offers an easy way to vary the pile roughness. It models small
relative pile-soil movement (slip), tension at the back of the pile when the pile is loaded
laterally, and the formation of a gap once the limiting tension exceeds some prescribed value.
The interface element was also quite successful in the stability analysis of surface footings.
Surprisingly, and in contrast to some of the earlier work done in this area (reported in Chapter
4), no numerical instabilities were found. In addition to the interface analysis, three-
dimensional finite element analysis was carried out to predict the ultimate bearing capacity
9-3
of a partially penetrated spudcan footing (Chapter 8). The analysis results help to explain the
applicability of the existing bearing capacity theories to conical footings.
9.4 Lateral load carrying capacity and response of single piles and pile
groups
In Chapters 5, 6 and 7, it was demonstrated that two-dimensional idealisation of a laterally
loaded pile or a pile group was quite suitable for studying the lateral reaction mechanisms in
both cohesive and frictional soils. The new material model introduced (in Chapter 7) to
represent the constitutive behaviour of the frictional soils successfully models dilation and the
post-peak behaviour of the soils. The simplification (two-dimensional idealisation) of the
problem was adopted to study pile-soil-pile interactions in closely spaced offshore pile groups.
It adds significantly to our understanding of the overlapping plastic zones within the group.
A method of pile group analysis was then suggested which combines plane strain finite
element analysis for predicting the pile head response with a p-y curve approach for
determining the bending moments and deflections along the piles. The p-y curves used are
theoretical curves obtained by a single plane strain finite element analyses along the length
of the piles.
9.5 Future research
The research work undertaken in this thesis is related to undrained loading of soils. As
demonstrated and described in Chapter 5, a significant amount of saving in the design of pile
foundations can be made if an effective stress design approach is followed. The consolidation
analysis performed in Chapter 5, in which the increase in shear strength of the soil around
the pile (taken from the logarithmic distribution of the shear strength after consolidation
9-4
obtained by Randolph et al. (1979)) was approximately incorporated in the finite element
analysis based on the Cam-Clay model (Roscoe and Schofield (1963)), is not sufficiently
refined. The deficiency could be made up by doing a two-dimensional large displacement
coupled consolidation analysis using a sophisticated work hardening soil model (modified
Cam-Clay model (Roscoe and Burland (1968)). This would be a significant leap forward
towards the development of effective stress analysis and economical design of pile
foundations.
Since, for the present, a complete three-dimensional analysis is prohibitively expensive, a new
versatile analysis approach is suggested in Chapter 6. However, another approach for pile
groups whose one dimension is significantly larger than the other, is to carry out a single two-
dimensional finite element analysis by developing a two-dimensional beam element
compatible with the soil element to represent the piles. The problem of rigid boundaries in
the finite element mesh could be resolved by developing an infinite element. In this case we
would have to compromise on the effect of the stiffness of the pile element and, instead, look
for a balanced out pile-soil stiffness representing the whole group, as for example suggested
by Stewart et al. (1993). Interface analysis could be helpful in determining active and passive
soil pressures around the piles. Large displacement formulations could also be useful since,
once a gap has opened up at the back of the pile, a small displacement is unlikely.
Consolidation analysis could also be coupled together. This approach, too, we can see as an
approximation of what is in reality a three-dimensional problem. It would not be possible to
capture the true collapse mechanism in the soil using this approach, simply because of its
inability to model overlapping plastic zones around individual piles in a group.
Another possibility would be to make use of the axial symmetry of the problem in the finite
element analysis and apply the loading in terms of Fourier harmonics in order to study the
lateral response of single piles and pile groups. This work too could be extended to coupled
consolidation analysis, as mentioned earlier.
9-5
An immediate work, as proposed in Chapter 7 is to carry out work similar to the one
described in Chapters 5 and 6, using the Matsuoka model (2) presented in Chapter 7.
9-6
Bibliography
BIBLIOGRAPHY
Aldridge, T.R. and Erbrich, C. (1990)
"A new metJwd for the prediction of lateral pile group behaviour - FUGROUP",
Proc. 2nd European Conf. on Numerical Methods in Geotechnical Engineering,
Santander: pp 349-353
Al Tabbaa, A. and Wood, D.M. (1989)
"An experimentally based 'bubble' model for clay",
Proc. 3rd Int. Symposium on Numerical Models in Geomechanics, Niagra, Canada:
pp 91-99
Arnesen, K., Dahlberg, R., Helge, K. and Carlsen, C.A. (1988)
"Soil-structure-interaction aspects for jackup platforms",
BOSS, Trondheim
Atkinson, J.H. and Bransby, P.L. (1982)
"The Mechanics of Soils",
The McGraw-Hill Publishing Company Ltd., London, UK: pp 310-312
Atkinson, J.H. and Stallebrass, S.E. (1991)
"A model for recent history and non-linearity in the stress-strain behaviour of
overconsolidated soil",
Proc. 7th lnt. Conf. of the Int. Association for Computer Methods and Advances in
Geomechanics, Cairns, Queensland, Australia: pp 555-560
Atkinson, J.H. (1993)
"A note on modelling small strain stiffness in Cam Clay",
Proc. Wroth Memorial Symposium, Oxford, UK: pp 628-643
Baguelin, F., Frank, R. and Said, Y.H. (1977)
"Theoretical study of lateral reaction mechanism of piles",
Geotechnique, VoL 27, No.3: pp 405-434
Baguelin, F. and Frank, R. (1980)
"Theoretical studies of piles using the finite element method",
Proc. lnst. Civ. Engrs. Conf. on Numerical Methods in Offshore Piling, London: Paper
No. 11
Bibliography
Baligh, M.M. (1986)
"Undrained deep penetration, I: shear stresses (and) II: pore pressures",
Geotechnique, Vol. 36, No.4: pp 471-502
Banerjee, P.K. and Driscoll, P.M. (1976)
"Three dimensional analysis of raked pile groups",
Proc. Inst. Civ. Engrs. Part 2, Vol. 61: pp 653-671
Banerjee, P.K. and Davies, T.G. (1978)
"The behaviour of axially and laterally loaded single piles embedded in non-
homogeneous soils",
Banerjee, P.K. and Mustoe, G.G.W. (1978)
"The boundary element method for two-dimensional problems of elasto-plasticity",
Proc. Int. Conf. on Recent Advances in Boundary Element Method, Pen tech Press
London: pp 309-326
Banerjee, P.K. and Davies, T.G. (1980)
"Analysis of some reported case histories of laterally loaded pile groups",
Proc. lnst. Civ. Engrs. Conf. on Numerical Methods in Offshore Piling, London: pp
83-91
Banerjee, P.K. and Butterfield, R. (1980)
"Boundary element methods in engineering science",
McGraw-Hill, London, UK
Banerjee, P.K. and Dargush, G.P. (1988)
"Progress in BEM applications in geomechanics via examples",
Proc. Conf. on Numerical Methods in Geomechanics, Innsbruck: pp 13-22
Barden, L. Ismail, H. and Tong, P. (1969)
"Plane strain deformation of granular material at low and high pressures",
Barton, Y.O. (1982)
"Laterally loaded model piles in sand: centrifuge tests and finite element analyses",
Ph.D Thesis submitted to the Department of Engineering, University of Cambridge.
Bell, R.W. (1991)
"The analysis of offshore foundations subjected to combined loading",
M.Sc. Thesis, University of Oxford, UK
Bettes, P. (1977)
"Infinite elements",
Int. Journal for Numerical Methods in Engineering, Vol. 11: pp 53-64
Bibliography
Bogard, D. and Matlock, H. (1983)
"Procedures for analysis of laterally loaded pile groups",
Proc. Conf. on Geotechnical Practice in Offshore Engineering, Austin, Texas: pp 499-
535
Bolton, M.D. (1979)
"A guide to soil mechanics",
The Macmillan Press Ltd., London, UK
Bolton, M.D. (1986)
"The strength and dilatancy of sands",
Bond, A.J. and Jardine, RJ. (1991)
"Effects of installing displacement piles in a high OCR clay",
Brocklehurst, c.J. (1993)
"Finite element studies of reinforced and unreinforced two-layer soil systems",
D.PhiL Thesis, University of Oxford, U.K.
Broms, RR (1964a)
"Lateral resistance of piles in cohesive soils",
Proc. ASCE, Journal of the Soil Mechanics and Foundation Engineering Division, 90
(SM3): pp 27-63
Broms, RB. (1964b)
"Lateral resistance of piles in cohesionless soils",
(SM3): pp 123-156
Brown, D.A. and Shie, c.-F. (1990a)
"Numerical experiments into group effects on the response ofpiles to lateral loading",
Computers and Geotechnics, Vol 10: pp 211-230
Brown, D.A. and Shie, c.-F. (1990b)
"Three dimensional finite element model of laterally loaded piles",
Computers and Geotechnics, Vol 10: pp 59-79
Burd, H.J. (1986)
"A large displacement finite element analysis of a reinforced unpaved road",
D.Phil. Thesis, University of Oxford, UK
Burd, H.J. and Houlsby, G.T. (1990)
"Finite element analysis of two cylindrical expansion problems involving nearly
incompressible material behaviour".
Int. Journal for Numerical and Analytical Methods in Geomechanics, Vol. 14: pp 351-
366
Bibliography
Burd, H.J. (1991)
"Report submitted to Edmund Hambly Ltd"
Burd, H.J. and Brocklehurst, C.J. (1991)
"Parametric studies of a soil reinforcement problem using finite element analysis",
Proc. 7th Int. Conf. of the Int. Association for Computer Methods and Advances in
Geomechanics, Cairns, Queensland, Australia, Vol. 3: pp 1783-1788
Burd, H.J. and Houlsby, G.T. (1992)
"Personal communication"
Butterfield, R. and Banerjee, P.K. (1971)
"The elastic analysis of compressible piles and pile groups",
Burland, J.B. (1973)
"Shaft friction of piles in clay - a simple fundamental approach",
Ground Engineering, Vol. 6, No.3: pp 30-34
Calladine, c.R. (1985)
"Plasticity for Engineers",
Ellis Horwood Ltd., Chichester, West Sussex, UK
Carter, J.P., Randolph, M.E and Wroth, c.P. (1979)
"Stress and pore pressure changes in clay during and after the expansion of a
cylindrical cavity",
322
Chen and Poulos (1993)
"Analysis of pile-soil interaction under lateral loading using infinite and finite
elements",
Computers and Geotechnics, Vol. 15: pp 189-220
Cheung, Y.K., Tham, L.G. and Guo, DJ. (1988)
"Analysis of pile group by infinite layer method",
Cheung, Y.K., Lee, P.K.K. and Zhao, W.B. (1991)
"Elastoplastic analysis of soil-pile interaction",
Computers and Geotechnics, Vol. 12: pp 115-132
Chow, Y.K. (1986)
"Analysis of vertically loaded pile groups",
Int. Journal for Numerical and Analytical Methods in Geomechanics, Vol. 10, No.1:
pp 59-72
Bibliography
Chow, Y.K. (1987)
"Three dimensional analysis of pile groups",
Proc. ASCE, Journal of the Geotechnical Engineering Division, 113 (GT6): pp 637 -651
Chow, L. (1994)
"The prediction of surface settlements due to tunnelling in soft ground",
M.Sc. Thesis, University of Oxford, UK
Clausen, C.S.F., Aas, P.M. and Hasle, E. (1984)
"SPLICE - a computer programfor analysing structure/pile/soil interaction programs",
Norwegian Geotechnical Institute Publication 152 Oslo
Coop, M.R and Wroth, c.P. (1989)
"Field studies of an instrumented model pile in clay",
Coyle, H.M. and Reese, L.c. (1966)
"Load transfer for axially loaded piles in clay",
(SM2): pp 1-26
Dean, E.T.R., James, RG., Schofield, A.N., Tan, F.S.C. and Tsukamoto, Y. (1993)
"The bearing capacity of conical footings on sand in relation to the behaviour of
spudcan footings of jackups",
De Borst, R. and Vermeer, P.A. (1984)
"Possibilities and limitations of finite elements for limit analysis",
De Josselin de Jong, G. (1971)
"The double sliding free rotating model for granular assemblies",
Geotechnique, Vol. 21: pp 155-163
De Josselin de Jong, G. (1988)
"Elastic-plastic version of the double sliding model in undrained simple shear tests",
Desai, C.S., Zaman, M.M., Lightner, J.G. and Siriwardane, H.J. (1984)
"Thin-layer elements for interfaces and joints",
43
De Santa Maria, P.E.L. (1988)
"Behaviour of footings for offshore structures under combined loads",
Bibliography
Drucker, D.C. and Prager, W. (1952)
"Soil mechanics and plastic analysis or limit design",
Quarterly Journal of Applied Mathematics, Vol. 10, No.2: pp 157-165
Duncan, J.M. and Chang, C.Y. (1970)
"Nonlinear analysis of stress and strain in soils",
(SM5): pp 1629-1653
Evangelista, A. and Viggiani, C. (1976)
"Accuracy of numerical solutions for laterally loaded piles in elastic half space",
Proc. 2nd Int. Conf. Numerical methods in Geomechanics, Blacksburg, Vol 3: pp
1367-1370
Fleming, W.G.K., Weltman, AJ., Randolph, M.F. and Elson, W.K. (1992)
"Piling Engineering",
Chapman and Hall, Glasgow, 2nd edition.
Frank, R, Guenor, A. and Humbert, P. (1982)
"Numerical analysis of contacts in geomechanics",
Proc. 4th Int. Conf. on Numerical Methods in Geomechanics, Edmonton, Canada, Vol.
1: pp 37-45
Gens, A., Carol, I. and Alonso, E.E. (1989)
"An interface element formulation for the analysis of soil-reinforcement interaction",
Computers and Geotechnics, Vol. 7, No 1-2: pp 133-151
Georgiadis, M. and Butterfield, R. (1982)
"Laterally loaded pile behaviour",
Proc. ASCE, Journal of the Geotechnical Engineering Division, 108 (OT1): pp 155-165
Ghaboussi, J., Wilson, E.L. and Isenberg, J, (1973)
"Finite element for rock joints and interfaces",
(SMlO): pp 833-848
Gibson, R.E. and Anderson, W.F. (1961)
"In situ measurement of soil properties with the pressuremeter",
Civil Engineering and Public Works Review, Vol. 56, No. 658: pp 615-618
Goodman, RE., Taylor, RL. and Brekke, T.L. (1968)
"A model for the mechanics of jointed rock",
(SM3): pp 637-659
Bibliography
Griffiths, D.V. (1985)
"Numerical modelling of interfaces using conventional finite elements",
Proc. 5th Int. Conf. on Numerical methods in Geomechanics, Nagoya, Japan: pp 837-
844
Guo, D.J . Tham. L.G. and Cheung, Y.K. (1987)
"Infinite layer for the analysis of a single pi/e",
Computers and Geotechnics, Vol. 3, No.4: pp 229-249
Hambly, E.C. (1992)
"Jackup spudcan sliding/bearing resistance on sand",
BOSS, London
Handel, E., Schweiger, H.F. and Yeo, K.C. (1990)
itA simple thin-layer element to model soil-geotextile interaction",
Int. Conf. on Reinforced Soil, Strathclyde, UK: pp 10-12
Hansen, J.B. (1961 a)
"The ultimate resistance of rigid piles against transversal forces",
Geoteknisk Institut. Bull. No. 12, Copenhagen.
Hansen, J.B. (1961b)
"A general formula for bearing capacity".
Danish Geotechnical Institute, Copenhagen, Bull. No. 11
Hansen, J.B. (1970)
"A revised and extended formula for bearing capacity",
Danish Geotechnical Institute, Copenhagen, Bull. No. 28: pp 38-46
Herrmann, L.R. (1965)
"Elasticity equation for incompressible and nearly incompressible materials by a
variational theorem",
Journal of American Institute of Aeronautics and Astronautics, 3, 1896-1900
Herrmann, L.R. (1978)
"Finite element analysis of contact problems",
Proc. ASCE, Journal of the Engineering Mechanics Division, 104 (EMS): pp 1043-
1057
Hetenyi, M. (1946)
"Beams on elastic foundations",
Univ. of Michigan Press, Michigan.
Heydinger, A.G. and O'Neil, M.W. (1986)
"Analysis of axial pile-soil interaction in clay",
381
Bibliography
Hibbit, Karlsson and Sorensen (1989)
"ABAQUS - Theory Manual, version 4.8"
Hill, R. (1950)
"The Mathematical Theory of Plasticity",
Oxford University Press, London
Houlsby, G.T. (1982)
"Theoretical analysis of the fall cone test",
Houlsby, G.T. and Wroth, C.P. (1982)
"Direct solution of plasticity problems in soils by the method of characteristics",
Proc. 4th Int. Conf. on Numerical Methods in Geomechanics, Edmonton, Canada, Vol.
3: pp 1059-1071
Houlsby, G.T. and Wroth, C.P. (1983)
"Calculations of stresses on shallow penetrometers and footings",
Proc. Int. Union of Theoretical and Applied MechanicslInt. Union of Geodesy and
Geophysics Symposium on Seabed Mechanics, Newcastle upon Tyne, pp 107-112
"How the dilatancy of soils affects their behaviour",
Report No. OUEL 1888/91, Engineering Science Department, University of Oxford,
UK
Houlsby, G.T. and Martin, C.M. (1993)
"Modelling of the behaviour of foundations of jack up units on clay",
Proe. Wroth Memorial Symposium, Oxford, UK: pp 235-253
Houlsby, G.T. and Chow, L. (1994)
"Modelling the variable stiffness of undrained clay using multiple yield surfaces",
Report No. OUEL 2024/94, Engineering Science Department, University of Oxford,
UK
Ingra, T.S. and Baecker, G.B. (1983)
"Uncertainty in bearing capacity of sands",
Proe. ASCE, Journal of the Geotechnical Engineering Division, 109 (G17): pp 899-914
Jardine, R.J. and Potts, D.M. (1988)
"Hutton tension leg plat/orm foundations: prediction of driven pile behaviour",
Bibliography
Jardine, RJ., Potts, D.M., St. John, H.D. and Hight, D.W. (1991)
"Some applications of a non-linear ground mode!",
Proc. 10th European Conf. on Soil Mechanics and Foundation Engineering, Florence,
Vol. 1: pp 223-228
Katona, M.G. (1983)
"A simple contact-jriction interface element with applications to buried culverts",
Int. Journal for Numerical and Analytical Methods in Geomechanics, Vol. 7, No.3:
pp 371-384
Kodikara, J.K. and Moore, LD. (1991)
"A versatile interaction analysis for nonlinear bodies",
Proc. 7th lnt. Conf. of the Int. Association for Computer Methods and Advances in
Geomechanics, Cairns, Queensland, Australia, Vol 2: pp 1171-1176
Kraft, L.M., Jr., Ray, RP. and Kagawa, T. (1981)
"Theoretical t-z curves",
Proc. ASCE, Journal of the Geotechnical Engineering Division, 197 (GTl1): pp 1543-
1561
Kooijman, A.P. and Vermeer, P.A. (1988)
"Elasto-plastic analysis of laterally loaded piles",
Kooijman, A.P. (1989)
"Comparison of an elasto-plastic Quasi three dimensional model for laterally loaded
piles with field tests",
Proc. 3rd. Int. Symposium on Numerical Models in Geomechanics, Niagara, Canada:
pp 675-682
Kuhlemeyer, RL. (1979)
"Static and dynamic laterally loaded floating piles",
Proc. ASCE, Journal of the Geotechnical Engineering Division, 105 (GT2): pp 289-304
Lane, P.A. and Griffiths, D.V. (1988)
"Computation of the ultimate pressure of a laterally loaded pile in frictional soil",
Langen, H. and Vermeer, P.A. (1990)
"Finite element analysis of a pile penetration problem in clay",
Proc. 2nd European Conf. on Numerical Methods in Geotechnical Engineering,
Santander: pp 18-20
Leung, C.F. and Chow, Y.K. (1987)
"Response of pile groups subjected to lateral loads",
lnt. Journal for Numerical and Analytical Methods in Geomechanics, VoL 11, No.3:
pp 307-314
Bibliography
Mandolini, A. and Viggiani, C. (1993)
"Settlement predictions for piled foundations from loading tests on single piles",
Matlock, H. and Reese, L.C. (1960)
"Generalised solutions for laterally loaded piles",
(SM5): pp 63-91
Matlock, H. (1970)
"Correlation for design of laterally loaded piles in soft clay",
Proc. 2nd Offshore Technology ConL, Houston, Vol 1: pp 577-594
Matlock, H., Ingram, W.B., Kelly, A.E. and Bogard, D. (1980)
"Field tests on the lateral behaviour of pile groups in soft clay",
Proc. 12th Offshore Technology Conf., Houston, Vol 1: pp 163-174
Matsuoka, H. (1976)
"On the significance of the spatial mobilised plane",
Soils and Foundations, Vol. 16, No.1: pp 91-100
Mattes, N.S. (1969)
"The influence of radial displacement compatibility on pile settlements",
Geotechnique, VoL 19: pp 157-159
Mattes, N.S. and Poulos, H.G. (1969)
"Settlement of single compressible pile",
(SMl): pp 189-207
Medina, F. (1981)
"An axisymmetric infinite element",
Meyerhof, G.G. (1953)
"The bearing capacity of foundations under eccentric and inclined loads",
Proc. 3rd Int. Conference on Soil Mechanics, Zurich, Vol 1: pp 440-445
Mindlin, R.D. (1936)
"Force at a point in the interior of a semi-infinite solid",
Physics, Vol 7: pp 195-202
Moore, I.D. and Booker, J.R. (1989)
"Geometrically non-linear analysis of buried cylinders",
Proc. 3rd Int. Symposium on Numerical Models in Geomechanics, Niagra, Canada:
pp 716-723
Bibliography
Mroz, Z., Norris, V.A. and Zienkiewicz, O.c. (1979)
"Application of an anisotropic hardening model in the analysis of elasto-plastic
deformation of soils",
Muqtadir, A. and Desai, C.S. (1986)
"Three-dimensional analysis of pile group foundation",
lnt. Journal for Numerical and Analytical Methods in Geomechanics, Vol. 10: pp
41-58
Nagtegaal, J.e., Parks, D.M. and Rice, J.R. (1974)
"On numerically accurate finite element solutions in the fully plastic range",
Computer Methods in Applied Mechanics and Engineering, Vol. 4: pp 153-177
Naylor, D.J. (1974)
"Stresses in nearly incompressible materials by finite elements with application to the
calculation of excess pore pressures",
Nutt, N.R.F. (1993)
"Development of the cone pressuremeter",
O'Neill, M.W., Ghazzaly, OJ. and Ha, H.B. (1977)
"Analysis of three dimensional pile groups with non-linear soil response and pile-soil-
pile interaction",
Proc. 9th Offshore Technology Conf, Houston, Texas, Vol 2: pp 245-256
Osborne, J.l., Trickey, J.c., Houlsby, G.T. and James, R.G. (1991)
"Findings from a joint industry study on foundation fixity of jackup units",
Offshore Technology Conf, Houston, Texas: paper OTC6615
Ottaviani, M. (1975)
"Three dimensional finite element analysis of vertically loaded pile groups",
Pasternak, P.L. (1954)
"On a new method of analysis of an elastic foundation by means of two foundation
constants" ,
Moscow: GosJzd.Lit.Po Stroit-i-Arch (Russian).
Potts, D.M. and Martin, J.P. (1982)
"The shaft resistance of axially loaded piles in clay",
Poulos, H.G. and Davis, E.H. (1968)
"The settlement behaviour of single axially loaded incompressible piles and piers",
Geotechnique, Vol 18: pp 351-371
Bibliography
Poulos, H.G. (l971a)
"Behaviour of laterally loaded piles: I-Single piles",
(SM5): pp 711-731
Poulos, H.G. (l971b)
"Behaviour of laterally loaded piles: II-Pile Groups",
(SM5): pp 733-751
"Elastic Solutions for Soil and Rock Mechanics",
Wiley: New York
Poulos, H.G. (1976)
"Behaviour of laterally loaded piles near a cut or slope",
Journal of Australian Geomechanics, Vol. G6, No. 1: pp 6-12
Poulos, H.G. (1980)
"An approach for the analysis of Offshore pile groups",
Proc. lnst. Civ. Engrs. Conf. on Numerical Methods in Offshore Piling, London: pp
119-126
Poulos, H.G. (1980)
"User's guide to DEFPIG - deformation analysis of pile group",
School of Civil Engineering, University of Sydney
"Pile foundation analysis and design",
Wiley New York.
Poulos, H.G. (1988)
"Marine Geotechnics",
Unwin Hyman, London, UK
Poulos, H.G. (1989)
"Pile behaviour - theory and application",
Randolph, M.F. (1977)
"A theoretical study of the performance of the piles",
Ph.D Thesis submitted to the Department of Engineering, University of Cambridge.
Randolph, M.F. and Wroth, c.P. (1978)
"Analysis of deformation of vertically loaded piles",
Proc. ASCE, Journal of the Geotechnical Engineering Division, 104 (GT12): pp 1465-
1488
Bibliography
Randolph, M.F. and Wroth, c.P. (l979a)
"An analysis of the vertical deformation of pile groups",
Randolph, M.F. and Wroth, c.P. (1979b)
"An analytical solution for the consolidation around a driven pile",
lnt. Journal for Numerical and Analytical Methods in Geomechanics, Vol. 3, No.3:
pp 217-229
Randolph, M.F., Carter, J.P. and Wroth, C.P. (1979)
"Driven piles in clay - The effects of installation and subsequent consolidation",
"PIGLET, a computer program for the analysis and design of piles under general
loading conditions",
Cambridge University Research Report, CUED/D-Soils TR91
"The response of flexible piles to lateral loading",
Randolph, M.F. and Houlsby, G.T. (1984)
"The limiting pressure on circular pile loaded laterally on cohesive soil".
"RATZ - load transfer analysis of axially loaded piles",
Rep. No. Geo 86033, Perth: Dept. Civil Engineering, University of Western Australia
Randolph, M.F. and Springman, S.M. (1991)
"Analysis of pile response due to external loads and soil movement",
Proc. 10th European Conf. on Soil Mechanics and Foundation Engineering, Florence.
Reese, L.C., Cox, W.R. and Koop, F.D. (1974)
"Analysis of laterally loaded piles in sand",
Proc. 6th Offshore Tech. Conf., Houston, Texas, paper OTC 2080: pp 473-483
Reese, L.C., Cox, W.R. and Koop, ED. (1975)
"Field testing and analysis of laterally loaded piles in stiff clay",
Proc. 7th Offshore Tech. Conf., Houston, Texas, paper OTC 2312: pp 671-690
Reese, L.C. (1975)
"Laterally loaded piles",
Proc. ASCE, Geotech. Group and Continuing Education Committee, San Francisco
section, Univ. of Calif., Berkeley
Bibliography
Reese, L.c. (1977)
"Laterally loaded piles: Program Documentation",
Proc. ASCE, Journal of the Geotechnical Engineering Division, 103 (GT4): pp 287-305
Roscoe, K.H. and Schofield, A.N. (1963)
"Mechanical behaviour of an idealised 'wet clay"',
Proc. 2nd European Conf. Soil Mechanics, Wiesbaden, Vol. 1: pp 47-54
Roscoe, K.H. and Burland, J.B. (1968)
"On the generalised stress-strain behaviour of 'wet clay"',
Engineering Plasticity, Cambridge University Press: pp 535-609
Rowe, P.W. (1962)
"The stress-dilatancy relation for static equilibrium of an assembly of particles in
contact".
Proc. Roy. Soc., Series A, VoL 269: pp 500-527
Sagaseta, C. (1984)
"Quasi-static undrained expansion of a cylindrical cavity in clay in the presence of
shaft friction and anisotropic initial stresses",
Unpublished Research Report, Engineering Science Department, University of Oxford,
UK
Schofield, A.N. and Wroth, C.P. (1968)
"Critical State Soil Mechanics",
The McGraw-Hill Publishing Company Ltd., London, UK: pp 32-37
Scott, R.F. (1981)
"Foundation Analysis",
Prentice-Hall, London, UK
SeInes, P.B. (1982)
"Geotechnical problems in offshore earthquake engineering",
Norwegian Geotechnical Institute Publication 140 Oslo
Simpson, B. (1993)
"Development and application of a new soil model for prediction of soil movements",
Proc. Wroth Memorial Symposium. Oxford, UK: pp 628-643
Sloan, S.W. (1981)
"Numerical analysis of incompressible and plastic solids using finite elements",
Ph.D. Thesis, University of Cambridge, UK
Sloan, S.W. and Randolph, M.P. (1982)
"Numerical prediction of collapse loads using finite element methods",
76
B ibliograph y
Smith, I.M. and Griffiths, D.V. (1988)
"Programming the Finite Element Method",
John Wiley and Sons, Chichester, New York, 2nd edition.
Springman, S.M. (1984)
"Lateral loading on piles due to embankment construction",
M.Phil. Thesis, University of Cambridge, UK
Springman, S.M. (1989)
"Lateral loading of piles due to simulated embankment construction",
Stroud, M.A. (1971)
"The behaviour of sand at low stress levels in the simple shear apparatus",
Sullivan, W. R., Reese, L. C. and Fenske, C. W. (1980)
"Unified method for analysis of laterally loaded piles in Clay",
Proc. Inst. Civ. Engrs. Conf. on Numerical Methods in Offshore Piling, London: pp
135-146
Swedlow, J.L. and Cruse, T.A. (1971)
"Formulation of Boundary Integral Equations for three-dimensional Elasto-plastic
flow",
Int. Journal for Solids and Structures, Vol. 7: pp 1673-1683
Tan, F.S.C. (1990)
"Centrifuge and theoretical modelling of conical footings on sand",
Taylor, D.W. (1948)
"Fundamentals of Soil Mechanics",
Wiley, New York
Telles, J.C.F. and Brebbia, C.A. (1981)
"New Developments in Elasto-plastic Analysis",
Proceedings of the Third International Seminar, Irvine, California: pp 350-370
Timoshenko, S.P. and Goodier, J.N. (1970)
"Theory of Elasticity",
The McGraw-Hill Book Company, Inc., New York
Verruijt, A. and Kooijman, A.P. (1989)
"Laterally loaded piles in a layered elastic medium",
Bibliography
Vesic, A.S. (1975)
"Bearing capacity of shallow foundations",
"Foundation Engineering Handbook", Ed. by H. F. Winterkorn and H. Y. Fang, Van
Nostrand Reinhold, New York, USA: pp 121-147
Vlasov, V. Z., and Leontiev, N.N. (1966)
"Beams, plates, shells on elastic foundations",
translated from Russian by Israel Program for Scientific translations, NTIS No.
N67-14238.
Wroth, C.P. and Bassett, R.H. (1965)
itA stress-strain relationship for the shearing behaviour of sand",
Wroth, C.P. and Houlsby, G.T. (1985)
"Soil mechanics - property characterization and analysis procedures",
Proc. 11th lnt. Conf. on Soil Mechanics and Foundation Engineering, San Francisco,
USA, Vol. 1: pp 1-57
Yan, L. and Byrne, P.M. (1992)
"Lateral pile response to monotonic pile head loading",
Canadian Geotechnical Journal, Vol. 29: pp 955-970
Yegian, M. and Wright, S.G. (1973)
"Lateral soil resistance - displacement relationships for pile foundations in soft clays",
Proc. 5th Annual Offshore Technology Conference, Paper OTC1893, Vol. 2: pp 663-
676
Yu, H-S. (1990)
"Cavity expansion theory and its application to the analysis of pressuremeters",
Zaman, M.M. (1985)
"Evaluation of 'thin layer element' and modelling of interface behaviour in soil-
structure interaction",
Proc. 5th lnt. Conf. on Numerical methods in Geomechanics, Nagoya, Japan: pp 1797-
1803
Zienkiewicz, D.C., Best, B., Dullage, C. and Stagg, K.G. (1970)
"Analysis of non-linear problems in rock mechanics with particular reference to
jointed rock systems",
Proc. 2nd lnt. Congo on Rock Mechanics, Belgrade, Vol. 2: pp 501-509
Zienkiewicz, D.C. (1977)
"The finite element method",
Third edition, McGraw-Hill, London, UK

Pile Soil Interaction

Uploaded by

Copyright:

Available Formats

Pile Soil Interaction

Uploaded by

Document Information

Copyright

Available Formats

Share this document

Share or Embed Document

Sharing Options

Did you find this document useful?

Is this content inappropriate?

Copyright:

Available Formats

Pile Soil Interaction

Uploaded by

Copyright:

Available Formats

What is the thesis about?

What methods are used in the analysis?

ST A TIC PILE-SOIL-PILE INTERACTION IN

OFFSHORE PILE GROUPS

You might also like