Return to computing page for the first course APMA0330
Return to computing page for the second course APMA0340
Return to computing page for the fourth course APMA0360
Return to Mathematica tutorial for the first course APMA0330
Return to Mathematica tutorial for the second course APMA0340
Return to Mathematica tutorial for the fourth course APMA0360
Return to the main page for the course APMA0330
Return to the main page for the course APMA0340
Return to the main page for the course APMA0360
Introduction to Linear Algebra with Mathematica

Glossary

Linear Equations of Elastic Media

The Lamé equations of elastodynamics (also called the Navier–Lamé equations) are the governing partial differential equations for the displacement field in a linear, homogeneous, isotropic elastic solid. They follow from combining Hooke's law with the equations of mechanical equilibrium:

\[ \mu\,\nabla^2 \mathbf{u} + \left( \lambda + \mu \right) \nabla \left( \nabla \cdot \mathbf{u} \right) + \mathbf{f} = \rho\,\ddot{\bf u} , \]
where If we apply the Helmholtz decomposition
\[ \mathbf{u} = \nabla \phi + \nabla \times \Psi , \qquad \nabla \cdot \Psi = 0, \]
the Lamé vector equation is separated into longitudinal and transverse components. In dynamics, these components correspond to the familiar P-waves and S-waves, with
\[ c_p = +\sqrt{\frac{\lambda + 2\mu}{\rho}} , \qquad c_s = \sqrt{\frac{\mu}{\mu}} \]
Material     cp (km/s)     cs (km/s)
Steel     ≈ 5.8     ≈ 3.1
Copper     ≈ 4.76     ≈ 2.28

Substituting the speed (longitudinal and transverse) definitions gives
\[ \frac{\partial^2 \mathbf{u}}{\partial t^2} = \left( c_p^2 - c_s^2 \right) \nabla \left( \nabla \cdot \mathbf{u} \right) + c_s^2 \nabla^2 \mathbf{u} , \]
or, using the vector identity
\[ \Delta \mathbf{u} = \nabla^2 \mathbf{u} = \nabla \left( \nabla \cdot \mathbf{u} \right) - \nabla \times \nabla \times \mathbf{u} , \]
we get
\[ \frac{\partial^2 \mathbf{u}}{\partial t^2} = c_p^2 \, \nabla \left( \nabla \cdot \mathbf{u} \right) - c_s^2 \,\nabla \times \left( \nabla \times \mathbf{u} \right) . \]
In component form:
\[ \ddot{u}_i = c_s^2 \,\nabla^2 u_i + \left( c_p^2 - c_s^2 \right) \frac{\partial}{\partial x_i} \left( \nabla \cdot \mathbf{u} \right) , \]
or explicitly in 3D,
\[ \ddot{u}_x = c_p^2 \,u_{x,xx} + c_s^2 \left( u_{x,yy} + u_{x,zz} \right) + \left( c_p^2 - c_s^2 \right) \left( u_{y,xy} + u_{z,xz} \right) , \]
with analogous equations for uy​, uz​.

The physical interpretation is especially clean: cp controls compressional (P) waves, while cs controls shear (S) waves. In a homogeneous isotropic medium, the equation separates naturally into these two wave modes.

For static elasticity with no body forces, the Lamé vector equation reduces to

\begin{equation} \label{EqKM.5} \left( c_p^2 - c_s^2 \right) \nabla \left( \nabla \cdot \mathbf{u} \right) + c_s^2 \nabla^2 \mathbf{u} = 0 . \end{equation}
The principal symbol of the Lamé operator is
\[ \sigma (L) = -\left( \lambda + \mu \right) \begin{pmatrix} \xi_1^2 & \xi_1 \xi_2 & \xi_1 \xi_3 \\ \xi_1 \xi_2 & \xi_2^2 & \xi_2 \xi_3 \\ \xi_1 \xi_3 & \xi_2 \xi_3 & \xi_3^2 \end{pmatrix} - \mu \left( \xi_1^2 + \xi_2^2 + \xi_3^2 \right)\begin{pmatrix} 1 &0&0 \\ 0&1 &0 \\ 0&0& 1 \end{pmatrix} . \]
We rewrite this symbol in concise way:
\[ \sigma (L) = \mu\, |\xi |^2 \mathbf{I} + \left( \lambda + \mu \right) \xi \otimes \xi . \]

