R

Finite Elements in Geotechnical Engineering

a (very) short introduction

Gioacchino (Cino) Viggiani

Laboratoire 3S-R (Sols, Solides, Structures - Risques), Grenoble, France
Université Joseph Fourier

Course “Computational Geotechnics” 1/21

 Contents
R

• Introduction

• Finite Elements (displacement) analysis

Elements for two-dimensional analysis
Displacement interpolation
Strains
Constitutive equation
Element stiffness matrix
Global stiffness matrix
Solution of the global stiffness equations

Course “Computational Geotechnics” 2/21

 three preliminary, (trivial yet) important remarks
R

1. while the FEM has been used in many fields of engineering practice
for over 30 years, it is only recently that it has begun to be widely
used for analyzing geotechnical problems. This is probably because
there are many complex issues which are specific to geotechnical
engineering and which have been resolved relatively recently

2. when properly used, this method can produce realistic results

which are of value to practical soil engineering problems

3. a good analysis, which simulates real behavior, allows the engineer

to understand problems better. While an important part of the
design process, analysis only provides the engineer with a tool to
quantify effects once material properties and loading conditions
have been set

Course “Computational Geotechnics” 3/21

 three remarks: short version
R

1. geotechnical engineering is complex. It is not because you’re

using the FEM that it becomes simpler

2. the quality of a tool is important, yet the quality of a result

also (mainly) depends on the user’s understanding of both
the problem and the tool

3. the design process involves considerably more than analysis

Course “Computational Geotechnics” 4/21

 introduction: the Finite Element Method
R

the FEM is a computational procedure that may be used to obtain an

approximate solution to a boundary value problem

the governing mathematical equations are approximated by a series

of algebraic equations involving quantities that are evaluated at
discrete points within the region of interest. The FE equations are
formulated and solved in such a way as to minimize the error in the
approximate solution

this lecture presents only a basic outline of the method

the discussion is restricted to:
 linear elasticity
 two-dimensional plane strain
attention is focused on the "displacement based" FE approach

Course “Computational Geotechnics” 5/21

 introduction: the Finite Element Method R

the first stage in any FE analysis is to generate a FE mesh

Footing
width = B

Node

Gauss point

a mesh consists of elements connected together at nodes

Course “Computational Geotechnics” 6/21

 examples: embankment R

Course “Computational Geotechnics” 7/21

 examples: multi-anchored diaphragm wall R

Course “Computational Geotechnics” 8/21

 the Finite Elements Method in one single slide! R

the nodes are the points where values Footing

width = B

of the primary variables (displacements)

are calculated
Node

the values of nodal displacements are

interpolated within the elements to give Gauss point

algebraic expressions for displacement

and strain throughout the complete mesh

a constitutive law is then used to relate strains to stresses and this

leads to the calculation of forces acting at the element nodes

the nodal forces are related to the nodal displacements by equations

which are set up and solved to find values of the nodal displacements

Course “Computational Geotechnics” 9/21

 the FEM involves the following steps R

Elements discretization
this is the process of modeling the geometry of the problem under
investigation by an assemblage of small regions, termed finite
elements. These elements have nodes defined on the element
boundaries, or within the elements

Primary variable approximation

a primary variable must be selected (e.g., displacements) and rules
as how it should vary over a finite element established. This
variation is expressed in terms of nodal values

 a polynomial form is assumed, where the order of the polynomial

depends on the number of nodes in the element
 the higher the number of nodes (the order of the polynomial), the
more accurate are the results (the longer takes the computation!)

Course “Computational Geotechnics” 10/21

 the FEM involves the following steps R

Element equations
use of an appropriate variational principle (e.g., minimum potential
energy) to derive element equations:

where is the element stiffness matrix, is the vector of nodal

displacements and is the vector of nodal forces
Global equations
combine element equations to form global equations

Boundary conditions
formulate boundary conditions and modify global equations. Loads
affect P, while displacements affect U
Solve the global equations
to obtain the displacements at the nodes. From nodal displacements
secondary quantities (stresses, strain) are evaluated

