Introduction

Civil engineering problems (stress analysis, seepage, flow in channels, soil-structure interaction) are governed by differential equations that rarely have closed-form solutions for real geometry and loading. Numerical techniques replace the continuous problem by a set of algebraic equations. The main techniques are FEM, FDM, BEM, DEM and mesh-free methods (SPH).

1. Finite Element Method (FEM)

The body is divided into elements joined at nodes; displacement inside an element is interpolated by shape functions and the equations $[K]\{U\}=\{F\}$ are formed from energy or virtual work (e.g. a beam with $[k]$ from Hermite functions).

• Advantages: handles irregular geometry, mixed materials, non-linearity, any boundary condition; systematic and general; large software base.
• Disadvantages: whole domain must be meshed; large equations; poor for infinite domains, cracks and large deformation (mesh distortion).
• Suitability: structures (frames, trusses, plates, shells), foundations, dams, soil stress analysis.

2. Finite Difference Method (FDM)

Derivatives are replaced by differences between grid-point values, e.g. $\dfrac{\partial u}{\partial x}\approx\dfrac{u_{i+1}-u_{i-1}}{2\Delta x}$.

• Advantages: simple concept and programming; efficient on regular grids; explicit and implicit time stepping.
• Disadvantages: difficult for curved or irregular boundaries; grid refinement is global; stability limits on $\Delta t$.
• Suitability: groundwater flow, seepage, open-channel and pipe flow, heat flow, flow nets.

3. Boundary Element Method (BEM)

The governing PDE is converted into an integral equation on the boundary only; only the boundary is discretised.

• Advantages: reduces dimension by one; small data; good for infinite or semi-infinite domains and stress concentration.
• Disadvantages: full, unsymmetric matrices; needs a fundamental solution; poor for non-linear and heterogeneous media.
• Suitability: soil-structure interaction, tunnels, fracture, potential flow.

4. Discrete Element Method (DEM)

Material is modelled as an assembly of separate blocks or particles; contact forces are found from force-displacement laws and Newton's second law is integrated in time.

• Advantages: models large displacement, sliding and separation naturally.
• Disadvantages: costly; contact parameters hard to calibrate.
• Suitability: granular soils, jointed rock, masonry, rockfall.

Comparison

FEM is the most widely used because it suits both solids and fluids and complex geometry.

Techniques of solving civil engineering problems

1. Analytical (exact, closed-form) methods
2. Numerical methods: finite element (FEM), finite difference (FDM), boundary element (BEM), discrete element (DEM), finite volume, mesh-free (SPH)
3. Experimental and model studies (laboratory and physical models)
4. Empirical and semi-empirical (code-based) methods

Finite element method

The body is divided into small finite elements connected at nodes. Inside each element the unknown (displacement $u$) is assumed as

$$u = [N]\{d\}$$

where $[N]$ are shape functions and $\{d\}$ are nodal values. Strain is $\{\varepsilon\}=[B]\{d\}$ and stress $\{\sigma\}=[D]\{\varepsilon\}$. Using minimum potential energy or virtual work, the element stiffness is

$$[k]=\int [B]^T[D][B]\,dV$$

Element matrices are assembled into $[K]\{U\}=\{F\}$; supports are applied; the equations are solved for nodal displacements; strains and stresses follow from $[B]$ and $[D]$.

Finite difference method

The domain is covered with a regular grid of points spacing $\Delta x$, $\Delta y$. Derivatives in the governing differential equation are replaced by differences of nodal values, for example

$$\frac{d^2u}{dx^2}\approx\frac{u_{i-1}-2u_i+u_{i+1}}{\Delta x^2}$$

Writing this at every grid point gives one algebraic equation per point. Boundary values are inserted, and the simultaneous equations are solved (directly or iteratively). For time-dependent problems the time derivative is also differenced, using an explicit or implicit scheme.

 FEM: elements + nodes      FDM: grid points
 o----o----o                 o  o  o  o
 | \  |  / |                 o  o  o  o
 o----o----o                 o  o  o  o

Necessity of computational techniques

• Real structures and soils have complex geometry, loading and material behaviour, so closed-form solutions do not exist.
• Classical hand methods (slope-deflection, moment distribution) become impractical for large indeterminate structures with hundreds of unknowns.
• Computers solve thousands of equations quickly and accurately, allowing many load cases and design iterations.
• Non-linear, dynamic, time-dependent (seepage, consolidation, flood routing) and soil-structure problems can only be solved numerically.
• They reduce the cost of physical model tests and allow parametric studies and safer, more economical design.

