Chapter 5
Analytical techniques and solutions for
linear elastic solids
5.4 Solutions to 3D static problems in
linear elasticity
The
field equations of linear elasticity are much more difficult to solve in 3D
than in 2D. Nevertheless, several important
problems have been solved. In this
section, we outline a common representation for 3D problems, and give solutions
to selected 3D problems.
5.4.1 Papkovich-Neuber Potential
representations for 3D solutions for isotropic solids
In this section we outline a general technique for solving 3D
static linear elasticity problems. The
technique is similar to the 2D Airy function method, in that the solution is
derived by differentiating a potential, which is governed by a PDE. Many other potential representations are used
in 3D elasticity, but most are simply special cases of the general
Papkovich-Neuber representation.
Assume
that
The
solid has Young’s modulus E, mass
density and Poisson’s ratio .
The
solid is subjected to body force distribution (per unit mass)
Part of
the boundary is subjected to prescribed displacements
A second
part of the boundary is subjected to prescribed tractions
The Papkovich-Neuber
procedure can be summarized as follows:
1. Begin by finding a vector function and scalar function which satisfy
as
well as boundary conditions
2.
Calculate
displacements from the formula
3.
Calculate
stresses from the formula
HEALTH WARNING: Although the displacements and stresses that solve a
linear elasticity problem are unique, the Papkovich-Neuber potentials that
generate a particular solution are not.
Consequently, if you find several different sets of potentials in the
literature that claim to solve the same problem, don’t panic. It is likely that they really do solve the
same problem.
5.4.2 Demonstration that the
Papkovich-Neuber solution satisfies the governing equations
We
need to show two things:
1. That the displacement field satisfies the equilibrium
equation (See sect 5.1.2)
2.
That the stresses
are related to the displacements by the elastic stress-strain equations
To
show the first result, differentiate the formula relating potentials to the
displacement to see that
Substitute this result into
the governing equation to see that
Finally, substitute the
governing equations for the potentials
and
simplify the result to verify that the governing equation is indeed satisfied.
The second result can be derived by substituting the formula for displacement
into the elastic stress-strain equations and simplifying.
5.4.3 Point force in an infinite solid
The
displacements and stresses induced by a point force acting at the origin of a large (infinite)
elastic solid with Young’s modulus E
and Poisson’s ratio are generated by the Papkovich-Neuber potentials
where .
The displacements, strains and stresses follow as
5.4.4 Point force normal to the surface
of an infinite half-space
The
displacements and stresses induced by a point force acting normal to the surface of a
semi-infinite solid with Young’s modulus E
and Poisson’s ratio are generated by the Papkovich-Neuber
potentials
where
The displacements and stresses follow as
5.4.5 Point force tangent to the surface
of an infinite half-space
The
displacements and stresses induced by a point force acting tangent to the surface of a
semi-infinite solid with Young’s modulus E
and Poisson’s ratio are generated by the Papkovich-Neuber potentials
The displacements and
stresses can be calculated from these potentials as
5.4.6 The Eshelby Inclusion Problem
The Eshelby problem
(Eshelby, Proc Roy Soc A 241 p. 376)
is posed as follows:
1. Consider an infinite, isotropic, linear elastic solid,
with (homogeneous) Young’s modulus and Poisson’s ratio .
2.
The solid is
initially stress free, with displacements, strains and stresses .
3. Some unspecified external agency then induces a
uniform `transformation strain inside an ellipsoidal region, with semi-axes centered at the origin. The `transformation strain’ can be visualized
as an anisotropic thermal expansion if the ellipsoidal region were separated from
the surrounding elastic solid, it would be stress free, and would change its
shape according to the strain tensor .
4. Because the ellipsoid is encapsulated within the
surrounding elastic solid, stress, strain and displacement fields are induced
throughout the elastic solid. These
fields must be defined carefully because the initial configuration for the solid
could be chosen in a number of different ways.
In the following, will denote the displacement of a material
particle from the initial, unstressed configuration, as the transformation
strain is introduced. The total strain
is defined as
Inside the ellipsoid, the total strain consists of the
transformation strain (which does not
induce stress); together with an additional elastic strain .
Outside the ellipsoid, .
The stress in the solid is related to the elastic part of the strain by the
usual linear elastic equations
The
Eshelby solution gives full expressions for these fields. It has proved to be one of the most important
solutions in all of linear elasticity: it is of some interest in its own right,
because it provides some insight into the mechanics of phase transformations in
crystals. More importantly, a number of
very important boundary value problems can be solved by manipulating the
Eshelby solution. These include (i) the
solution for an ellipsoidal inclusion embedded within an elastically mismatched
matrix; (ii) the solution for an ellipsoidal cavity in an elastic solid; (iii)
solutions for circular and elliptical cracks in an elastic solid. In addition, the Eshelby solution is used
extensively in theories that provide estimates of elastic properties of
composite materials.
The displacement field is
generated by Papkovich-Neuber potentials
where
the integral is taken over the surface of the ellipsoid, denotes the components of a unit vector
perpendicular to the surface of the ellipsoid (pointing outwards); ,
and
is the transformation stress (i.e. the stress that would be
induced by applying an elastic strain to the inclusion that is equal to the
transformation strain). The stresses
outside the inclusion can be calculated using the standard Papkovich-Neuber
representation given in Section 5.4.1.
To calculate stresses inside the inclusion, the formula must be modified
to account for the transformation strain, which gives
For
the general ellipsoid, the expressions for displacement and stress can be
reduced to elliptic integrals, with the results
Solution inside the ellipsoid: Remarkably, it turns out that the stresses and strains
are uniform inside the ellipsoid. The
displacements, strains and stresses can be expressed as:
(i) Displacement
(ii) Strain
(iii) Stress
Here,
is a constant called the `Eshelby Tensor,’ and
is a second (anonymous) constant tensor. These
tensors can be calculated as follows. Choose the coordinate system so that .
Define
where and
are elliptic integrals of
the first and second kinds. Then
The
remaining components of can be calculated by cyclic permutations of
(1,2,3). Any components that cannot be
obtained from these formulas are zero: thus ,
,
etc. Note that has many of the symmetries of the elastic
compliance tensor (e.g ), but does not have major symmetry .
For certain special shapes
the expressions given for break down and simplified formulas must be
used
Oblate spheroid
Prolate spheroid
Sphere .
In this case the Eshelby tensor can
be calculated analytically
Additional terms follow
from the symmetry conditions . The remaining terms are zero.
Cylinder .
For this case the Eshelby tensor
reduces to
Additional terms follow
from the symmetry conditions . The remaining terms are zero.
Solution outside the ellipsoid: The solution outside the ellipsoid can also be
expressed in simplified form:
J.
D. Eshelby “The Elastic Field Outside an Ellipsoidal Inclusion,” Proceedings of
the Royal Society of London.
A252, pp. 561-569, (1959) shows that
the displacement can be obtained from a single scalar potential . For actual calculations only the derivatives
of the potential are required, which can be reduced to
where is the greatest positive root of ,
with
and . Additional derivatives can be computed using
the relations
The displacements follow as
where
The
remaining displacement components can be calculated by cyclic permutations of
(1,2,3), and strains and stresses can be calculated by differentiating the
displacements appropriately. The results
are far too complicated to write out in full, and in practice the algebra can
only be done with the aid of a symbolic manipulation program. Some special results can be reduced to a
tractable form, however:
Displacements far from the ellipsoid
Solution outside a spherical inclusion: For this case the Papkovich-Neuber potentials can be
reduced to
The displacements and
stresses follow as
where and
5.4.7 Elastically mismatched ellipsoidal
inclusion in an infinite solid subjected to remote stress

|
The
figure shows an ellipsoidal inclusion, with semi-axes . The inclusion is made from an isotropic,
elastic solid with Young’s modulus and Poisson’s ratio . It is embedded in an infinite, isotropic
elastic matrix with Young’s modulus and Poisson’s ratio . The solid is loaded at infinity by a uniform
stress state ,
strains and displacements .
The
solution is constructed by superposing the Eshelby solution to the uniform
stress state. To represent the Eshelby
solution, we introduce:
1. The Eshelby transformation strain
2. The Eshelby tensor
3. The displacement induced by the Eshelby transformation
4. The stresses induced by the Eshelby transformation
The
functions and can be calculated using the results given in
Section 5.4.6 (the elastic properties of the matrix should be used when
evaluating the formulas).
The
solution for the solid containing the inclusion follows as
where
the transformation strain is calculated by solving
for
. Here,
is the stiffness of the
matrix, with a similar expression for the stiffness of the inclusion.
5.4.8 Spherical cavity in an infinite
solid subjected to remote stress
The
figure shows a spherical cavity with radius a
in an infinite, isotropic linear elastic solid. Far from the cavity, the
solid is subjected to a tensile stress ,
with all other stress components zero.
The solution is generated
by potentials
The displacements and
stresses follow as
Derivation: This solution can be derived by superposing two
solutions:
- A uniform state of
stress ,
which can be generated from potentials ,
- The Eshelby solution
for a sphere with transformation stress .
The unknown coefficients A and B must be chosen to
satisfy the traction free boundary condition on the surface of the hole . Noting that and working through some tedious algebra shows
that
Substituting
back into the Eshelby potentials and simplifying yields the results given. The same approach can be used to derive the
solution for a rigid inclusion in an infinite solid subjected to remote stress,
as well as the solution to an elastically mismatched spherical inclusion in an
infinite solid.
5.4.9 Flat ended cylindrical indenter in
contact with an elastic half-space
A
rigid flat ended cylindrical punch with radius a is pushed into the surface of an elastic half-space with Young’s
modulus E and Poisson’s ratio by a force P. The indenter sinks into the surface by a
depth h. The interface between the
contacting surfaces is frictionless
The
load is related to the displacement of the punch by
The solution can be generated
from Papkovich-Neuber potentials
where ,
and denotes the imaginary part of z.
The displacements and stresses follow as
A symbolic manipulation
program can handle the complex arithmetic in these formulas without
difficulty. If you wish to find
analytical formulas for the displacement or stress, the following expressions
are helpful
where
Important features of these
results include:
1. Contact
pressure: The pressure exerted by the
indenter on the elastic solid follows as
2.
Surface displacement: The vertical displacement of the surface is
3. Contact
stiffness: the stiffness of a contact
is defined as the ratio of the force acting on the indenter to its displacement
,
and is of interest in practical applications.
The stiffness of a 3D contact is well defined (unlike 2D contacts
discussed in Section 5.4) and is given by . This turns out to be a universal relation for
any axisymmetric contact with contact radius a.
5.4.10 Frictionless contact between two
elastic spheres

|
This
solution is known as the `Hertz contact problem’ after its author. The figure illustrates the problem to be
solved. Two elastic spheres, with radii and elastic constants ,
initially meet at a point, and are pushed into
contact by a force P. The two spheres deform so as to make contact
over a small circular patch with radius ,
and the centers of the two spheres approach one another by a distance h.
The
solution is conveniently expressed in terms of an effective modulus and radius
for the contact pair:
Relations between : The force P,
approach of distant points h and contact area a are related by
Contact pressure: The two solids are subjected
to a repulsive pressure within the contact area. The maximum contact pressure is related to
the load applied to the spheres by
Contact stiffness: the stiffness of a contact is defined as the ratio of
the force acting on the indenter to its displacement ,
and is of interest in practical applications.
The stiffness of a 3D contact is well defined (unlike 2D contacts
discussed in Section 5.4) and is given by . This turns out to be a universal relation for
any axisymmetric contact with contact radius a.
Stress field The two spheres are subjected to the same contact
pressure, and are both assumed to deform like a half-space (with a flat
surface). Consequently, the stress field
is identical inside both spheres, and can be calculated from formulas derived
by Hamilton, Proc I. Mech E. 197C
(1983) p. 53.
where
The stresses on r=0 must
be computed using a limiting process, with the result
Conditions to initiate yield: The material under the contact yields when the maximum
von-Mises effective stress reaches the uniaxial tensile yield stress Y.
The location of the maximum von-Mises stress can be found by plotting
contours of as a function of . For the maximum value occurs at and has value . Yield occurs when .
5.4.11 Contact area, pressure, stiffness
and elastic limit for general non-conformal contacts
A
non-conformal contact has the following properties: (i) the two contacting
solids initially touch at a point or along a line; (ii) both contacting solids
are smooth in the neighborhood of the contact, so that their local geometry can
be approximated as ellipsoids, (iii) the size of the contact patch between the
two solids is much smaller than either solid.
Complete
solutions for such contacts can be found in Bryant and Keer Journal of Applied
Mechanics, 49, p. 345 (1982) or Sackfield
and Hills Journal of Strain Analysis 18
p107, p. 195 (1983). These papers also
account for the effects of friction under sliding contacts. The results are lengthy. Here, we give formulas that predict the most
important features of frictionless nonconformal contacts.
Contact Geometry: The geometry of the contacting solids is characterized as follows:
- The principal radii of
curvature of the two solids at the point of initial contact are denoted by
,
. The radii of curvature are positive if
convex and negative if concave.
- The angle between the principal directions of
curvature of the two solids is denoted by .
Note that while labels 1 and 2 can
be assigned to the radii of curvature of the two surfaces arbitrarily, must specify the angle between the two
planes containing the radii and .
- Define the principal relative contact radii as
- Introduce an effective contact radius
Elastic constants: The two contacting solids are isotropic, with
Young’s modulus and Poisson’s ratio . Define the effective modulus

