the kelvin problem

arXiv:1202.1719v1 [math-ph] 8 Feb 2012On the Kelvin Problem Antonino Favata Dipartimento di Ingegneria Civile, Universit` a di Roma Tor Vergata1 Abstract The Kelvin problem of an isotropic elastic space subject to a concentrated load is solved in a manner that exploits the problem’s built- in symmetries so as to determine in the first place the unique balanced and comp atible stress field. Keywords: Kelvin Problem, Concentrated Loads, Unlimited Domains, Li n- ear Elasticity 1 Introduction This paper deals with a classical problem in linear elastici ty, whose solution in terms of displacements was given by Lord Kelvin in a short p aper dated 1848 [12]. The problem (see Figure 1) consists in finding the e quilibrium Figure 1: The Kelvin problem. state of a linearly elastic, isotropic material body occupy ing the whole space 1Via Politecnico 1, 00133 Rome, Italy. Email:[email protected] 1 and being subject to a point load. Kelvin’s is the fundamental solution (in modern terminology, the Green function) later used by other XIX century elasticians to determine a number of strain nuclei (such as e.g. centers of dilatation and rotation), that is to say, of solutions to p roblems involv- ing unbounded elastic domains subject to force doublets (ordipoles) with or withoutmoment, invarious combinations, all inducingsing ular stress states. Singular stress states are also encountered when modeling t he strongly lo- calized residual stress fields due to, say, gas trapped in sma ll cavitites or the thermal mismatch associated with an inclusion. More recent ly, Kelvin solu- tion has been widely used as a basic ingredient in integral-e quation methods for solving boundary-value problems in elastostatics (e.g ., in boundary ele- ment methods [13]). Accounts of Kelvin problem, of different levels of clarity, co mpleteness and purposes, are found in many expositions of the mathemati cal theory of classic elasticity: we recall the books by Love [6], Sokolni koff [10], Malvern [7], Benvenuto [2], Davis & Selvadurai [3], Barber [1], and S add [9], but there are others. All these authors tackle the problem by sol ving the Navier equations for a displacement field having the representatio n constitutively implied by a tentative representation of the stress field in t erms of a given scalar potential . We hereproposeadifferent approach, wherethebasicequa- tions of an elasticity problem — namely, the equilibrium, co nstitutive, and compatibility equations — are not combined preliminarily i nto one partial differential equation for displacement, as Navier did, but in stead used se- quentially, in a semi-heuristic fashion, guided by the disp lacement and stress symmetries intrinsic to the problem at hand. In the two subsections to follow, we briefly sum up first the cla ssical resolution technique, then ours; certain not completely ob vious technical steps of the former, that are not found in any of the elasticit y books quoted above, are contained in Appendix 1. Displacement and stress symmetries are discussed in Section 2; the balanced and compatible stre ss field that solves Kelvin problem is constructed in Section 3, the bulk o f this paper; the easy task of deriving the associated strain and displace ment fields is carried out in Section 4. 1.1 Kelvin solution according to Love In his paper [12], Kelvin does not disclose the details of the technique he adoptedtosolvetheproblemnamedafterhim: inlittlemoret hantwopages, he just proposes the solution, stressing the similarity wit h the problem in the theory of thermal conduction, where a point heat source parallels the 2 role of a strain source under form of a concentrated load. We now give a quick account of Kelvin problem and its solution modeled aft er the account given by Love [6]. Love’s procedure is essentially the same a dopted by those authors who solved problems of the same type of Kelvin’s afte r him; while the flow is clear, the mathematical developments are at times skipped. We warn the reader that the notation we use is rather different fro m the original — but coherent with the rest of the paper — and that some slight changes in presentation have been found either necessary or simply c onvenient. Navier equation for the displacement vector field u=u(x) is: (λ+µ)∇(Divu)+µ∆u+d=0, (1) whereλandµare the Lam´ e constants and dis the body force. With a view toward solving this equation, both the displacement and the force field are given a Helmhholtz representation in terms of potential pai rs: u=∇ϕ+curlw,d=∇ψ+curlb,with Div w= Divb=0.2(2) The distance force is taken null in H \ Bρ, where His the whole space andBρ⊂ Hdenotes a sphere of radius ρabout the point owhere the concentrated load fis applied. On recalling that ∆ ∇(·) =∇∆(·) and ∆curl( ·) = curl∆( ·), equation (1.1) can be written as: ∇/parenleftbig (λ+2µ)∆ϕ+ψ/parenrightbig +curl/parenleftbig µ∆w+∇ψ+b/parenrightbig =0. (3) A solution of this equation can be obtained by solving the two equations: (λ+2µ)∆ϕ+ψ= 0, µ∆w+b=0,3(4) 2We recall that, given any sufficiently smooth field uover a bounded regular region R, itsHelmhholtz representation consists in the pair of a scalar field ϕand a divergenceless vector field woverRsuch that u=∇ϕ+curlwwith Div w=0, (foru∈C(¯R)∩CM(R),M≥1, bothϕandware of class CM(R)). A straightforward application of the identities curl ∇ϕ=0and div∇ϕ= ∆ϕ, yields: curlu= curlcurl wand div u= ∆ϕ. 3Note that (1.1) is recovered by taking the gradient of the firs t equation and the curl of the second, and then summing up the results. 3 both set over H\Bρ; the solutions ϕρandwρof (1.1) are: ψρ(x) =−1 4π/integraldisplay Bρd(ξ)·∇ξ(γ−1(x,ξ))dv(ξ), bρ(x) =−1 4π/integraldisplay Bρd(ξ)×∇ξ(γ−1(x,ξ))dv(ξ). where, forγ(x,ξ):=|x−ξ|, (4πγ)−1is the Green kernel of the laplacian. The potential pair yielding the representation (1.1)1of the solution to (1.1) is found by taking the limits of ϕρandwρforρ→0: ϕ=1 8πf·∇ρ,w=1 8πf×∇ρ, ρ(x):=|x−o|. (5) Finally, the displacement vector is determined by insertin g (1.1) in (1.1)1. Fore2,e3two unit vectors completing with the direction e1of the applied load an orthonormal triplet, one finds: u1=−(λ+µ)f 8πµ(λ+2µ)∂2ρ ∂x2 1+f 4πµρ=f 16πG(1−ν)/parenleftbigg2(1−2ν) ρ+1 ρ+x2 1 ρ3/parenrightbigg , u2=(λ+µ)f 8πµ(λ+2µ)∂2ρ ∂x1x2=fx1x2 16πG(1−ν), u3=(λ+µ)f 8πµ(λ+2µ)∂2ρ ∂x1x3=fx1x3 16πG(1−ν), (6) whereui:=u·ei. 1.2 Our method As shown in the previous subsection, it was Kelvin’s concern to find a repre- sentation of the displacement field that would make Navier eq uations (1.1) manageable; with his procedure, the elastic state {u,E,S}is determined sequentially, starting from the knowledge of u. We here propose a different, in fact inverse path. Firstly, we find the unique stress field SinH\{o}that is: (i)balanced, in that it satisfies the balance equation DivS=0, (ii)compatible , in that it satisfies the compatibility equation ∆S+1 1+ν∇∇(trS) =0, (7) 4 so that, in particular, ∆(trS) = 0; (8) (iii)balances the applied load , in the sense that /integraldisplay ∂NSnda+f=0, (9) for every full-dimensional neighborhood N ⊂ Hof the point owhere theforce fisapplied(here ndenotes theouter normaltotheboundary ∂NofN). Once the stress field is found, the deformation field Ecan be computed by using the constitutive equation: E=1 2G/parenleftBig S−ν 1+ν(trS)I/parenrightBig , (10) whereGandνare the shear and Poisson’s moduli, and Iis the identity tensor. Finally, the determination of the displacement fiel d is accomplished by solving the differential equations: ∇u+∇uT= 2E. 2 Symmetries For a problem like Kelvin’s, our intuition captures and rest itutes a num- ber of symmetries that the displacement and contact-intera ction fields must possess. As we shall see, recognizing such symmetries has th ree direct and useful consequences: (i) it reduces the representation of the stress tensor to a mo re tractable form; (ii) it makes it simpler to satisfy the compatibility equati on; (iii) it furnishes a condition that is by all means equivalen t to a bound- ary condition even in the present case of an infinite domain wi th no boundary. Different kinds of symmetries can manifest, depending on the b ody part we consider, either a right circular cylinder whose axis is dir ected as the load and passes through its application point (Figure 2) or a ball centered at the load. 5 Figure 2: Cartesian axes and cylindrical coordinates. 2.1 Symmetries of the displacement field We expect the displacement field that solves Kelvin problem t o have cylin- drical symmetry, that is, when represented in the basis {e,h,h′}, to besuch that bothuzandurare independent of ϕand, in addition, uϕis null: u(z,r,ϕ) =uz(z,r)e1+ur(z,r)h(ϕ). (11) Transforming this expectation into a parametric represent ation of the can- didate solution has two relevant consequences. 1.Consequences on the strain and stress fields. In view of the definition of the strain measure E, we have that 2Eϕr=uϕ,r+r−1(ur,ϕ−uϕ), 2Eϕz=uϕ,z+r−1uz,ϕ. Accordingly, as it is not difficult to see, whenever (2.1) prev ail, both strain components EϕrandEϕzturn out to be everywhere null. But then, given that His supposedto becomprised of an isotropic, linearly elastic material, both stress components SϕrandSϕzmust also be everywhere null: Sϕr≡0 andSϕz≡0. Therefore, the class of stress fields of interest has the redu ced repre- sentation: S=σ1e1⊗e1+σ2h⊗h+σ3h′⊗h′+σ4(e1⊗h+h⊗e1),(12) 6 parameterized by the four scalar-valued mappings /hatwideσi, with /hatwideσ1(z,r) =Szz=S·e1⊗e1, /hatwideσ2(z,r) =Srr=S·h⊗h, /hatwideσ3(z,r) =Sϕϕ=S·h′⊗h′, /hatwideσ4(z,r) =Szr=S·e1⊗h. 2.Consequences on the compatibility condition. It follows from (2.1) that ∇u=uz,ze1⊗e1+ur,zh⊗e1+r−1urh′⊗h′+uz,re1⊗h+ur,rh⊗h; consequently, Ezz=uz,z, Err=ur,r, Eϕϕ=r−1ur, Ezr=1 2/parenleftBig uz,r+ur,z/parenrightBig , (13) and hence, in particular, Err= (rEϕϕ),r. (14) In view of the inverse constitutive equation (1.2), this con sistency con- dition of representation (2.1) can be rewritten as: σ2−(rσ3),r+νrα,r= 0. (15) In addition to (2.1), another symmetry condition prevails, namely, lim r→0+ur(z,r) = 0; (16) asweshallseelater on, thislast conditionworksinallresp ectsasaboundary conditions. Remark. Condition (2) is related to the compatibility equation curlcurlE=0. In fact, given the noted symmetries in the displacement and d eformation fields, it is not difficult to see that the e1⊗e1−component of curlcurl Eis given by: curlcurlE·e1⊗e1=:(curlcurl E)zz= 0 =r−1Err,r−Eϕϕ,rr−2r−1Eϕϕ,r, or rather, equivalently, Err,r= (rEϕϕ),rr, a direct consequence of (2). 7 2.2 Symmetries on the traction Recall, to begin with, thefollowing relationships between thecylindrical and spherical coordinates of a given point: z=ρcosϑ, r=ρ|sinϑ|;ρ2=z2+r2,|tanϑ|=r z (Figure 3). Figure 3: Cartesian axes and spherical coordinates. When written as functions of the spherical coordinates, the stress com- ponents /hatwideσiwill be denoted by /tildewideσi, with /hatwideσi(z,r) =/hatwideσi(ρcosϑ,ρ|sinϑ|) =:/tildewideσi(ρ,ϑ), a set of functions assumed to be such that /tildewideσi(ρ,ϑ) =/tildewideσi(ρ,−ϑ)∀ρ>0. (17) Consider now a body part under form of a ball Bρcentered at the point of application of the load. On writing /hatwiden(ρ,ϑ,ϕ) = cosϑe1+|sinϑ|h(ϕ),(ϑ,ϕ)∈(0,+π)×(0,2π) for the outer normal to the boundary ∂BρofBρ, the traction on that surface has the expression: Sn= (cosϑσ1+|sinϑ|σ4)e1+(|sinϑ|σ2+cosϑσ4)h. According to (iii), the equilibrium condition is: fe1=−/integraldisplay ∂BρSnda=−ρ2/integraldisplayπ 0/parenleftBig/integraldisplay2π 0Sndϕ/parenrightBig |sinϑ|dϑ =−2πρ2/parenleftBig/integraldisplayπ 0(cosϑ/tildewideσ1(ρ,ϑ)+|sinϑ|/tildewideσ4(ρ,ϑ))|sinϑ|dϑ/parenrightBig e1,(18) 8 whateverρ >0. Therefore, for the right side to remain finite when ρis chosen arbitrarily big, it is necessary that /tildewideσ1(ρ,ϑ) =ρ−2/tildewideτ1(ϑ)+o(ρ−2); /tildewideσ4(ρ,ϑ) =ρ−2/tildewideτ4(ϑ)+o(ρ−2). This result suggests the Ansatz: /tildewideσi(ρ,ϑ) =ρ−2/tildewideτi(ϑ),/tildewideτi(ϑ) =/tildewideτi(−ϑ) (i= 1,...,4),(19) which takes into account the parity conditions (2.2). 3 The stress field Given that, when using cylindrical coordinates, DivS= (Se1),z+(Sh),r+r−1(Sh′),ϕ+r−1Sh, a stress field in the class (1) is balanced for null distance fo rces if 0=r(σ1e1+σ4h),z+/parenleftbig r(σ2h+σ4e1)/parenrightbig ,r+(σ3h′),ϕ =/parenleftbig rσ1,z+(rσ4),r/parenrightbig e1+/parenleftbig rσ4,z+(rσ2),r−σ3/parenrightbig h, inH\{0}. The choice of the parameter mappings is therefore restrict ed to those satisfying the following partial differential equatio ns: rσ1,z+(rσ4),r= 0, (20) rσ4,z+(rσ2),r−σ3= 0. (21) Moreover, it is not difficult to see that, under the present cir cumstances, the vectorial compatibility condition (ii) is equivalent to a s ystem of four scalar equations, namely, ∆σ1+α,zz= 0, ∆σ2−2r−2(σ2−σ3)+α,rr= 0, ∆σ3+2r−2(σ2−σ3)+r−1α,r= 0, ∆σ4−r−2σ4+α,zr= 0,(22) where α:= (1+ν)−1trS= (1+ν)−1(σ1+σ2+σ3) (23) must be harmonic: ∆α= 0, (24) to satisfy (ii). To solve (3), we propose the following seque ntial procedures: 9 1. to determine, a multiplicative constant apart, an approp riate solution α=/hatwideα(z,r) of the Laplace equation (3); 2. to integrate (3)1forσ1, in the form: σ1=−∆−1[α,zz]+c1α, (25) where ∆−1denotes the integral operator that formally inverts the laplacian, and c1is a constant to be determined. 3. to integrate (3)4forσ4; 4. to determine the fields /hatwideσ2and/hatwideσ3, by solving the system of (3)2and (2): σ2+σ3= ∆−1[α,zz]+c2α, σ2−(rσ3),r+νrα,r= 0,(26) wherec2is a second constant, such that c1+c2= 1+ν.4(27) Remark. For an alternative procedure, note the following conseque nce of (3) and (3): ∇(rσ4) =−((rσ2),r−σ3)e1−rσ1,zh. (28) Were step 4. taken right after steps 1. and 2., the field /hatwideσ4could be de- termined by integrating equation (3) along any regular curv eC, arc-length parameterized, with tangent t, and going from a fixed point x0to the vari- able pointx: (rσ4)|x x0=−/integraldisplay C((rσ2),r−σ3)e1+rσ1,zh)·tds. (29) Remarkably, given the extreme points, the choice of the join ing curve is irrelevant, due to a well known result in the theory of differen tial forms: for Ropen and star-shaped, and for ω:=ω1e1+ω2e2+ω3e3a vector field of classC1(R), the differential form ω=ω1dx1+ω2dx2+ω3dx3is exact if and only if curl ω=0inR. In our case, the differential form and associated vector field under scrutiny are: ω=/parenleftbig (rσ2),r−σ3/parenrightbig dz+rσ1,zdr,ω=/parenleftbig (rσ2),r−σ3/parenrightbig e1+rσ1,zh; 4This condition is arrived at by adding (3) and (3)1and by taking into account of (3). 10 consequently, curlω=/parenleftbig σ1,zz−((rσ2),r−σ3),r/parenrightbig h′. Now, by differentiating (3) with respect to zand (3) with respect to r, and by subtracting the resulting relations, we obtain: σ1,zz−((rσ2),r−σ3),r= 0, which allows us to conclude that curl ω=0, and hence that ωis exact. 3.1 Determination of α Relations (2.2) imply a preliminary representation for tr S: trS=ρ−2/tildewideτ(ϑ),/tildewideτ(ϑ) =/tildewideτ(−ϑ). For such a field to be harmonic: ∆(ρ−2/tildewideτ(ϑ)) = 0, function /tildewideτmust satisfy the following ordinary differential equation: sinϑτ′′+cosϑτ′+2sinϑτ= 0, whose only regularsolution is /tildewideτ(ϑ) =τ0cosϑ. All in all, we are induced to choose trS=τ0ρ−2cosϑ, (30) withτ0aconstant proportional tobeappliedload, to bedetermined later on (a method to arrive at (3.1) without using (2.2) is expounded in Appendix 2). At this point, on recalling (3), we have: α=α0ρ−2cosϑ,withα0= (1+ν)−1τ0. 3.2 Determination of σ1 We start by showing how to solve the equation (3)1, that we here reproduce: ∆σ1+α,zz= 0, 11 with the use of the information we gathered so far: that the un known field has the preliminary representation (2.2)1: σ1=/tildewideσ1(ρ,ϑ) =ρ−2/tildewideτ1(ϑ),/tildewideτ1(ϑ) =/tildewideτ1(−ϑ), and that α=/tildewideα(ρ,ϑ) =α0ρ−2cosϑ. The equation we have to solve for /tildewideτ1is: τ′′ 1+cotϑτ′ 1+2τ1+3α0cosϑ(2cos2ϑ−3sin2ϑ) = 0.(31) A particular solution of (3.2) can be obtained by the method o f parameter variation: /tildewideτ(p) 1(ϑ) =3 2α0cos3ϑ; moreover, the even solutions of the homogeneous equation as sociated with (3.2) are: /tildewideτ(h) 1(ϑ) =β0cosϑ, β0= a constant , in conclusion, /tildewideσ1(ρ,ϑ) =ρ−2/parenleftbig3 2α0cos3ϑ+β0cosϑ/parenrightbig . (32) 3.3 Determination of σ4 We have to solve (3)4for a field /tildewideσ4that as the representation (2.2)4, namely, σ4=/tildewideσ4(ρ,ϑ) =ρ−2/tildewideτ4(ϑ),/tildewideτ4(ϑ) =/tildewideτ1(−ϑ), The function /tildewideτ4is determined by the ordinary differential equation τ′′ 4+cotϑτ′ 4+(1−cot2ϑ)τ4+3α0|sinϑ|(4cos2ϑ−sin2ϑ) = 0, whose particular solutions are: /tildewideτ(p) 4(ϑ) =3 2α0cos2ϑ|sinϑ|, (33) while the even solutions of the homogeneous associated equa tion are: /tildewideτ(h) 4(ϑ) =γ0|sinϑ|, γ0= a constant . We conclude that: /tildewideσ4(ρ,ϑ) =ρ−2/parenleftbig3 2α0cos2ϑ|sinϑ|+γ0|sinϑ|/parenrightbig . (34) 12 As to the constants β0andγ0, we note that the fields (3.2) and (3.3) satisfy the balance equation (3) only if β0=γ0; with this choice, the balance equation (3) is also satisfied. Furthermore, since the part-wise balance condition (2.2) can now be writt en as f=−2π/integraldisplayπ 0(cosϑ/tildewideτ1+|sinϑ|/tildewideτ4)|sinϑ|dϑ, we find that α0+2β0=−f 2π. (35) Remark. Were we dealing with another classic problem in linear elas ticity, Boussinesq’s, the problem of a half space with a concentrate d load perpen- dicular to the boundary, we could employ so far the same solvi ng procedure. At this point, though, we would have had to satisfy the bounda ry condition, and this requires that β0= 0; the reader is referred to [8] for details. Remark. On using cylindrical coordinates, we have that /hatwideσ1(z,r) =3 2α0z3 (z2+r2)5/2+β0z (z2+r2)3/2, /hatwideσ4(z,r) =3 2α0z2r (z2+r2)5/2+β0r (z2+r2)3/2, /hatwideα(z,r) =α0z (z2+r2)3/2. A consequence of the first and third of these relations is: /hatwideσ1(z,r) =3 2α0z3 (z2+r2)5/2+β0 α0/hatwideα(z,r), Since ∆−1[α,zz] =−3 2α0z3 (z2+r2)5/2, (36) a comparison with (3) shows that c1=β0/α0. (37) 13 3.4 Determination of σ2andσ3 The functions σ2andσ3we seek are to solve system (3), which, with the use of (3), (3.3), and (3.3), can be rewritten as: σ2+σ3=−3 2α0z3 ρ5+(α0(1+ν)−β0)z ρ3, σ2−(rσ3),r−3α0νrz ρ5= 0. This system can be solved sequentially for σ3andσ2; it can be checked that solutions have the form: /hatwideσ3(z,r) =−/parenleftbig α0(1−2ν)−β0/parenrightbigz ρ3−/parenleftbig α0(1−2ν)−2β0/parenrightbigz3 2r2ρ3+g(z) r2, /hatwideσ2(z,r) =−3 2α0z3 ρ5+(α0(1+ν)−β0)z ρ3+ +/parenleftbig α0(1−2ν)−β0/parenrightbigz ρ3+/parenleftbig α0(1−2ν)−2β0/parenrightbigz3 2r2ρ3−g(z) r2, parameterized by an arbitrary function g(z).5 3.5 Wrapping up At this point, we have satisfied the equilibrium and compatib ility equations, and we have balanced the load. However, we have not yet comple tely found the solution, because the form of function g(z) and the value of the constant α0arestillwanted: weshallfindthisinformationbyexploitin gthesymmetry property (2.1) of the displacement field. By using the inverse constitutive law (1 .2)2, we find that Eϕϕ=−1 2G/parenleftbigg/parenleftbig α0(1−ν)−β0/parenrightbigz ρ3+/parenleftbig α0(1−2ν)−2β0/parenrightbigz3 2r2ρ3−g(z) r2/parenrightbigg , whence ur=rEϕϕ=−1 2G/parenleftbigg/parenleftbig α0(1−ν)−β0/parenrightbigzr ρ3+/parenleftbig α0(1−2ν)−2β0/parenrightbigz3 2rρ3−g(z) r/parenrightbigg . (38) With problemsformulated onadomainwithnoboundary,aninc omplete, no matter if inessential, specification of thesolution is to be expected, dueto an 5Recall that the constant β0can be determined by means of (3.3) in terms of α0and the magnitude of the applied load. 14 inevitable deficiency of conditions. This is not the case for the Kelvin prob- lem, where the deficit is covered by the symmetry condition (2 .1), namely, lim r→0+ur(z,r) = 0, which plays the role of a Dirichlet boundary condition. Inse rting (3.5) in it, we deduce that: α0(1−2ν)−2β0= 0, g(z) = 0. The first of these relations, together with (3.3), yields the values of the constantsα0andβ0: α0=−f 4π(1−ν), β0=−f(1−2ν) 8π(1−ν). (39) We are now in a position to write the expression of the stress fi eld: /hatwideσ1(z,r) =−f 8π(1−ν)/parenleftbigg 3z3 ρ5−(1−2ν)z ρ3/parenrightbigg , /hatwideσ2(z,r) =−f 8π(1−ν)/parenleftbigg 3zr2 ρ5+(1−2ν)z ρ3/parenrightbigg , /hatwideσ3(z,r) =f(1−2ν) 8π(1−ν)z ρ3, /hatwideσ4(z,r) =−f 8π(1−ν)/parenleftbigg 3z2r ρ5−(1−2ν)z ρ3/parenrightbigg .(40) Remark. In case the procedure delineated in Remark 2 is applied, the dis- covery of condition (3.3) is delayed. It is convenient to cho ose the curve C as the union of two different curves C=C1∪C2: C1:={x(z,r,ϕ)∈ H|z≥z0,r=r0,ϕ=ϕ0},C2:={x(z,r,ϕ)∈ H|r≥r0,ϕ=ϕ0}, whose tangent vectors are, respectively, t1=e1,t1=h(ϕ0). The integral in (3) can be determined as: /integraldisplay Cω·t=/integraldisplay C(ωze1+ωrh)·t=/integraldisplayz z0ωz(t,r0)dt+/integraldisplayr r0ωr(z,t)dt, with ωz:= (rσ2),r−σ3ωr:=rσ1,z. 15 It is not difficult to see that /integraldisplay Cω·t=f 8π(1−ν)/parenleftbigg 3z2r2 ρ5−(1−2ν)zr ρ3−(1−2ν)r4 0−2(2−ν)z2 0r2 0 ρ5 0/parenrightbigg , whereρ0:=/hatwideρ(z0,r0). On invoking equation (3), we have that r/hatwideσ4(z,r)−r0/hatwideσ4(z0,r0) =−f 8π(1−ν)/parenleftbigg 3z2r2 ρ5−(1−2ν)zr ρ3−(1−2ν)r4 0−2(2−ν)z2 0r2 0 ρ5 0/parenrightbigg , whence, by choosing ( z0,r0) = (0,0), we recover (3.5)4. 4 The strain and displacement fields TodeducetheKelvinstrainfield, wehaverecoursetotheinve rseconstitutive equation (1.2). We find: Ezz=−f 16πG(1−ν)ρ5/parenleftbig 4(1+ν)z3+(1−4ν)zr2/parenrightbig , Err=f 16πG(1−ν)ρ5/parenleftbig z3−2zr2/parenrightbig , Eϕϕ=f 16πG(1−ν)z ρ3, Ezr=−f 16πG(1−ν)ρ5/parenleftbig 2(2−ν)z2r+(1−2ν)r3/parenrightbig .(41) As to the displacement field, given (3.5) and (3.5), we have th at ur=f 16πG(1−ν)zr ρ3; (42) moreover, on integrating the obvious combination of (2)1and (4), we obtain the following preliminary representation of uz: uz=f 16πG(1−ν)/parenleftbigg2(1−2ν) ρ+1 ρ+z2 ρ3/parenrightbigg +h(r). With this, (2)4, and (4), we find that h′(r) = 0 ⇔h(r) =h0e1, and we take h0= 0 (cf. (1.1)1). 16 Remark. Forh0/\e}atio\slash= 0, the vector h0represents an arbitrary translation of the whole space in the vertical direction. Such an indetermi nacy was to be expected, because the solution of an elasticity problem on a region without boundary is always determined to within a rigid displacemen t. In this spe- cific case, a vertical translation of an arbitrary amount is t he only type of rigid displacement compatible with the postulated symmetr ies. We regard disposingofsuchanarbitrarinessbysetting h0= 0asacovenient completion of the boundary conditions . Appendix 1. Developments complementing those in Subsection 1.1 From the classic point of view, a concentrated force is the li mit of a body force field having a shrinking support; according to the prec ise definition found in [11], a sequence {dn}of body force fields defined on an open neigh- borhood Rof a pointotends to the load fconcentrated at oif: (i)dn∈C2(R); (ii)dn=0suR\Brn(o), where {Brn(o)}is a sequence of spheres of radius rnsuch thatrn→0 whenn→ ∞; (iii) lim n→∞/integraldisplay Rdn=f; (iv) the sequence {/integraltext R|dn|}is bounded. Given the vector d, the determination of the fields ψandbsatisfying the Helhmoltz’s decomposition (1.1)2can be achieved in two steps. 1. On applying the divergence operator to the equation (1.1) , we obtain the Poisson equation ∆ψ= Divd, whose solution has the well-known representation (see e.g. [4]): ψ(x) =−/integraldisplay BrG(x,ξ)Divd(ξ)dv(ξ), x∈ H, where G(x,ξ) =/parenleftbig 4πρ(x,ξ)/parenrightbig−1, ρ(x,ξ):=|x−ξ|. 17 On recalling the identity Div(ϕv) =ϕDivv+v·∇ϕ, forϕandvsmooth scalar-valued and vector-valued fields, and on using the divergence theorem, we obtain: ψr(x) =−1 4π/integraldisplay Brd(ξ)·∇ξ(ρ−1(x,ξ))dv(ξ). 2. On taking the curl of d, we obtain: curld= curlcurl b=∇(Divb)−∆b=−∆b. Thus, we have to solve another Poisson equation: −∆b= curld, whose solution is: b(x) =1 4π/integraldisplay BrG(x,ξ)curld(ξ)dv(ξ). An application of Stokes theorem combined with the identity curl(ϕv) =ϕcurlv+∇ϕ×v yields: br(x) =−1 4π/integraldisplay Brd(ξ)×∇ξ(ρ−1(x,ξ))dv(ξ). We now computethelimits of ψrandbrforr→0, underthe assumption that lim r→0/integraldisplay Brd(x)dv(x) =f. We find: ψ(x) =−1 4πf·∇(r−1),b(x) =−1 4πf×∇(r−1), r(x):=|x−o|, or rather, since 2 ∇(r−1) = ∆∇r, ψ(x) =−1 8π∆/parenleftbig f·∇r/parenrightbig ,b(x) =−1 8π∆/parenleftbig f×∇r/parenrightbig . We can now write system (1.1) as follows: ∆/parenleftbigg ϕ−1 8πf·∇r/parenrightbigg = 0, ∆/parenleftbigg w−1 8πf×∇r/parenrightbigg =0, and read out the particular solution (1.1). 18 Appendix 2. An alternative way to obtain (3.1) A way to find solutions independent of ϕto the Laplace equation in cylin- drical coordinates: ∆α=α,zz+α,rr+r−1α,r= 0 consists in looking for solutions, if any, having the form /hatwideα(z,r) =ρazbrc. It is not dfficult to see that, for such an Ansatz to be successfu l, the expo- nentsa,b,cmust satisfy the following condition: a(a+2b+2c+1)+b(b−1)ρ2z−2+c2ρ2r−2= 0∀z,r>0, or rather, equivalently, they must be chosen so as to satisfy the following three algebraic conditions: c= 0, b(b−1) = 0, a(a+2b+2c+1) = 0. There are four possibilities: (i) b= 0 anda= 0; (ii)b= 0 anda= 1; (iii) b= 1 anda= 0; (iv)b= 1 anda=−3; accordingly, the desired field must have the following form: /hatwidea(z,r) =α0z ρ3+α11 ρ+α2z+α3. (43) This result can be applied when, as it happens when dealing wi th Kelvin problem, the trace of a stress field that solves Kelvin proble m is found to be proportional to a harmonic field that, when expressed in cylindrical coordinates, does not depend on ϕ. In such an application, one starts from the representation (4) and quickly sets both constants α2andα3to zero, so as to comply with the physical palusibility requirement tha t the stress field – and hence its trace – vanishes at infinity. Next, one does the same forα1, this time on the basis of an application of a result due to A. Si gnorini, that we now recall in a version appropriate to our present context (cf. e.g. Sect. 18 of [5]). Signorini’s Lemma . LetSbe a stress field that balances the distance and contact forces dandcacting on a domain Rwith boundary ∂R: divS+d=0inR,Sn=con∂R. 19 Moreover, let wbe a conveniently smooth vector field on R∪∂R. Then, /integraldisplay R(∇w)Sdv=/integraldisplay Rw⊗ddv+/integraldisplay ∂Rw⊗cda. (44) We specialize (4) for R ≡ B ρ, a ball of radius ρcentered at the point of application of the load, d≡0, andw=x; by taking the trace of the resulting identity, we obtain: /integraldisplay BρtrSdv=/integraldisplay ∂Bρρn·Snda. (45) Now, /integraldisplay ∂BρρSn·nda=O(ρ), because we have from (iii) that /integraldisplay ∂BρSnda=O(1). On the other hand, as to the left side of (4), we have that /integraldisplay BρtrSdv= 2π(1+ν)/integraldisplayρ 0/integraldisplayπ 0(α0cosϑ+sα1)dϑds; for it to be O(ρ) as well,α1must be set equal to 0. Acknowledgements The author gratefully acknowledges the valuable suggestio ns he received from Prof. Paolo Podio-Guidugli in the course of a number of e xtensive discussions. References [1] J.R. Barber, Elasticity . 2nd Ed., Kluwer (2002). [2] E. Benvenuto, La Scienza delle Costruzioni e il suo Sviluppo Storico . Sansoni (1981). [3] R.O. Davis and A.P.S. Selvadurai, Elasticity and Geomechanics , Cam- bridge University Press (1996). 20 [4] L.C. Evans, Partial Differential Equations , American Mathematical So- ciety, 2 edition (2010). [5] M.E. Gurtin, The Linear Theory of Elasticity , in Handbuch der Physik, Vol. VIa/2, S. Fl¨ ugge Ed., Springer-Verlag (1972). [6] A.E.H. Love, A Treatise on the Mathematical Theory of Elasticity . 4th Ed., Dover (1927). [7] L.E. Malvern, Introduction to the Mechanics of a Continuous Medium , Prentice-Hall (1969). [8] P. Podio-Guidugli, A. Favata, Elasticity for Geotechnicians , (Forthcom- ing). [9] M.H. Sadd, Elasticity. Theory, Applications, and numerics , Academic Press (2009). [10] I.S. Sokolnikoff, Mathematical Theory of Elasticity , 2nd Ed., Tata McGraw-Hill (1956). [11] M.J. Turteltaub, E. Sternberg, 1968, On concentrated L oads and Green’s functions in Elastostatics, Arch. Rat. Mech. Anal. , 29, 193-240. [12] W. Thompson(Lord Kelvin), 1848, Note on theintegratio n of theequa- tions of equilibrium of an elastic solid, Cambr. Dubl. Math. J.,3, 87-89. [13] Z.Q. Yue, H. T. Xiao, Generalized Kelvin solution based boundary element method for crack problems in multilayered solids, E ng. Anal. Bound. Elem., 26, 691-705 (2002). 21