Algorithm of the finite element method

1. Discretisation: divide the body into finite elements joined at nodes (select element type and mesh).
2. Select interpolation (shape) functions: assume displacement within an element as $u = [N]\{d\}$.
3. Strain-displacement and stress-strain relations: $\{\varepsilon\} = [B]\{d\}$, $\{\sigma\} = [D]\{\varepsilon\}$.
4. Element stiffness matrix and load vector: $[k] = \int [B]^T[D][B]\,dV$ and equivalent nodal loads $\{f\}$.
5. Assembly: combine element matrices into the global equations $[K]\{U\} = \{F\}$.
6. Boundary conditions: apply supports (known displacements) and loads.
7. Solve for nodal displacements.
8. Post-processing: compute strains, stresses, reactions and plot results.

Advantages of FEM

• Handles complex geometry, irregular boundaries and any support condition.
• Handles different materials, non-homogeneous and non-linear behaviour in the same model.
• Mesh can be refined where stress changes rapidly.
• Applicable to many fields: structures, soil, seepage, heat and fluid flow.
• Algorithm is systematic and easy to program; wide software availability.

Numerical computation in civil engineering grew with the available tools:

• Before 1900: hand analysis; Euler, Lagrange and Navier created the analytical basis. Graphical and relaxation ideas appeared.
• Early 20th century: Richardson (1910) introduced finite differences; Southwell's relaxation method and Cross's moment distribution (1930) allowed hand analysis of frames.
• 1940s: Hrennikoff and Courant proposed lattice and variational discretisation ideas; electronic computers appeared.
• 1950s-60s: matrix methods of structural analysis; Turner, Clough, Martin and Topp (1956) developed the direct stiffness method for aircraft structures; Clough named the "finite element method" (1960).
• 1960s-70s: Zienkiewicz and Cheung extended FEM to general field problems; boundary element method (Brebbia) and Cundall's discrete element method (1971) were developed.
• 1980s onward: personal computers, graphical pre- and post-processors, commercial packages (SAP, STAAD, ANSYS, ABAQUS, PLAXIS, MODFLOW) made computational analysis routine in design offices.
• Today: non-linear, dynamic, 3D and mesh-free analysis (SPH), parallel computing.