Course “Computational Geotechnics” 11/21

 displacement interpolation R

two-dimensional analysis of continua is generally based on the use of

either triangular or quadrilateral elements

the most used elements are based on an iso-parametric approach

Course “Computational Geotechnics” 12/21

 displacement interpolation R

primary unknowns: values of the nodal displacements

displacement within the element: expressed in terms of the nodal values

how does it work?  polynomial interpolation

shape functions

= nodal displacements

Course “Computational Geotechnics” 13/21

 illustration for the six-noded triangular element R

3 y quadratic interpolation

v x
6 5

1 2 12 coefficients, depending on the values of the 12

4 nodal displacements

strains may be derived within the element using the standard definitions

B-matrix

Course “Computational Geotechnics” 14/21

 constitutive relation (elasticity) R

elasticity: one-to-one relationship between stress and strain

in a FE context, stresses s and strains e are written in vector form
the stress-strain relationship is then expressed as: s = D e

material stiffness matrix

linear isotropic elasticity in plane strain

 
1  v v 0 
D
E  v 1 v 0 
(1  2v)(1  v)  1  2v 
 0 0 
 2 

in this case the coefficients of the matrix are constants, which means
that (for linear kinematics) the resulting F.E. equations are linear

Course “Computational Geotechnics” 15/21

 what happens with inelastic constitutive relations? R

advantage with elasticity: the coefficients of the matrix are constants,

the resulting F.E. equations are linear

the problem may be solved by applying all of the external loads

in a single calculation step

soils usually do not behave elastically

D e
s  D
with D depending on the current and past stress history

it is necessary to apply the external load in separate increments

and to adopt a suitable non-linear solution scheme

Course “Computational Geotechnics” 16/21

 element stiffness matrix R

body forces and surface tractions applied to the element may be

generalized into a set of forces acting at the nodes (vector of nodal forces)

 P1x 
nodal forces may be related  
to the nodal displacements by: 3  P1 y 
 
 P2 x 
K e Ue  Pe 6 5
e
 P2 y 
P  
P1x 1  
Ke element stiffness matrix 4
2  
P1y  
 P6 x 
Ke
 B T DBdv P 
 6y 

D material stiffness matrix

recall
B matrix relating nodal displacements to strains

Course “Computational Geotechnics” 17/21

 Gauss points R


Ke B T DBdv

to evaluate Ke, integration must be performed for each element

a numerical integration scheme must be employed (Gaussian integration)

essentially, the integral of a function is replaced by a weighted sum

of the function evaluated at a number of integration points

Course “Computational Geotechnics” 18/21

 global stiffness matrix (1) R

the stiffness matrix for the complete mesh is evaluated by combining

the individual element stiffness matrixes (assembly)

this produces a square matrix K of dimension equal to the number of

degrees-of-freedom in the mesh

the global vector of nodal forces P is obtained in a similar way by

assembling the element nodal force vectors

the assembled stiffness matrix and force vector are related by:

KU  P

where vector U contains the displacements at all the nodes in the mesh

Course “Computational Geotechnics” 19/21

 global stiffness matrix (2) R

if D is symmetric (elasticity), then Ke and hence K will be symmetric

the global stiffness matrix generally contains many terms that are zero
if the node numbering scheme is efficient then all of the non-zero
terms are clustered in a band along the leading diagonal

assembly
schemes for storage
solution

take into account its sym and

banded structure

number of dofs

Course “Computational Geotechnics” 20/21

 solution of the global stiffness equations R

once the global stiffness equations have been established

(and the boundary conditions added), they mathematically
form a large system of symultaneous (algebraic) equations

KU  P

these have to be solved to give values for the nodal displacements

it is advantageous to adopt special techniques to reduce

computation time (e.g. bandwidth and frontal techniques)

detailed discussion of such techniques is beyond the scope of

this lecture

Course “Computational Geotechnics” 21/21