Hankel transformation method for solving the Westergaard problem for point, line and distributed loads on elastic half-space

The Hankel transformation method was used in this work to determine the normal and shear stress distributions due to point, line and distributed loads applied to the surface of an elastic media. The elastic media considered in this study was assumed to be inextensible in the horizontal directions, and the only non-vanishing displacement component is the vertical component. Such materials were first considered by Westergaard as models of elastic half-space with alternating layers of soft and stiff materials with the stiff materials of negligible thickness and so closely spaced that the composite material characteristics is idealised as isotropic and homogeneous. The study adopted a displacement formulation. The Hankel transformation was applied to the governing Cauchy – Navier differential equation of equilibrium to reduce the problem to a second order ordinary differential equation (ODE) in terms of the deflection function in the Hankel transform space. Solution of the ODE subject to the boundedness condition yielded bounded deflection function in the Hankel transform space. Equilibrium of the internal vertical forces and the external applied load, was used to obtain the constant of integration, for the deflection function in the Hankel transform space. Inversion yielded the deflection in the physical domain variable. The stress displacement equations were then used to determine the Cauchy stresses. The vertical stress distributions due to line and distributed load over a rectangular area were also determined by using the point load solution for vertical stresses as Green functions, and then performing integration along the line and over the rectangular area of the load. The results obtained for vertical stresses due to point, line and distributed loads were determined in terms of dimensionless influence coefficients which were presented. The results obtained for the deflection, normal and shear stresses due to point, line and distributed load agreed with the solutions originally presented by Westergaard who used a stress function method.