The basic equations of elastostatics consist of equilibrium equations of stresses, strain-displacement relations, and Hooke’s law that relates stresses and strains. In plane elasticity (plane strain and plane stress), the equilibrium equations are (body forces are absent)

\begin{equation} \label{EqKM.1} \frac{\partial \sigma_{xx}}{\partial x} + \frac{\partial \sigma_{xy}}{\partial y} + \frac{\partial \sigma_{xz}}{\partial z} = 0, \qquad \frac{\partial \sigma_{xy}}{\partial x} + \frac{\partial \sigma_{yy}}{\partial y} + \frac{\partial \sigma_{yz}}{\partial z} = 0 , \qquad \frac{\partial \sigma_{xz}}{\partial x} + \frac{\partial \sigma_{yz}}{\partial y} + \frac{\partial \sigma_{zz}}{\partial z} = 0 , \end{equation}
where σₓₓ, σᵧᵧ, and σₓᵧ are stresses and (x, y, z) are Cartesian coordinates.

Special Case: Plane Stress vs. Plane Strain

Depending on whether your 2D case represents a thin plate (Plane Stress) or a long cylinder (Plane Strain), Lamé's first parameter λ must be adjusted using Young's Modulus (E) and Poisson's ratio (ν):

Parameter Plane Strain Plane Stress
μ (Shear Modulus)    $$\frac{E}{2(1+\nu)}$$ $$\frac{E}{2(1+\nu)}$$
λ (Lamé's First)    $$\frac{E\nu}{(1+\nu)(1-2\nu)}$$ $$\frac{E\nu}{1-\nu^2}$$

The basic equations of elastostatics consist of equilibrium equations of stresses, strain-displacement relations, and Hooke’s law that relates stresses and strains. In plane elasticity (plane strain and plane stress), the equilibrium equations are (body forces are absent)

\begin{equation} \label{EqKM.1} \frac{\partial \sigma_{xx}}{\partial x} + \frac{\partial \sigma_{xy}}{\partial y} + \frac{\partial \sigma_{xz}}{\partial z} = 0, \qquad \frac{\partial \sigma_{xy}}{\partial x} + \frac{\partial \sigma_{yy}}{\partial y} + \frac{\partial \sigma_{yz}}{\partial z} = 0 , \qquad \frac{\partial \sigma_{xz}}{\partial x} + \frac{\partial \sigma_{yz}}{\partial y} + \frac{\partial \sigma_{zz}}{\partial z} = 0 , \end{equation}
where σₓₓ, σᵧᵧ, and σₓᵧ are stresses and (x, y, z) are Cartesian coordinates. The strains and displacements are related by
\begin{equation} \label{EqKM.2} \varepsilon_{xx} = \frac{\partial u_x}{\partial x} , \quad \varepsilon_{yy} = \frac{\partial u_y}{\partial y} , \quad \varepsilon_{xy} = \frac{1}{2} \left( \frac{\partial u_x}{\partial y} + \frac{\partial u_y}{\partial x} \right) , \end{equation}
where εₓₓ, εᵧᵧ, and εₓᵧ are tensorial strain components, and uₓ and uᵧ are displacements in x and y directions, respectively. The stress–strain relations are given by
\begin{equation} \label{EqKM.3} \begin{split} \sigma_{xx} &= \lambda^{\ast} \left( e_{xx} + e_{yy} \right) + 2\mu\,e_{xx} , \\ \sigma_{yy} &= \lambda^{\ast} \left( e_{xx} + e_{yy} \right) + 2\mu\,e_{yy} , \\ \sigma_{xy} &= 2\mu\, e_{xy} , \end{split} \end{equation}
or inversely
\begin{equation} \label{EqKM.4} \begin{split} \varepsilon_{xx} &= \frac{1}{2\mu} \left[ \sigma_{xx} - \frac{\lambda^{\ast}}{2\left( \lambda^{\ast} + \mu \right)} \left( \sigma_{xx} + \sigma_{yy} \right) \right] , \\ \varepsilon_{yy} &= \frac{1}{2\mu} \left[ \sigma_{yy} - \frac{\lambda^{\ast}}{2\left( \lambda^{\ast} + \mu \right)} \left( \sigma_{xx} + \sigma_{yy} \right) \right] , \\ \varepsilon_{xy} &= \frac{1}{2\mu} \,\sigma_{xy} , \end{split} \end{equation}
where μ is the Lamé constant (the shear modulus) and
\[ \lambda^{\ast} = \frac{3 - \kappa}{\kappa -1} , \qquad \kappa = \begin{cases} 3 - 4\nu \qquad&\mbox{for plain strain}, \\ \frac{3 - \nu}{1 + \nu} \qquad&\mbox{for plain stress}. \end{cases} \]
Here ν is Poisson's ratio. The compatibility equation of strains can be obtained from Eq.\eqref{EqKM.2} by eliminating the displacements as follows:
\[ \frac{\partial^2 e_{xx}}{\partial y^2} + \frac{\partial^2 e_{yy}}{\partial y^2} = 2\,\frac{\partial^2 e_{xy}}{\partial x\,\partial y} . \]
By using the stress–strain relations Eq.\eqref{EqKM.4} together with the equations of equilibrium Eq.\eqref{EqKM.1}, the compatibility condition can be expressed in terms of stresses as
\[ \nabla^2 \left( \sigma_{xx} + \sigma_{yy} \right) = 0 \qquad \iff \qquad \Delta \left( \sigma_{xx} + \sigma_{yy} \right) = 0 , \]
where Δ is Laplace's operator.

In a 2D polar coordinate system (r, θ), the displacement vector is split into a radial displacement ur and a circumferential (tangential) displacement uθ. Then 2D Lamé equations in polar coordinates become a system of two coupled partial differential equations:

\[ \begin{split} c_p^2 \,\frac{\partial e}{\partial r} - \frac{2\,c_s^2}{r}\,\frac{\partial \omega_z}{\partial \theta} + f_r &= 0 , \\ c_p^2\,\frac{1}{r}\,\frac{\partial e}{\partial\theta} + 2\,c_s^2 \,\frac{\partial \omega_z}{\partial \theta} + f_{\theta} &= 0 . \end{split} \]
To fully expand these equations, we define the dilatation (e, representing volumetric strain) and the rotation (ωz, representing rigid-body rotation around the z-axis):
  1. Volumetric Strain / Dilatation (e)
    \[ e = \varepsilon_{rr} + \varepsilon_{\theta\theta} = \frac{\partial u_r}{\partial r} + \frac{u_r}{r} + \frac{1}{r}\,\frac{\partial u_{\theta}}{\partial\theta} . \]
  2. Rotation component (ωz)
    \[ \omega_z = \frac{1}{2} \left[ \frac{1}{r}\,\frac{\partial \left( r\,u_{\theta} \right)}{} - \frac{1}{r}\,\frac{\partial u_r}{\partial\theta} \right] = \frac{1}{2} \left[ \frac{\partial u_{\theta}}{\partial r} + \frac{u_{\theta}}{r} - \frac{1}{r}\,\frac{\partial u_r}{\partial\theta} \right] \]
  3. Radial strain (εᵣᵣ): \( \displaystyle \quad \varepsilon_{rr} = \frac{\partial u_r}{\partial r} . \)
  4. Circumferential strain (εθθ): \( \displaystyle \quad \varepsilon_{\theta\theta} = \frac{u_r}{r} + \frac{1}{r}\,\frac{\partial u_{\theta}}{\partial\theta} . \)
  5. Shear strain (γrθ): \( \displaystyle \quad \gamma_{r\theta} = \frac{1}{r}\,\frac{\partial u_r}{\partial\theta} + \frac{\partial u_{\theta}}{\partial r} - \frac{u_{\theta}}{r} . \)
  6. Material and external constants:
    • λ, μ: are the Lamé Parameters (μ is also the shear modulus G).
    • fr, fθ are components of the body force.
Note that 2D polar Laplacian is \( \displaystyle \quad \Delta = \nabla^2 = \frac{\partial^2}{\partial r^2} + \frac{1}{r}\,\frac{\partial}{\partial r} + \frac{1}{r^2}\,\frac{\partial^2}{\partial \theta^2} . \) The 2D Lamé system in polar coordinates:
\[ \begin{split} \mu \left( \nabla^2 u_r - \frac{u_r}{r^2} - \frac{2}{r^2}\,\frac{\partial u_{\theta}}{\partial\theta} \right) + \left( \lambda + \mu \right) \frac{\partial}{\partial r} \left( \frac{\partial u_r}{\partial r} + \frac{u_r}{r} + \frac{1}{r}\,\frac{\partial u_{\theta}}{\partial\theta} \right) + f_r &= 0, \\ \mu \left( \nabla^2 u_{\theta} - \frac{u_{\theta}}{r^2} + \frac{2}{r^2}\,\frac{\partial u_{r}}{\partial\theta} \right) + \left( \lambda + \mu \right) \frac{\partial}{\partial \theta} \left( \frac{\partial u_r}{\partial r} + \frac{u_r}{r} + \frac{1}{r}\,\frac{\partial u_{\theta}}{\partial\theta} \right) + f_{\theta} &= 0 . \end{split} \]

Kolosov–Muskhelishvili formalism

Among various mathematical methods in plane elasticity, the complex potential function method by Kolosov and Muskhelishvili are one of the powerful and convenient methods to treat two-dimensional problems. In the complex poten- tial method, stresses and displacements are expressed in terms of analytic functions of complex variables. The problem of obtaining stresses and displacements around a crack tip is converted to finding some analytic functions subjected to appropriate boundary conditions. A brief introduction of the general formulation of the Kolosov and Muskhelishvili complex potentials is given in this section.

Two‑dimensional isotropic elasticity possesses a remarkable analytic structure: the displacement field can be expressed entirely in terms of two analytic functions. This reduction is unique to 2D and has no analogue in 3D.

The formalism originates from the works of Gury Kolosov (1909) and Nikoloz Muskhelishvili (1930–1953) and forms the backbone of complex variable methods in elasticity, fracture mechanics, and wedge singularity theory. Kolosov introduced key ideas connecting plane elasticity problems with analytic functions of a complex variable, while Muskhelishvili systematically developed and generalized the method, particularly through his influential 1910s–1950s work on boundary-value problems. The resulting formulation expresses stresses, displacements, and boundary conditions in terms of a small number of analytic complex potentials.

The formalism became one of the foundational methods of plane elasticity, especially for problems involving holes, cracks, inclusions, and other geometrically complicated boundaries. Its importance grew substantially through Muskhelishvili's monograph Some Basic Problems of the Mathematical Theory of Elasticity, which helped establish complex-variable methods as a powerful alternative to direct solution of the governing partial differential equations. The same mathematical framework later became closely associated with classical analyses of stress concentrations and fracture mechanics.

In a Cartesian coordinate system (x, y), the complex variable z and its conjugate z✶ are defined as

\[ z = x + \mathbf{j}\,y , \qquad z^{\ast} = \overline{z} = x - \mathbf{j}\,y , \]
respectively, where j or ⅉ is the imaginary unit on complex plane ℂ, so ⅉ² = −1. They can also be expressed in polar coordinates (r, θ) as
\[ z = r\left( \cos\theta + \mathbf{j}\,\sin \theta \right) = r\,e^{\mathbf{j}\theta} , \qquad z^{\ast} = r\left( \cos\theta - \mathbf{j}\,\sin \theta \right) = r\,e^{-\mathbf{j}\theta} , \]
respectively. Consider a function, f(z), of the complex variable z. The derivative of f(z) with respect to z is by definition
\[ f' (z) = \frac{{\text d}f}{{\text d}z} = \lim_{\Delta z \to 0} \ \frac{f(z + \Delta z) - f(z)}{\Delta z} . \]
If f(z) has a derivative at point z₁ and also at each point in some neighborhood of z₁, then f(z) is said to be holomorphic (sometime analytic) at z₀. The complex function f(z) can be expressed in the form
\[ f(z) = u(x,y) + \mathbf{j}\,v(x,y) \]
where u and v are real-valued functions. If f(z) is holomorphic, we have
\[ \frac{\partial}{\partial x}\,f(z) = f'(z)\,\frac{\partial z}{\partial x} = f' (z) \]
and
\[ \frac{\partial}{\partial y}\,f(z) = f'(z)\,\frac{\partial z}{\partial y} = \mathbf{j}\,f' (z) \]
where a prime stands for differentiation with respect to z. Thus,
\[ \frac{\partial}{\partial x}\,f(z) = -\mathbf{j}\,\frac{\partial}{\partial y}\,f(z) \]
or
\[ \frac{\partial u}{\partial x} = + \mathbf{j}\,\frac{\partial v}{\partial x} = \frac{\partial v}{\partial y} - \mathbf{j}\,\frac{\partial u}{\partial y} . \]
From this equation, we obtain the Cauchy-Riemann equations:
\[ \frac{\partial u}{\partial x} = \frac{\partial v}{\partial y} , \qquad \frac{\partial u}{\partial y} = - \frac{\partial v}{\partial x} . \]
These equations can also be shown to be sufficient for f(z) to be holomorphic.

From the Cauchy-Riemann equations it is easy to derive the following:

\[ \nabla^2 u = \nabla^2 v = 0 , \]
that is, the real and imaginary parts of an holomorphic function are harmonic.

The Kolosov–Muskhelishvili (KM for short) representation expresses stresses and displacements in terms of two analytic functions, Φ(z) and Ψ(z). Namely, we have stress representation

\begin{align*} \sigma_{xx} + \sigma_{yy} &= 4\,\Re\left9 \Phi' (z) \right] , \\ \sigma_{yy} - \sigma_{xx} + 2\mathbf{j}\, \sigma_{xy} &= 2 \left[ z\,\Phi'' (z) + \Psi' (z) \right] , \end{align*}
and displacement representation
\[ 2\mu \left( u_x + \mathbf{j}\,u_y \right) = \kappa\,\Phi (z) -z\,\Phi' (z) - \Psi (z) , \]
where \( \displaystyle \ \kappa = \frac{3-\nu}{1+\nu} \ \) (plane stress) or κ = 3 − 4ν (plane strain).

Traction boundary conditions become

\[ \sigma_{n} + \mathbf{j}\,\tau_{n} = 2\.\Re \left\{ \Phi' (z) \right\} + e^{-\mathbf{j}\theta} \left[ z\,\Phi'' (z) + \Psi' (z) \right] . \]

In polar coordinates z = rejθ, we have

\begin{align*} \sigma_{rr} + \sigma_{\theta\theta} &= 4\,\Re \left[ \Phi' (z) \right] , \\ \sigma_{\theta\theta} - \sigma_{rr} + 2\mathbf{j}\,\sigma_{r\theta} &= 2\,e^{-\mathbf{j}\theta} \left[ z\,\Phi'' (z) + \Psi' (z) \right] , \end{align*}
and displacement
\[ 2\mu \left( u_r + \mathbf{j}\.u_{\theta} \right) = e^{-\mathbf{j}\theta} \left[ \kappa\,\Phi (z) - z\,\Phi' (z) - \Psi (z) \right] . \]

Boundary conditions on a circular boundary r = R are

\[ \sigma_{rr} + \mathbf{j}\,\sigma_{r\theta} = 2\,\Re \left[ \Phi' (z) \right] + e^{-\mathbf{j}\theta} \left[ z\,\Phi'' (z) + \Psi' (z) \right] . \]
This form of boundary conditions are used in

KM representations in elliptic coordinates μ, ν (confocal ellipses).

\[ x = c\,\cosh\mu \,\cos\nu , \qquad y = c\,\sinh\mu\,\sin\nu . \]
Then complex coordinates become
\[ z = x + \mathbf{j}\,y = c\left( \cosh\mu \,\cos\nu + \mathbf{j}\,\sinh\mu\,\sin\nu \right) . \]
The KM formulas in elliptic coordinates remain formally identical, but the conversion to elliptic components requires the metric coefficients:
\[ h_{\mu} = h_{\nu} = c\,\sqrt{\sinh^2 \mu + \sin^2 \nu} . \]
Thus,
\[ \sigma_{\mu\mu} = \sigma_{xx} \left( \partial_{\mu} x \right)^2 + \sigma_{yy} \left( \partial_{\mu} y \right)^2 + 2\,\sigma_{xy} \left( \partial_{\mu} x \right) \left( \partial_{\mu} y \right) , \]
etc.

Boundary conditions on a confocal ellipse μ = μ₀:

\[ \sigma_{\mu\mu} + \mathbf{j}\,\sigma_{\mu\nu} = \left( \sigma_{xx} + \mathbf{j}\,\sigma_{xy} \right) \left( \partial_{\mu} x + \mathbf{j}\,\partial_{\mu} y \right) + \left(\sigma_{xy} + \mathbf{j}\,\sigma_{yy} \right) \left( \partial_{\mu} x - \mathbf{j}\,\partial_{\mu} y \right) . \]
This form of boundary conditions is used for    
Example 1: We consider the domain Ω and its boundary ∂Ω = Γ, so
\[ \Omega = \left\{ (x,y) \in \mathbb{R}^2 \ : \ y > 0 \right\} , \qquad \Gamma = = \partial\Omega = \left\{ y =0 \right\} . \]
Suppose that on the boundary y = 0 the traction is prescribed:
\[ \sigma_{yy} (x,0) = p(x) , \qquad \sigma_{xy} (x,0) = q(x) . \]
Define the complex traction
\[ T(x) \ := \sigma_{yy} (x,0) + \mathbf{j}\,\sigma_{xy} (x,0) = p(x) + \mathbf{j}\,q(x) . \]
On y = 0, the KM formulas give a relation of the form
\[ T(x) = 2\,\Re \left[ \Phi' (x) \right] + x\,\Phi'' (x) + \Psi' (x) , \]
where Φ, Ψ are boundary values of analytic functions in the upper half-plane. Now KM actually asks to determine these holomorphic potentials Φ and Ψ in half-plane y > 0. This could be achieved via Cauchy integrals or Hilbert transform. However, we take one Fourier mode:
\[ p(x) = \Re \left[ P\, e^{\mathbf{j}kx} \right] , \qquad q(x) = \Re \left[ Q\, e^{\mathbf{j}kx} \right] , \qquad k > 0 . \]
Then we seek potentials in the form: .
\[ \Phi (z) = A\, e^{\mathbf{j}kz} , \qquad \Psi (z) = B\, e^{\mathbf{j}kz} , \]
with z = x + ⅉy, so \( \displaystyle \quad e^{\mathbf{j}kz} = e^{\mathbf{j}kx} \cdot e^{-ky} \quad \) gives decay into the half-plane. Compute derivatives:
\[ \Phi' (z) = \mathbf{j}kA\,e^{\mathbf{j}kz} , \quad \Phi'' (z) = -k^2 A\,e^{\mathbf{j}kz} , \quad \Psi'(z) = \mathbf{j}kB\,e^{\mathbf{j}kz} \]
On y = 0, z = x, the KM stress formulas give
\begin{align*} \sigma_{xx} + \sigma_{yy} &= 4\,\Re \left[ \mathbf{j}kA\,e^{\mathbf{j}kx} \right] , \\ \sigma_{yy} - \sigma_{xx} + 2\mathbf{j}\,\sigma_{xy} &= 2 \left[ x \left( -k^2 A\,e^{\mathbf{j}kx} \right) + \mathbf{j}kB\,e^{\mathbf{j}kx} \right] . \end{align*}
From these equations, one can solve algebraically for A, B in terms of P, Q, k, obtaining a mode that matches the prescribed traction \( \displaystyle \quad T(x) = P\, e^{\mathbf{j}kx} \quad \) on y = 0.

By superposition (Fourier transform), a general nonconstant traction

\[ T(x) = \int_{-\infty}^{\infty} \,\hat{T} (k)\,e^{\mathbf{j}kx}\ {\text d}k \]
is represented by
\[ \Phi (z) = \int_{-\infty}^{\infty} \,A (k)\,e^{\mathbf{j}kz}\ {\text d}k , \qquad \Psi (z) = \int_{-\infty}^{\infty} \,B (k)\,e^{\mathbf{j}kz}\ {\text d}k , \]
with A(k), B(k) determined from the same algebraic relations mode-by-mode. ;    ■
End of Example 1
   
Example 2: We consider the half-plane y > 0, with boundary y = 0. We consider the rigid plane in half-plane, loaded by flat punch on |x| < 𝑎, normal pressure p(x), zero traction outside.

In KM form, one uses potentials with logarithmic behavior:

\[ \Phi (z) = \frac{P}{2\pi}\,\ln (z) , \qquad \Psi (z) = - \frac{P}{2\pi}\,\ln (z) \]
for a concentratied load P. For a flat punch of width 2𝑎, one builds
\[ \Phi (z) = \frac{P}{2\pi}\,\int_{-a}^a , \frac{p(t)}{t-z}\,{\text d} t, , \qquad \Psi (z) = - \Phi (z) , \]
so that the boundary traction on y = 0 matches p(x) on |x| < 𝑎, zero outside. For the classical flat punch with
\[ p(x) = \frac{P}{\pi\,\sqrt{a^2 - x^2}} , \qquad |x| < a . \]
This gives the standard Hertzian type singular pressure and corresponding stress field in closed form via φ, ψ.

In contact mechanics, a Hertzian-type singular pressure distribution refers to a classical mathematical solution where the contact pressure spikes to infinity (∞) at the boundaries or sharp edges of an indenter.

This stands in direct contrast to a standard Hertzian contact (like a sphere or a cylinder), where the pressure smoothly tapers down to zero at the edges of the contact zone. Instead, "Hertzian-type singular" profiles mimic the inverted geometric structure of the original Hertz integral equations, causing an algebraic stress singularity.

The most prominent example of a Hertzian-type singular pressure occurs when a rigid, flat-ended cylindrical indenter (of radius 𝑎) is pressed normally into an elastic half-space.    ■

End of Example 2
All structural elements have discontinuities, either created unintentionally during fabrication or developed during service conditions under repeated load to which the structure is subjected. Discontinuities will affect the strength of a structural element and may lead to fracture under a particular state of stress induced by static or fatigue loading.

Harold Westergaard and Nikoloz Mushkelishvili have developed independently general methods for stress function solutions which are well suited to solve plane problem of elasticity for cracked body. The Westergaard method is more particularly applicable to crack problems in infinite plate.

Westergaard's approach is a clever, simplified "semi-inverse" method that reduces the problem to a single complex function, but only works for specific crack geometries and symmetric loading conditions. It is a simplified version of the Kolosov–Muskhelishvili approach that is a comprehensive, mathematically rigorous method that can solve any 2D isotropic elasticity problem using two complex functions.    
Example 3: Suppose an infinite plane has a crack along     −𝑎 < x < 𝑎,     y = 0.    . In elasticity and fracture mechanics, remote tension (σ∞, spoken as "sigma infinity") refers to a uniform tensile stress applied far away from any geometric disturbance, such as a hole, crack, notch, or boundary. It is the standard, undisturbed "baseline" stress that exists in a material before it encounters a defect that causes stress concentrations.

KM/Westergaard potensials:

\[ \Phi (z) = \frac{\sigma_{\infty}}{2}\, z , \qquad \Psi (z) = \frac{\sigma_{\infty}}{2}\,\sqrt{z^2 - a^2} , \]
with the branch cut along the cruck.

From there

\[ \sigma_{yy} - \sigma_{xx} + 2\mathbf{j}\,\sigma_{xy} = 2 \left[ z\,\Phi'' (z) + \Psi' (z) \right] = \sigma_{\infty}\,\sqrt{z^2 - a^2} , \]
and the stress-intensity factor
\[ K =\sigma_{\infty}\,\sqrt{a\pi} . \]
   ■
End of Example 3
   
Example 4: We consider an infinite plate with a circular hole under uniform tension. Using KM potentials Φ(z), Ψ(z), the displacement and stress fields around a circular hole of radius 𝑎 under remote tension σ∞ are given by analytic functions:
\[ \Phi (z) = \frac{\sigma_{\infty}}{2}\,z + \frac{\sigma_{\infty}\,a^2}{2z} , \qquad \Psi (z) = - \frac{\sigma_{\infty}\,a^2}{2z} . \]
From these formulas, one obtains the famous stress concentration
\[ \sigma_{\theta\theta} (a, \theta ) = \sigma_{\infty} \left( 1 - 2\,\cos 2\theta \right) . \]
   ■
End of Example 4
   
Example 5: Let a circular inclusion of radius 𝑎 have Lamé constants (λ₁, μ₁) embedded in a matrix with (λ₀, μ₀).

KM potentials allow one to write the displacement field as

\[ \Phi (z) = Az + \frac{B}{z} , \qquad \Psi (z) = Cz + \frac{D}{z} , \]
with constants A, B, C, D determined by continuity of displacement and traction.    ■
End of Example 5
   
Example 6: Let the domain be a wedge 0 < θ < α. KM potentials allow one to write the solution as
\[ \Phi (z) = z^{\lambda} , \qquad \Psi (z) =z^{\lambda} , \]
where λ is determined by boundary conditions. This leads to the Kondrat'ev exponents for elastostatic problem in angular domain. Pri    ■
End of Example 6
   
Example 1:    ■
End of Example 1

 

 

  1. Grisvard, P., Elliptic Problems in Non-Smooth Domains,
  2. V. A. Kondrat’ev (1967) В. А. Кондратьев, «Краевые задачи для эллиптических уравнений в областях с углами» (Boundary value problems for elliptic equations in domains with corners) Trudy Moskov. Mat. Obshch., 16 (1967), 209–292 Transactions of the Moscow Mathematical Society
  3. V. Kozlov, V. Maz'ya, and J. Rossmann, Elliptic Boundary Value Problems in Domains with Point Singularities,
  4. Krainer, T., (2005) Elliptic boundary problems on manifolds with polycylindrical ends, arXiv:math/0508516 [math.AP]
  5. Krupchyk, K. & Tuomela, J., (2006) The Shapiro–Lopatinskij Condition for Elliptic Boundary Value Problems, London Mathematical Society
  6. Lopatinskii, Ya.B., (1953) On a method of reducing boundary problems for a system of differential equations of elliptic type to regular integral equations Journal: Ukrain. Mat. Zh., 5 (1953), 123–151. Лопатинский Я. Б. (1953) «Об одном способе приведения краевых задач для систем дифференциальных уравнений эллиптического типа к регулярным интегральным уравнениям» (Украинский математический журнал)
  8. S. Nazarov and B. Plamenevsky, Elliptic Equations in Polyhedral Domains,
  9. Shapiro, Z.Ya,, (1953) On general boundary problems for equations of elliptic type, Journal: Izvestiya Akad. Nauk SSSR. Ser. Mat., 17 (1953), 539–562.

 

Return to Mathematica page
Return to the main page (APMA0340)
Return to the Part 1 Matrix Algebra
Return to the Part 2 Linear Systems of Ordinary Differential Equations
Return to the Part 3 Non-linear Systems of Ordinary Differential Equations
Return to the Part 4 Numerical Methods
Return to the Part 5 Fourier Series
Return to the Part 6 Partial Differential Equations
Return to the Part 7 Special Functions