The boundary element method (BEM) is a numerical method in which the governing differential equation is converted into an integral equation over the boundary of the domain, using a fundamental solution (Green's function). Only the boundary is divided into elements, so a 3D problem becomes a 2D surface problem.

Procedure

1. Write the governing PDE (e.g. Laplace or Navier equation) and its weighted-residual form.
2. Apply Green's theorem to move all derivatives to the weight function, which is the fundamental solution.
3. Obtain the boundary integral equation, then discretise the boundary into elements.
4. Evaluate integrals and form $[H]\{u\}=[G]\{t\}$, where $u$ is boundary displacement and $t$ traction.
5. Apply boundary conditions, solve for the unknown boundary values, and then compute interior values.

Advantages: fewer unknowns and easy data preparation; accurate for stress concentration; suited to infinite and semi-infinite domains.

Disadvantages: matrices are full and unsymmetric; a fundamental solution is needed; non-linear and non-homogeneous problems are difficult.

Applications: soil-structure interaction, tunnels, foundations, fracture mechanics, potential flow and seepage.

The discrete element method (DEM), introduced by Cundall (1971), models a material as an assembly of separate rigid or deformable blocks or particles that interact at contacts. It is suited to discontinua where displacement and rotation are large and contacts open and close.

Working

1. Define the particles/blocks and their initial positions.
2. Detect contacts between neighbours.
3. Compute contact forces from a force-displacement law, for example normal force $F_n=k_n\,\delta_n$ and a friction-limited shear force.
4. Sum forces on each particle and apply Newton's second law $m\ddot{x}=\sum F$.
5. Integrate in time with an explicit scheme, update positions and repeat (the time step must be small).

Advantages: handles fracture, sliding, rotation and large displacements naturally; no continuity assumption.

Disadvantages: computationally expensive; contact parameters (stiffness, friction, damping) need calibration.

Applications: granular soils and silo flow, jointed rock slopes, masonry and rockfill, landslides and rockfall.

Smoothed particle hydrodynamics (SPH) is a mesh-free, Lagrangian particle method originally developed for astrophysics (Lucy, Gingold and Monaghan, 1977) and now used for fluids and solids with large deformation. The medium is represented by particles carrying mass, density, velocity and pressure; there is no grid.

Principle: value of a field $f$ at a point is the weighted sum of values at neighbouring particles inside a smoothing length $h$:

$$f(x_i)\approx\sum_j \frac{m_j}{\rho_j}\,f_j\,W(x_i-x_j,h)$$

where $W$ is the smoothing (kernel) function. Gradients are obtained by differentiating $W$, so the governing equations (continuity and momentum) become ordinary differential equations for each particle and are integrated in time.

Advantages: no mesh distortion; free surfaces and large deformation are handled easily; simple conservation of mass.

Disadvantages: costly neighbour search; boundary treatment and tensile instability are difficult; lower accuracy than FEM for small deformation.

Applications: dam break, wave impact, flood flow, landslide and debris flow, soil flow.

The finite element method replaces a continuum by a finite number of elements and solves the resulting algebraic equations. The basic steps are:

1. Discretisation: divide the body into finite elements joined at nodes (select element type and mesh).
2. Select interpolation (shape) functions: assume displacement within an element as $u = [N]\{d\}$.
3. Strain-displacement and stress-strain relations: $\{\varepsilon\} = [B]\{d\}$, $\{\sigma\} = [D]\{\varepsilon\}$.
4. Element stiffness matrix and load vector: $[k] = \int [B]^T[D][B]\,dV$ and equivalent nodal loads $\{f\}$.
5. Assembly: combine element matrices into the global equations $[K]\{U\} = \{F\}$.
6. Boundary conditions: apply supports (known displacements) and loads.
7. Solve for nodal displacements.
8. Post-processing: compute strains, stresses, reactions and plot results.

 Body -> Mesh -> [k] each -> [K]{U}={F} -> U -> stress

FEM is applied to any problem governed by a differential equation with boundary conditions. The main types are:

1. Structural / solid mechanics (equilibrium problems): stress, strain and displacement in trusses, frames, plates, shells, dams, retaining walls and foundations. This includes linear static analysis, plane stress, plane strain and axisymmetric problems.
2. Eigenvalue problems: natural frequencies and mode shapes (vibration) and buckling loads of structures, giving $[K]\{\phi\}=\lambda[M]\{\phi\}$.
3. Propagation (transient) problems: time-dependent behaviour such as dynamic response to earthquake and wind, consolidation of soil, heat conduction with time.
4. Field problems: steady seepage through dams and soils, heat conduction, electric potential, torsion of non-circular shafts governed by Laplace/Poisson equations.
5. Fluid flow: open channel, pipe and groundwater flow.
6. Non-linear problems: material non-linearity (plasticity, cracking of concrete), geometric non-linearity (large displacement) and contact.
7. Soil-structure interaction and construction-stage analysis.

Examples: stress in a gravity dam, seepage below a weir, settlement under a footing, modal analysis of a bridge.

Computer programs let engineers solve FEM and FDM problems on realistic models. They normally contain a pre-processor (geometry, mesh, loads), a solver and a post-processor (results display).

FEM software

• SAP2000, ETABS, STAAD.Pro: frames, trusses, buildings and bridges (structural analysis and design).
• ANSYS, ABAQUS, NASTRAN, LS-DYNA: general-purpose, linear and non-linear stress, thermal and dynamic analysis.
• PLAXIS, GeoStudio (SIGMA/W, SEEP/W), Geo5: geotechnical, slopes, excavations and seepage.
• Midas Civil, CSiBridge, SAFE: bridges and slabs/foundations.
• OpenSees, FreeFEM, MATLAB (own coding): research and teaching.

FDM software

• MODFLOW (USGS): groundwater flow, using a finite difference grid.
• HEC-RAS, MIKE 11/MIKE 21, SWMM, EPANET: river, flood and pipe network flow (finite difference/implicit schemes).
• FLAC (Itasca): explicit finite difference for soil and rock mechanics.
• MATLAB/Python scripts for heat conduction, kinematic wave, unsteady flow.

In each package the user builds a model, assigns materials and boundary conditions, runs the analysis and checks results (displacement, stress, head, water level).

Discretisation means dividing a continuous body or domain into a finite number of smaller parts (elements, or grid points) so that a problem with infinite unknowns becomes a problem with a finite number of unknowns, usually at nodes. The continuous differential equation is then replaced by a set of algebraic equations.

Example (FEM): a rectangular plate with a hole under tension is divided into triangular and quadrilateral elements joined at nodes. Smaller elements are used near the hole where stress changes rapidly.

 +--+--+--+--+--+
 |\ |\ |  |  |  |
 +--+--o--+--+--+     o = hole
 |  |/ |\ |  |  |
 +--+--+--+--+--+

Unknown displacements are found at the nodes only; displacement inside an element is obtained from shape functions.

Example (FDM): for a bar $u''=f(x)$ the length is divided into points $x_0,x_1,\dots,x_n$ with spacing $h$ and the equation is written at each point as

$$\frac{u_{i-1}-2u_i+u_{i+1}}{h^2}=f_i$$

A finer discretisation gives better accuracy but needs more computation.

Finite element method (concept)

The structure is divided into elements joined at nodes. The displacement in an element is expressed by shape functions, $u=[N]\{d\}$, and the element stiffness $[k]=\int [B]^T[D][B]\,dV$ is formed using minimum potential energy. Assembling gives $[K]\{U\}=\{F\}$.

• Applications: stress analysis of trusses, frames, plates, dams, foundations, seepage, vibration.
• Advantages: complex geometry and boundaries; different materials; non-linearity; local refinement.
• Disadvantages: whole domain must be meshed; large computing time; poor for infinite domains and very large deformation.

Finite difference method (concept)

Derivatives of the governing equation are replaced by differences of values at regular grid points, e.g. $\dfrac{\partial^2 h}{\partial x^2}\approx\dfrac{h_{i-1}-2h_i+h_{i+1}}{\Delta x^2}$, giving one algebraic equation per point.

• Applications: groundwater and seepage (flow net grids), unsteady open-channel and pipe flow, heat conduction, kinematic wave routing.
• Advantages: simple, easy to program; good for regular domains; both explicit and implicit schemes.
• Disadvantages: irregular or curved boundaries are difficult; stability limits for explicit schemes; mesh is not locally refined easily.

Comparison with other methods

FEM is the most flexible; FDM is the simplest for regular grids and fluid problems.

Computational techniques

• Finite element method (FEM)
• Finite difference method (FDM)
• Boundary element method (BEM)
• Discrete element method (DEM)
• Finite volume method and mesh-free methods (SPH)

Why FEM predominates

• Handles any shape of boundary and any combination of loads and supports.
• Handles non-homogeneous, anisotropic and non-linear materials.
• Same formulation covers structures, soil, seepage, heat and dynamics.
• Mesh can be refined locally and element types can be mixed.
• Stiffness matrices are sparse and symmetric, so solution is efficient.
• Many mature commercial programs and a strong mathematical basis (energy principles) exist.

Steps in FEM

1. Discretisation: divide the body into finite elements joined at nodes (select element type and mesh).
2. Select interpolation (shape) functions: assume displacement within an element as $u = [N]\{d\}$.
3. Strain-displacement and stress-strain relations: $\{\varepsilon\} = [B]\{d\}$, $\{\sigma\} = [D]\{\varepsilon\}$.
4. Element stiffness matrix and load vector: $[k] = \int [B]^T[D][B]\,dV$ and equivalent nodal loads $\{f\}$.
5. Assembly: combine element matrices into the global equations $[K]\{U\} = \{F\}$.
6. Boundary conditions: apply supports (known displacements) and loads.
7. Solve for nodal displacements.
8. Post-processing: compute strains, stresses, reactions and plot results.

Foundation of the finite element method

FEM rests on the idea that a complex continuum can be approximated by a collection of simple pieces (elements) whose behaviour is known. Its mathematical basis is:

1. Variational / energy principle: the correct displacement field makes the total potential energy $\Pi=U-W$ stationary (minimum), or equivalently satisfies the virtual work equation. Weighted residual (Galerkin) methods can also be used.
2. Piecewise approximation: displacement in each element is $u=[N]\{d\}$ with shape functions that are continuous across elements.
3. Element equations: from $\{\varepsilon\}=[B]\{d\}$ and $\{\sigma\}=[D]\{\varepsilon\}$, $[k]=\int [B]^T[D][B]\,dV$.
4. Assembly and solution: compatibility at nodes and equilibrium give $[K]\{U\}=\{F\}$.

Why FEM is less suitable for large deformation

• Elements are attached to the material (Lagrangian); under very large strain the mesh distorts, elements become tangled or negative in volume, and accuracy falls.
• Small-strain formulation $[B]$ and constant $[K]$ assume small displacement; for large deformation the stiffness changes with geometry and needs non-linear incremental solution and re-meshing.
• Material separation, fracture and flow are not represented by continuous elements. Methods such as DEM, SPH or ALE formulations are better.

Choosing a numerical method

Selection depends on geometry, material behaviour (linear or non-linear, continuum or discontinuum), boundary conditions, required accuracy and computing resources.