INTRODUCTION 1.1 Background
The determination of the normal and shear stress distributions due to point, line and distributed loads applied to the surface of linear elastic media are problems of the classical mathematical theory of elasticity (Ike, 2006(Ike, , 2018a(Ike, , 2018b)).Such problems are extensively encountered and applied in the elastic stress and settlement analysis of solids and structural footings (Rocscience, 2018;Westergaard, 1938).The foundational principles for the formulation of such problems, as indeed all elasticity problems are the simultaneous consideration of a system of fifteen differential equations in fifteen unknowns, namely: three differential equations of equilibrium, six equations representing the generalized Hooke's stress -strain relations and six equations expressing strains in terms of the displacement components.The system of fifteen equations are required to additionally satisfy the compatibility conditions as well as Latin American Journal of Solids and Structures, 2019, 16(1), e155 2/19 the traction and deformation boundary conditions.The governing equations are more complex when anisotropy and heterogeneity of the three dimensional elastic media as well as the non-linear considerations of large deformations are taken into account.However, even for the cases of simplifications in the material properties of the media introduced by the assumptions of isotropy and homogeneity and the consideration of small deformations, the governing equations of the mathematical theory of elasticity are unwieldy in their numbers and demand extensively intensive analytical demanding techniques.Three methods used in the formulation of the mathematical theory of elasticity problems are the displacement method, stress method and mixed (hybrid) methods (Ojedokun and Olutoge, 2012;Westergaard, 1964;Barbar, 2010;Sadd, 2014;Bowles, 1997).The displacement methods of the theory of elasticity presented by Navier and Lamé involve a reformulation of the system of fifteen governing equations of equilibrium, stress -strain and kinematics such that only the three components of displacement (in a three dimensional problem) become the primary unknown variables (Kachanov et al., 2003;Sitharam and Govinda Reju, 2017;Palaniappan, 2011;Hazel, 2015).The advantage of the displacement formulation is the obvious reduction in the number of equations to be solved from fifteen to three coupled equations in terms of the three displacement components.The other unknowns Cauchy stresses and strains are then found using the strain displacement, and stress -strain relations after the Navier -Lamé equations are solved for the Cartesian displacement components.
Stress methods, introduced by Michell and Beltrami involve the reformulation of the system of fifteen governing equations such that only the six Cauchy stress components (in a three dimensional formulation) become the primary variables (Ike, 2006(Ike, , 2018a(Ike, , 2018b;;Rocscience, 2018;Westergaard, 1938;Ojedokun and Olutoge, 2012;Westergaard, 1964;Barbar, 2010;Sokolnikoff, 1956;Timoshenko and Goodier, 1970;Green and Zerna, 1954).The advantage is the simplification inherent in the reduction of the number of governing equations from fifteen to six.The system of coupled Beltrami -Michell equations in terms of stresses are then solved to find the six stress components.Thereafter, the strains are found using the stress -strain relations, and the displacements obtained by integration of the straindisplacement equations subject to the displacement boundary conditions.
The mixed (hybrid) method, which is not commonly used in the theory of elasticity involve a reformulation of the governing equations in terms of some Cauchy stress components and some displacement components as the primary unknown variables.
The problem of point load applied at the origin (0, 0, 0) of a homogeneous, semi-infinite, linear elastic half-space is the Boussinesq problem, found in the theory of elasticity applied to soil/solid mechanics.The Boussinesq problem is a fundamental problem whose solutions are fundamental Green functions used in deriving solutions for distributed loads acting on given lines or given areas on the boundary of the semi-infinite linear elastic soil/solid media.
The Boussinesq stress distribution theory is very extensively used in the determination of the distributions of stresses in elastic media considered semi-infinite and three dimensional.The theory assumes that the semi-infinite medium is linear elastic, homogeneous and isotropic.However, most soils are not isotropic neither are they homogeneous.
Westergaard sought to solve the theory of elasticity problem of finding stress distributions at any point under a foundation that would cater for alternating horizontal layers of soft and stiff materials (Ike, 2006(Ike, , 2018a(Ike, , 2018b;;Rocscience, 2018;Anyaegbunam, et al., 2011).For such materials with alternating soft and stiff layers, the Boussinesq stress distribution theory cannot correctly give the stresses, since the fundamental assumptions for constructing the Boussinesq theory have been violated by the material properties of the half-space.
Westergaard presented solutions for the stress fields for such problems by making the following fundamental assumptions (Ike, 2006(Ike, , 2018a(Ike, , 2018b;;Rocscience, 2018): (i) The stiff layers are very closely spaced such that the composite material characteristics could be approximated or idealized as isotropic and homogeneous.
Latin American Journal of Solids and Structures, 2019, 16(1), e155 3/19 (ii) The stiff layers are of negligible thickness as compared to the thickness of the soft materials.
(iii) The stiff layer is inextensible, in the horizontal directions.Thus the horizontal components of displacement 0 xy uu  for both the x and the y -Cartesian coordinate directions.
(iv) (iv) The inextensible stiff layer restrains any horizontal strains in both the stiff and the soft layers.The horizontal normal strains in the x and y Cartesian coordinate directions are zero, 0 ( ).
xx yy
Generally, sedimentary soils such as natural clay strata highlight the anisotropic characteristics of soil masses (Bowles, 1997).Consequently, the Westergaard theory is a better representation of the elasticity problem of three dimensional soil media (Rocscience, 2018).However many geotechnical engineers prefer the use of Boussinesq theory to model elastic soils primarily because it gives higher values of vertical stresses and could be considered more conservative from the perspective of design.

The Hankel transformation method and why it is used/adopted in this study
Hankel transformations are integral transformations whose kernels are Bessel functions, (Everitt and Kalf, 2007;Duffy, 2004;Negero, 2016).They can also be referred to as Fourier -Bessel transforms.Hankel transforms are integral transforms applied to physical domains represented by infinite two dimensional plane with polar/cylindrical coordinate (r, θ, z) where 0 , r    (Tuteja et al., 2014;Negero, 2016;Malacka Zuzana, 2018).The radial coordinate term of the Laplacian is the object of the Hankel transform.The Hankel transforms arise naturally in solving boundary value problems formulated and expressed in terms of the cylindrical polar coordinates system (Tuteja et al., 2014;Andrews and Shivamoggi, 1999;Voegtle, 2017).
They are thus particularly useful in dealing with boundary value problems that are characterised by circular symmetry (Piessens, 2000;Yokoyama, 2014).The Hankel transform of order, v equal to zero is the two-dimensional Fourier transform of a circularly symmetric function of three Cartesian coordinate variables.The most important special cases of the Hankel transform correspond to the order zero and one (i.e.v = 0, and v = 1).
The conditions for the existence of the Hankel transformation for a function and its derivatives are similar to the conditions for the existence of the Fourier transform for the function and its derivatives, and these conditions are usually satisfied for the functions that/which describe most physical and engineering models.The Hankel transform and its inverse are linear operators.
Similar to the Laplace and Fourier transforms, it is required that for the function being transformed in the Hankel transform space, the function and its derivatives with respect to the radial coordinate variable must vanish as r approaches infinity (r → ∞).(Tuteja et al., 2014;Malacka Zuzana, 2018;Negero, 2016;Andrews and Shivamoggi, 1999).
In the elastic half-space problem, the stresses and displacements are required to vanish as r → ∞, and this satisfies a key requirement for the application of the Hankel transformation.Again elastic half-space problems involving point load acting at the origin are axisymmetric and the stresses and displacement fields are expected to show circular symmetry about the axis of symmetry which is the axis of the application of the point load.Being boundary value problems with circular symmetry, the elastic half-space problem can be solved using the Hankel transform method.
Liao and Wang (1998) derived and presented closed-form solutions for the stresses and displacement fields in transversely isotropic, elastic half-space under the action of a point load.They used the Fourier and Hankel transforms in a cylindrical coordinate system in deriving the closed-form solutions.They assumed the plane of transverse isotropy to be parallel to the horizontal surface, and the half-space to be bounded by the horizontal surface (z = 0).They obtained their closed-form solution by using the superposition of two infinite spaces; one with a point load acting in the interior, and the other being free loading.Their solutions were identical with the Mindlin and Boussinesq solutions if the transversely isotropic half-space becomes homogeneous, linear elastic and isotropic.
Tarn and Wang (1987) used the Fourier and Hankel transforms to solve the three dimensional (3D) elasticity problem of transversely isotropic, infinite space, and obtained solutions which reduced to the Kelvin solution for the case of isotropic media.Closed-form solutions for the stresses and displacement fields in transversely isotropic media due to the action of a point load have also been presented by Barden (1963), and Misra and Sen (1975).
The Hankel transform method was used by Ike (2018a) in a stress-based formulation to obtain general solutions for stress and displacement fields in homogeneous, isotropic, linear elastic semi-infinite soil subject to uniformly distributed axisymmetric load over a circular area on the surface.Ike (2018a) applied the Hankel transformation to the biharmonic stress compatibility equation to obtain bounded stress functions for the elastic half-space problem.He also Latin American Journal of Solids and Structures, 2019, 16(1), e155  4/19 applied Hankel transform to the Love stress functions to obtain the stresses and displacements in the Hankel transform space; and enforced the boundary conditions to determine the unknown constants of the Love stress function in the Hankel transform space.By inversion of the Hankel transform expressions, the stresses and displacements for the elastic half-space under uniformly distributed load on a circular area were found by Ike (2018a) in the physical domain variables.His solutions were identical with solutions in literature that used Boussinesq stress functions.Ike (2018b) used the Hankel transform method to obtain general solutions for stresses and displacements in semi-infinite, linear elastic isotropic media under axisymmetric load.He used Hankel transformation of the governing equations of the 3D axisymmetric elasticity problem in a stress-based formulation to obtain the Love stress function in the Hankel transform space.Ike (2018b) also applied the Hankel transformation to the Love stresses and displacement fields to obtain general solutions for the stresses and displacements in the Hankel transform space.The general solutions obtained were used to solve the particular axisymmetric elasticity problem originally solved using Boussinesq stress potential functions by Boussinesq.Ike (2018b) obtained the same expressions for stresses and displacements fields as obtained by Boussinesq even with the Hankel transform method in a Love stress function formulation.

Research aim and objectives
The research aim is to use the Hankel transformation method to solve the Westergaard problem for point, line and distributed loads on the surface of linear elastic, isotropic half-space.The specific objectives include: (i) to use the Hankel transformation method to obtain bounded solutions to the Cauchy -Navier displacement equations of equilibrium for the Westergaard problem for the case of point load Q0 acting at the origin of the semi-infinite elastic space.
(ii) to use the condition of equilibrium of the internal vertical normal stresses and vertical point load to determine the unknown constant of integration in the integration problem of the Cauchy -Navier displacement equation of equilibrium.
(iii) to apply the Hankel transformation to the vertical stress -displacement equation and obtain the vertical normal stresses in the Hankel transform space variable, β.
(iv) to obtain the vertical stress distribution w(r, z) in the physical domain space variables by Hankel inverse transformation of the solution obtained in the Hankel transform space.
(v) to obtain the vertical normal stresses for the Westergaard point load problem in the physical domain variables using the vertical normal stress -displacement relations; and express the stresses in terms of non-dimensional Westergaard vertical stress influence coefficients.
(vi) to use the solution for a point load as Green function and obtain by superposition theory the solution for vertical stresses in the elastic half-space due to uniformly distributed line load of intensity q0 and finite length.
(vii) to use the point load solution as Green's function and obtain by superposition principles, the vertical stress distribution at any arbitrary depth under a corner of rectangular foundation areas subject to uniformly distributed load of intensity, p0.

THEORETICAL FRAMEWORK
In the Westergaard problem, the horizontal components of the displacement field are assumed to be zero, and the only non-vanishing displacement component is the vertical displacement.The displacement field components in terms of the three dimensional Cartesian coordinates x, y, z components are Cauchy -Navier displacement formulation of the differential equation of equilibrium in the absence of body force then reduce to the partial differential equation (PDE): where ∇ 2 is the Laplacian operator, given by: G is the shear modulus of elasticity given in terms of the Young's modulus of elasticity E, and the Poisson's ratio, μ as: and λ is the Lamé's constant (parameter or modulus) given by: Thus, Navier displacement equation of equilibrium can be expressed as: The strain -displacement equations from the small deformation assumptions of the three dimensional (3D) theory of elasticity become, for the Westergaard problem: where εxx, εyy and εzz are normal strains, γxy, γxz and γyz are shear strains.The stress -displacement relations for the Westergaard problem are found form the stress -strain laws using the strain -displacement relations as: where σzz is the vertical normal stress τxz and τyz are shear stresses.

Hankel transformation of the governing equations
For the Westergaard problem involving a point load Q0 the cylindrical polar coordinate system can be adopted to express the Cauchy -Navier displacement equation of equilibrium as: (22) where α(μ) is a parameter, defined to depend on the Poisson's ratio, μ.
Applying the Hankel transformation to the Cauchy-Navier differential equation of equilibrium, we have: where rJ0(βr) is the kernel (nucleus) of the Hankel transformation, β is the Hankel transformation parameter, and J0(βr) is the Bessel function of the first kind of order zero.
Using the linearity property of the Hankel transformation, we have: is the Hankel transform of w(r, z). Then, Using differential (D) operator methods or the method of trial functions, the general solution for ( , ) where c1(β) and c2(β) are the two integration constants.
For bounded solutions, ( , )  Wz  is required to be finite and defined as z → ∞.Therefore, Then for bounded solution, The constant of integration c1(β) is found from the boundary condition.

Enforcement of boundary condition
For a point load Q0 acting at the origin of the half-space as shown in Figure 1, the Hankel transformation of the load is () q  where: where in Equation ( 33 The boundary condition is established from the equilibrium of internal vertical normal stresses and external forces as: In Hankel transform space, Using the boundary condition,

Vertical displacement in the physical (real) space coordinates
The uniformly distributed line load of intensity q0 and finite length is considered to act on the surface of the halfspace as shown in Figure 2. By inversion, the vertical displacement is obtained in the domain space variables as: Latin American Journal of Solids and Structures, 2019, 16(1), e155 9/19 where

Vertical stresses in the physical domain variables
This is found from Equation ( 44) as: The other normal stresses are computed from the stress -strain relation as: 00 0 2 1 ( , , ( )) where k0 is the coefficient of lateral stress at rest.
The shear stresses are found from the shear stress -displacement relations as: (66)

Uniformly distributed line load of intensity q0 and finite length
Here, the vertical stress is found by using the solution for the point load as a Green function, as: where

, ( ))
WL I m n  is the non-dimensional vertical stress influence coefficient for uniformly distributed line load of intensity q0 and finite length for a Westergaard elastic half-space.

Vertical stress distribution below a corner point of a rectangular area due to uniformly distributed load on a rectangular foundation (B × L)
A uniformly distributed load is considered to act on a rectangular area B × L on the surface of the elastic half-space as shown in Figure 3.The vertical stress for the Westergaard problem is obtained in this case as the double integration problem over the loaded rectangular foundation area:    ( , ; .) .  of values of the Westergaard vertical stress influence coefficient for μ = 0 for various values of r/z is shown in Table 1.The results agree with results obtained by Fadum (1948).The influence values are also calculated for various values of μ and r/z and presented in Table 2.

DISCUSSION
The Hankel transformation method has been successfully used in this work to solve the Westergaard problem in the theory of elasticity involving finding vertical stresses due to a point load acting at the origin of an elastic half-space for which all the horizontal components of the displacement are assumed to vanish and the only non-vanishing displacement component is the vertical displacement.A displacement based approach was adopted in the study.Accordingly, the governing partial differential equation (PDE) of equilibrium is the Cauchy -Navier (PDE) given as Equation ( 4), which is solved by the Hankel transformation technique.Other considerations used in the solution were the strain displacement relations of the small -displacement assumptions of the 3D theory of elasticity given for the Latin American Journal of Solids and Structures, 2019, 16(1), e155 16/19 Westergaard problem as Equations (9 -14); and the stress displacement equations given for the Westergaard problem as Equations (15 -17).
Due to the symmetrical nature of the point load problem about the axis of application of the point load, the Cauchy -Navier displacement equation of equilibrium to be solved was expressed using the cylindrical polar coordinates system as Equation ( 18) or in expanded form, Equations ( 20) or ( 21).The Hankel transformation was applied to obtain Equation ( 23).The linearity property of the Hankel transformation was used to express the Hankel transformed Cauchy -Navier equation as the second order ordinary differential equation (ODE) given as Equation ( 28) in terms of ( , ), Wz  the vertical displacement in the Hankel transform space.Equation (28) was solved using methods for solving ODEs to obtain the solution as Equation ( 30) which contained two integration constants.The requirements for the boundedness of the vertical displacement was used to obtain the integration constant c2(β) as Equation (31); yielding the bounded solution for ( , )   Wz  as Equation ( 32) which still contained one constant of integration.The equilibrium of internal vertical normal stresses and external load was used to obtain the boundary condition as Equation ( 36).Enforcement of the vertical stress boundary condition gives the unknown constant of integration c1(β) as Equation ( 41); and the fully determined bounded solution for ( , )  Wz  as Equation ( 42).By inversion of Equation ( 42), the vertical displacement in the physical domain space variables was found as Equation ( 44).The vertical stress distribution was obtained from Equation (44) using the vertical stress -displacement relation as Equation ( 48).Evaluation of the integral yielded the solution as Equation ( 49).The vertical stress was expressed in dimensionless form as Equation ( 54) or ( 55) where the Westergaard vertical stress influence coefficient is defined as Equation ( 56) or (58).The normal stress distribution in the horizontal directions σxx and σyy are found as Equation ( 59), while the shear stress distributions τxy, τzy, and τzx are found as Equations ( 60), ( 63) and ( 66).
The point load solution was used to express the vertical stresses due to uniformly distributed line load of intensity q0 and finite length L as Equation ( 67).Evaluation of the integral gives the vertical stress for the distributed line load as Equation ( 70).
Similarly, the point load solution was used to express the vertical stress at points below the corner of rectangular foundations (B × L) carrying uniformly distributed load of intensity p0 as Equation ( 67).The evaluation of Equation ( 75) gives the solution in dimensionless form as Equation ( 79) where the non-dimensional coefficient was defined as Equation ( 83).
For value of the Poisson ratio, μ = 0, the Westergaard vertical stress influence coefficient for point load was found as Equation (86); which was calculated and presented in tabulated form for values of r/z ranging from r/z =0 to r/z = 10 in Table 1.Values of the Westergaard vertical stress influence coefficients for various values of r/z and μ = 0, and μ = 0.40 are also shown in Table 2.
The deflection of points (x, y, 0) on the xy coordinate plane is given by substitution of z =0 into Equation (44) by: However, at x = y = z =0 which is the origin and the point of application of the load, the deflection function becomes undefined and indeterminate, and this is expectedly due to the singular property of the point load acting at the origin, which introduces a singularity in the deflection function at the point of application (0,0,0).The normal and shear stress fields are similarly undefined and singular at the origin.However, on the xy coordinate plane at points (x, y, z = 0) the normal and shear stress fields are found by substitution of z = 0 into the expressions for the stresses.
It is thus observed that at points on the xy Cartesian coordinate plane, the shear stresses τzy and τzx do not vanish but have definite values that depend on the coordinates of the point; and this is in line with observations made by Anyaegbunam et al. (2011).However, the non-vanishing of the shear stresses τzy and τzx on the xy Cartesian coordinate plane do not in any way vitiate the results obtained, since the shear stress free boundaries were never imposed as a condition for the formulation of the governing equations unlike the classical Boussinesq problem where the displacement functions were derived to satisfy apriori the shear stress free boundary conditions on the xy plane.

CONCLUSIONS
The following conclusions are made from the study: (i) The Hankel transformation method has been successfully used to solve the Westergaard problem of finding vertical stresses in the elastic half-space due to a point load at the origin, and the solution used as Green functions to obtain solutions for uniformly distributed loads on line and on a rectangular area.
(ii) The Hankel transformation of the governing Cauchy -Navier displacement equation of equilibrium formulated in cylindrical polar coordinates simplified the boundary value problem (BVP) to a second order ordinary differential equation which could be solved using methods for solving ODEs, to obtain bounded solutions for the displacement.
(iii) The displacement function obtained was found to be defined and single valued at all points in the elastic halfspace except at the origin which is the point of application of the concentrated load at which the displacement function becomes singular, (undefined or indeterminate).
(iv) The normal stresses at any point on the xy coordinate plane are singular, (undefined or indeterminate).
(v) The shear stress τxy on the xy plane, as well as on the entire elastic half-space is zero.This is due to the foundational assumption of longitudinal inextensibility which made u = v = 0, and thus the shear stress γxy would naturally vanish from the strain -displacement relations for all points in the elastic half-space.
(vi) The shear stresses τzy and τzx do not vanish at points on the xy plane for which z = 0, but are singular (undefined or indeterminate) at the origin.
y z  are respectively the x, y, and z Cartesian coordinate components of the displacement.Latin American Journal of Solids and Structures, 2019, 16(1), e155 5/19

Figure 1 :
Figure 1: Point load at the origin O of an elastic half space made of Westergaard material

Figure 2 :
Figure 2: Uniformly distributed line load of finite length acting on the surface of an elastic half-space  is the non-dimensional vertical stress influence coefficient for uniformly distributed load on a rectangular foundation area for points at any depth z under a corner point in a Westergaard elastic half-space.

Figure 3 ;For
Figure 3; Uniformly distributed load of intensity p0 on a rectangular foundation (B × L) for various value of r/z can be calculated and presented in Tables.Table  is the Westergaard vertical stress influence coefficient at the point (r, z) in soil with Poisson's ratio, μ.
w I r z

Table 1 :
Influence values for vertical, normal stresses due to point load at the origin (Hankel transform solution of theWestergaard problem)

Table 2 :
Westergaard vertical stress influence coefficient for point load applied at the origin on the surface of an elastic half-space