|
Contact area: The area of contact between the two solids is
elliptical, with semi-axes . The dimensions of the contact area may be
calculated as follows:
- Solve the following equation (numerically) for ,
with
where
and are complete elliptic integrals of the first
and second kind
- Calculate the contact area from
(The
limit is helpful )
- The dimensions of the contact patch follow as
Contact pressure: The contact pressure distribution is ellipsoidal, with the form
Approach of contacting solids: Points distant from the contact in the two solids
approach one another by a displacement
Contact stiffness: The contact stiffness is defined as and is given by
Elastic Limit: The stresses in both solids are identical, and
therefore yield occurs first in the solid with the lower yield stress. The figure shows the critical load required
to cause yield in a solid with Von-Mises yield criterion, and uniaxial tensile
yield stress Y, based on tabular
values on page 99 of Johnson “Contact Mechanics” CUP (1985)
5.4.12 Load-displacement-contact area
relations for arbitrarily shaped axisymmetric contacts
The
most important properties of general frictionless axisymmetric contacts can be
calculated from simple formulas, even when full expressions for the stress and
displacement fields cannot be calculated.
Assume that:
- The two contacting solids have elastic constants . Define an effective elastic constant as
- The surfaces of the two solids are axisymmetric
near the point of initial contact.
- When the two solids
just touch, the gap between them can be described by a monotonically
increasing function ,
where r is the distance from the
point of initial contact. For example, a cone contacting a flat surface
would have ,
where is the cone angle; a sphere contacting a
flat surface could be approximated using where D is the sphere diameter. In the following we will use
- The two solids are
pushed into contact by a force P. The solids deform so as to make contact
over a circular region with radius a,
and move together by a distance h
as the load is applied.
- The relationship between h and the contact radius a
will be specified by a functional relationship of the form h=H(a). The derivative of this function with
respect to its argument will be denoted by
These quantities are related
by the following formulas:
1. Approach as a function of contact radius
2. Applied force as a function of contact radius
3.
Contact
stiffness
4.
Contact pressure
distribution
Once
these formulas have been evaluated for a given contact geometry the results can
be combined to determine other relationships, such as contact radius or
stiffness as a function of load or approach h.