Phil Lucht Math & Physics Archive
Home / Math and Physics Files / Physics / gases and other fluids / Sound Waves

chap-3 2D waves MIT

PDF · 54 pages · 7.0 MB
Open PDF file

Chapter 3 of MIT course notes 1.138J/2.062J/18.376J Wave Propagation by C. C. Mei, Fall 2004. It derives reflection and transmission coefficients for a plane sound wave at a fluid interface, with Snell's law, critical angle and total reflection. It then begins the equations of elastic waves, covering Lamé constants, stress and strain, and the dilatation wave equation. The reflection section appears twice in the extracted text. Source is MIT course material, not Phil's own work.

AI-written summary; may contain errors.

Extracted text (machine-read; may contain errors)
3.1.Reflection of sound by an interface 1 1.138J/2.062J/18.376J, WAVE PROPAGATION Fall, 2004 MIT Notes by C. C. Mei CHAPTER THREE TWO DIMENSIONAL WAVES 1 Reflection and tranmission of sound at an inter- face Reference : Brekhovskikh and Godin §.2.2. The governing equation for sound in a honmogeneous fluid is given by (7.31) and (7.32) in Chapter One. In term of the the veloctiy potential defined by u=∇φ (1.1) it is 1 c2∂2φ ∂t2=∇2φ (1.2) wherecdenotes the sound speed. Recall that the fluid pressure p=−ρ∂φ/∂t (1.3) also satisfies the same equation. 1.1 Plane wave in Infinite space Let us first consider a plane sinusoidal wave in three dimensional space φ(x,t)=φoei(k·x−ωt)=φoei(kn·x−ωt)(1.4) Here the phase function is θ(x,t)=k·x−ωt (1.5) 3.1.Reflection of sound by an interface 2 The equation of constant phase θ(x,t)=θodescribes a moving surface. The wave number vector k=knis defined to be k=kn=∇θ (1.6) hence is orthogonal to the surface of constant phase, and represens the direction of wave propagation. The frequency is defined to be ω=−∂θ ∂t(1.7) Is (2.40) a solution? Let us check (2.38). ∇φ=/parenleftbigg∂ ∂x,∂ ∂y,∂ ∂z/parenrightbigg φ=ikφ ∇2φ=∇·∇φ=ik·ikφ=−k2φ ∂2φ ∂t2=−ω2φ Hence (2.38) is satisfied if ω=kc (1.8) 1.2 Two-dimensional reflection from a plane interface Consider two semi-infinite fluids separated by the plane interface along z= 0. The lower fluid is distinguished from the upper fluid by the subscript ”1”. The densities and sound speeds in the upper and lower fluids are ρ,candρ1,c1respectively. Let a plane incident wave arive from z>0 at the incident angle of θwith respect to the zaxis, the sound pressure and the velocity potential are pi=P0exp[ik(xsinθ−zcosθ)] (1.9) The velocity potential is φi=−iP0 ωρexp[ik(xsinθ−zcosθ] (1.10) The indient wave number vector is ki=(ki x,ki z)=k(sinθ,−cosθ) (1.11) 3.1.Reflection of sound by an interface 3 The motion is confined in the x,zplane. On the same (incidence) side of the interface we have the reflected wave pr=Rexp[ik(xsinθ+zcosθ)] (1.12) whereRdenotes the reflection coefficient. The wavenumber vector is kr=(kr x,kr z)=k(sinθ,cosθ) (1.13) The total pressure and potential are p=P0{exp[ik(xsinθ−zcosθ)] +Rexp[ik(xsinθ+zcosθ)]} (1.14) φ=−iP0 ρω{exp[ik(xsinθ−zcosθ)] +Rexp[ik(xsinθ+zcosθ)]} (1.15) In the lower medium z<0 the transmitted wave has the pressure p1=TP 0exp[ik1(xsinθ1−zcosθ1)] (1.16) whereTis the transmission coefficient, and the potential φ1=−iP0 ρ1ωTexp[ik1(xsinθ1−zcosθ1)] (1.17) Along the interface z= 0 we require the continutiy of pressure and normal velocity, i.e., p=p1,z = 0 (1.18) and w=w1=0,z=0, (1.19) Applying (2.54), we get P0/braceleftbig eikxsinθ+Reikxsinθ/bracerightbig =TP 0eik1xsinθ1,−∞<x< ∞. Clearly we must have ksinθ=k1sinθ1 (1.20) or, sinθ c=sinθ1 c1(1.21) 3.1.Reflection of sound by an interface 4 With (2.56), we must have 1+R=T (1.22) Applying (2.55), we have iP0 ρω/bracketleftbig −kcosθeiksinθ+Rkcosθeiksinθ/bracketrightbig =iP0 ρ1ω/bracketleftbig −k1cosθ1Teik1sinθ1/bracketrightbig which implies 1−R=ρk1cosθ1 ρ1kcosθT (1.23) Eqs (2.58) and (2.59) can be solved to give T=2ρ1kcosθ ρk1cosθ1+ρ1kcosθ(1.24) R=ρ1kcosθ−ρk1cosθ1 ρ1kcosθ+ρk1cosθ1(1.25) Alternatively, we have T=2ρ1c1cosθ ρccosθ1+ρ1c1cosθ(1.26) R=ρ1c1cosθ−ρccosθ1 ρ1c1cosθ+ρccosθ1(1.27) Let m=ρ1 ρ,n =c c1(1.28) where the ratio of sound speeds nis called the index of refraction. We get after using Snell’s law that R=mcosθ−ncosθ1 mcosθ+ncosθ1=mcosθ−n/radicalBig 1−sin2θ n2 mcosθ+n/radicalBig 1−sin2θ n2(1.29) The transmission coefficient is T=1+R=2mcosθ mcosθ+n/radicalBig 1−sin2θ n2(1.30) We now examine the physics. 1. Ifn=c/c 1>1, the incidence is from a faster to a slower medium, then Ris always real. For normal incidence θ=θ1=0 , R=m−n m+n(1.31) 3.1.Reflection of sound by an interface 5 is real. Ifm>n ,0<R< 1. Ifθ=π/2, R=−n n=−1 (1.32) Hence R lies on a segment of the real axis as shown in Figur 1.a. If m<n , then R< 0 for allθas shown in figure 1.b. 2. If however n<1 thenθ1>θ. There is a critical incidence angle δ, called Brewster’s angle and defined by sinδ=n (1.33) Whenθ→δ,θ1becomesπ/2. Below this critical angle ( θ<δ ), R is real. In particular, when θ= 0, (2.67) applies. At the critical angle R=mcosδ mcosδ=1 , as shown in figure 2.c for m>n and in 2.d. for m<n . Whenθ>δ , the square roots above become imaginary. We must then take cosθ1=/radicalbigg 1−sin2θ n2=i/radicalbigg sin2θ n2−1 (1.34) This means that the reflection coefficient is now complex R=mcosθ−in/radicalBig sin2θ n2−1 mcosθ+in/radicalBig sin2θ n2−1(1.35) with|R|= 1, implying complete reflection. As a check the transmitted wave is now given by pt=Texp/bracketleftbigg k1/parenleftbigg ixsinθ1+z/radicalBig sin2θ/n2−1/parenrightbigg/bracketrightbigg (1.36) so the amplitude attenuates exponentially in zasz→− ∞ . Thus the wave train cannot penetrate much below the interface. The dependence of Ron various parameters is best displayed in the complex plane R=/RfracturR+i/IfracturR. It is clear from (2.71 ) that /IfracturR< 0 so thatRfalls on the half circle in the lower half of the complex plane as shown in figure 2.c and 2.d. 3.1.Reflection of sound by an interface 6 Figure 1: Complex reflection coefficient. From Brekhovskikh and Godin §.2.2. 2 Reflection and tranmission of sound at an inter- face Reference : Brekhovskikh and Godin §.2.2. The governing equation for sound in a honmogeneous fluid is given by (7.31) and (7.32) in Chapter One. In term of the the veloctiy potential defined by u=∇φ (2.37) it is 1 c2∂2φ ∂t2=∇2φ (2.38) wherecdenotes the sound speed. Recall that the fluid pressure p=−ρ∂φ/∂t (2.39) also satisfies the same equation. 3.1.Reflection of sound by an interface 7 2.1 Plane wave in Infinite space Let us first consider a plane sinusoidal wave in three dimensional space φ(x,t)=φoei(k·x−ωt)=φoei(kn·x−ωt)(2.40) Here the phase function is θ(x,t)=k·x−ωt (2.41) The equation of constant phase θ(x,t)=θodescribes a moving surface. The wave number vector k=knis defined to be k=kn=∇θ (2.42) hence is orthogonal to the surface of constant phase, and represens the direction of wave propagation. The frequency is defined to be ω=−∂θ ∂t(2.43) Is (2.40) a solution? Let us check (2.38). ∇φ=/parenleftbigg∂ ∂x,∂ ∂y,∂ ∂z/parenrightbigg φ=ikφ ∇2φ=∇·∇φ=ik·ikφ=−k2φ ∂2φ ∂t2=−ω2φ Hence (2.38) is satisfied if ω=kc (2.44) 2.2 Two-dimensional reflection from a plane interface Consider two semi-infinite fluids separated by the plane interface along z= 0. The lower fluid is distinguished from the upper fluid by the subscript ”1”. The densities and sound speeds in the upper and lower fluids are ρ,candρ1,c1respectively. Let a plane incident wave arive from z>0 at the incident angle of θwith respect to the zaxis, the sound pressure and the velocity potential are pi=P0exp[ik(xsinθ−zcosθ)] (2.45) 3.1.Reflection of sound by an interface 8 The velocity potential is φi=−iP0 ωρexp[ik(xsinθ−zcosθ] (2.46) The indient wave number vector is ki=(ki x,ki z)=k(sinθ,−cosθ) (2.47) The motion is confined in the x,zplane. On the same (incidence) side of the interface we have the reflected wave pr=Rexp[ik(xsinθ+zcosθ)] (2.48) whereRdenotes the reflection coefficient. The wavenumber vector is kr=(kr x,kr z)=k(sinθ,cosθ) (2.49) The total pressure and potential are p=P0{exp[ik(xsinθ−zcosθ)] +Rexp[ik(xsinθ+zcosθ)]} (2.50) φ=−iP0 ρω{exp[ik(xsinθ−zcosθ)] +Rexp[ik(xsinθ+zcosθ)]} (2.51) In the lower medium z<0 the transmitted wave has the pressure p1=TP 0exp[ik1(xsinθ1−zcosθ1)] (2.52) whereTis the transmission coefficient, and the potential φ1=−iP0 ρ1ωTexp[ik1(xsinθ1−zcosθ1)] (2.53) Along the interface z= 0 we require the continutiy of pressure and normal velocity, i.e., p=p1,z = 0 (2.54) and w=w1=0,z=0, (2.55) Applying (2.54), we get P0/braceleftbig eikxsinθ+Reikxsinθ/bracerightbig =TP 0eik1xsinθ1,−∞<x< ∞. 3.1.Reflection of sound by an interface 9 Clearly we must have ksinθ=k1sinθ1 (2.56) or, sinθ c=sinθ1 c1(2.57) With (2.56), we must have 1+R=T (2.58) Applying (2.55), we have iP0 ρω/bracketleftbig −kcosθeiksinθ+Rkcosθeiksinθ/bracketrightbig =iP0 ρ1ω/bracketleftbig −k1cosθ1Teik1sinθ1/bracketrightbig which implies 1−R=ρk1cosθ1 ρ1kcosθT (2.59) Eqs (2.58) and (2.59) can be solved to give T=2ρ1kcosθ ρk1cosθ1+ρ1kcosθ(2.60) R=ρ1kcosθ−ρk1cosθ1 ρ1kcosθ+ρk1cosθ1(2.61) Alternatively, we have T=2ρ1c1cosθ ρccosθ1+ρ1c1cosθ(2.62) R=ρ1c1cosθ−ρccosθ1 ρ1c1cosθ+ρccosθ1(2.63) Let m=ρ1 ρ,n =c c1(2.64) where the ratio of sound speeds nis called the index of refraction. We get after using Snell’s law that R=mcosθ−ncosθ1 mcosθ+ncosθ1=mcosθ−n/radicalBig 1−sin2θ n2 mcosθ+n/radicalBig 1−sin2θ n2(2.65) The transmission coefficient is T=1+R=2mcosθ mcosθ+n/radicalBig 1−sin2θ n2(2.66) We now examine the physics. 3.1.Reflection of sound by an interface 10 1. Ifn=c/c 1>1, the incidence is from a faster to a slower medium, then Ris always real. For normal incidence θ=θ1=0 , R=m−n m+n(2.67) is real. Ifm>n ,0<R< 1. Ifθ=π/2, R=−n n=−1 (2.68) Hence R lies on a segment of the real axis as shown in Figur 1.a. If m<n , then R< 0 for allθas shown in figure 1.b. 2. If however n<1 thenθ1>θ. There is a critical incidence angle δ, called Brewster’s angle and defined by sinδ=n (2.69) Whenθ→δ,θ1becomesπ/2. Below this critical angle ( θ<δ ), R is real. In particular, when θ= 0, (2.67) applies. At the critical angle R=mcosδ mcosδ=1 , as shown in figure 2.c for m>n and in 2.d. for m<n . Whenθ>δ , the square roots above become imaginary. We must then take cosθ1=/radicalbigg 1−sin2θ n2=i/radicalbigg sin2θ n2−1 (2.70) This means that the reflection coefficient is now complex R=mcosθ−in/radicalBig sin2θ n2−1 mcosθ+in/radicalBig sin2θ n2−1(2.71) with|R|= 1, implying complete reflection. As a check the transmitted wave is now given by pt=Texp/bracketleftbigg k1/parenleftbigg ixsinθ1+z/radicalBig sin2θ/n2−1/parenrightbigg/bracketrightbigg (2.72) so the amplitude attenuates exponentially in zasz→− ∞ . Thus the wave train cannot penetrate much below the interface. The dependence of Ron various parameters is best displayed in the complex plane R=/RfracturR+i/IfracturR. It is clear from (2.71 ) that /IfracturR< 0 so thatRfalls on the half circle in the lower half of the complex plane as shown in figure 2.c and 2.d. 3.2.Equations for Elastic Waves 11 Figure 2: Complex reflection coefficient. From Brekhovskikh and Godin §.2.2. 3 Equations for elastic waves Refs: Graff: Wave Motion in Elastic Solids Aki & Richards Quantitative Seismology, V. 1. Achenbach. Wave Propagation in Elastic Solids Let the displacement vector at a point xjand timetbe denoted by ui(xj,t), then Newton’s law applied to an material element of unit volume reads ρ∂2ui ∂t2=∂τij ∂xj(3.1) whereτijis the stress tensor. We have neglected body force such as gravity. For a homogeneous and isotropic elastic solid, we have the following relation between stress and strain τij=λekkδij+2µeij (3.2) whereλandµare Lam´ e constants and eij=1 2/parenleftbigg∂ui ∂xj+∂uj ∂xi/parenrightbigg (3.3) 3.2.Equations for Elastic Waves 12 is the strain tensor. Eq. (3.2) can be inverted to give eij=1+ν Eτij−ν Eτkkδij (3.4) where E=µ(3λ+µ) λ+µ(3.5) is Young’s modulus and ν=λ 2(λ+µ). (3.6) Poisson’s ratio. Substituting (3.2) and (3.3) into (3.1) we get ∂τij ∂xj=λ∂ekk ∂xjδij+µ∂ ∂xj/parenleftbigg∂ui ∂xj+∂uj ∂xi/parenrightbigg =λ∂ekk ∂xi+µ∂2ui ∂xj∂xj+µ∂2uj ∂xi∂xj =(λ+µ)∂2ui ∂xixj+µ∇2ui In vector form (3.1) becomes ρ∂2u ∂t2=(λ+µ)∇(∇·u)+µ∇2u (3.7) Taking the divergence of (3.1) and denoting the dilatation by ∆≡ekk=∂u1 ∂x1+∂u2 ∂x2+∂u3 ∂x3(3.8) we get the equation governing the dilatation alone ρ∂2∆ ∂t2=(λ+µ)∇·∇ ∆+µ∇2∆=(λ+2µ)∇2∆ (3.9) or, ∂2∆ ∂t2=c2 L∇2∆ (3.10) where cL=/radicalBigg λ+2µ ρ(3.11) 3.2.Equations for Elastic Waves 13 Thus the dilatation propagates as a wave at the speed cL. To be explained shortly, this is a longitudinal waves, hence the subscript L. On the other hand, taking the curl of (3.7) and denoting by /vectorωthe rotation vector: /vectorω=∇×u (3.12) we then get the governing equation for the rotation alone ∂2/vectorω ∂t2=c2 T∇2/vectorω (3.13) where cT=/radicalbigg µ ρ(3.14) Thus the rotation propagates as a wave at the slower speed cT. The subscript Tindicates that this is a transverse wave, to be shown later. The ratio of two wave speeds is cL cT=/radicalBigg λ+µ µ>1. (3.15) Since µ λ=1 2ν−1 (3.16) it follows that the speed ratio depends only on Poisson’s ratio cL cT=/radicalbigg 2−2ν 1−2ν(3.17) There is a general theorem due to Helmholtz that any vector can be expressed as the sum of an irrotational vector and a solenoidal vector i.e., u=∇φ+∇×H (3.18) subject to the constraint that ∇·H= 0 (3.19) The scalar φand the vector Hare called the displacement potentials. Substituting this into (3.7), we get ρ∂2 ∂t2[∇φ+∇×H]=µ∇2[∇φ+∇×H]+(λ+µ)∇∇ · [∇φ+∇×H] 3.3. Free waves in infinite space 14 Since∇·∇φ=∇2φ,and∇·∇× H=0w eg e t ∇/bracketleftbigg (λ+2µ)∇2φ−ρ∂2φ ∂t2/bracketrightbigg +∇×/bracketleftbigg µ∇2H−ρ∂2H ∂t2/bracketrightbigg = 0 (3.20) Clearly the above equation is satisied if (λ+2µ)∇2φ−ρ∂2φ ∂t2= 0 (3.21) and µ∇2H−ρ∂2H ∂t2= 0 (3.22) Although the governing equations are simplified, the two potentials are usually coupled by boundary conditions, unless the physical domain is infinite. 4 Free waves in infinite space The dilatational wave equation admits a plane sinusoidal wave solution: φ(x,t)=φoeik(n·x−cLt)(4.1) Here the phase function is θ(x,t)=k(n·x−cLt) (4.2) which describes a moving surface. The wave number vector k=knis defined to be k=kn=∇θ (4.3) hence is orthogonal to the surface of constant phase, and represents the direction of wave propagation. The frequency is ω=kcT=−∂θ ∂t(4.4) A general solution is φ=φ(n·x−cLt) (4.5) Similarly the the following sinusoidal wave is a solution to the shear wave equation; H=Hoeik(n·x−cTt)(4.6) 3.4. Elastic waves in a plane 15 A general solution is H=H(n·x−cTt) (4.7) −−−−−−− Note: We can also write (4.5) and (4.9) as φ=φ(t−n·x cL) (4.8) and H=H(t−n·x cT) (4.9) where sL=n cL,sT=n cT(4.10) are called the slowness vectors of longitudinal and transverse waves respectively. −−−−− In a dilatational wave the displacement vector is parallel to the wave number vector: uL=∇φ=φ/primen (4.11) from (3.5), where φ/primeis the ordinary derivative of φwith repect to its argument. Hence the dilatational wave is a longitudinal (compression) wave. On the other hand in a rotational wave the displacement vector is perpendicular to the wave number vector, uT=∇×H=ex/parenleftbigg∂Hz ∂y−∂Hy ∂z/parenrightbigg +ey/parenleftbigg∂Hx ∂z−∂Hz ∂x/parenrightbigg +ez/parenleftbigg∂Hy ∂x−∂Hx ∂y/parenrightbigg =ex/parenleftbig H/prime zny−H/prime ynz/parenrightbig +ey(H/prime xnz−H/prime znx)+ez/parenleftbig H/prime ynx−H/prime xny/parenrightbig =n×H/prime(4.12) from (3.7). Hence a rotational wave is a transverse (shear) wave. 4 Elastic waves in a plane Refs. Graff, Achenbach, Aki and Richards : Quantitative Seismology, v.1 3.4. Elastic waves in a plane 16 Let us examine waves propagating in the vertical plane of x,y. All physical quantities are assumed to be uniform in the direction of z, hence∂/∂z = 0, then ux=∂φ ∂x+∂Hz ∂y,u y=∂φ ∂y−∂Hz ∂x,u z=−∂Hx ∂y+∂Hy ∂x(4.13) and ∂Hx ∂x+∂Hy ∂y= 0 (4.14) where ∂2φ ∂x2+∂2φ ∂y2=1 c2 L∂2φ ∂t2, (4.15) ∂2Hp ∂x2+∂2Hp ∂y2=1 c2 T∂2Hp ∂t2,p=x,y,z (4.16) Note thatuzis also governed by (4.16). Note that the in-plane displacements ux,uydepend only on φandHz, and not on Hx,Hy. Out-of-plane motion uzdepends on Hx,Hybut not on Hz. Hence the in- plane displacement components ux,uyare independent of the out-of-plane component uz. The in-plane displacements ( ux,uy) are associated with dilatation and in-plane shear, represented respectively by φandHz, which will be refered to as the P wave and the SV wave. The out-of-plane displacement uzis associated with HxandHy, and will be refered to as the SH wave. From Hooke’s law the stress components can be written τxx=λ/parenleftbigg∂ux ∂x+∂uy ∂y/parenrightbigg +2µ∂ux ∂x=(λ+2µ)/parenleftbigg∂ux ∂x+∂uy ∂y/parenrightbigg −2µ∂uy ∂y =(λ+2µ)/parenleftbigg∂2φ ∂x2+∂2φ ∂y2/parenrightbigg −2µ/parenleftbigg∂2φ ∂y2−∂2Hz ∂y∂x/parenrightbigg (4.17) τyy=λ/parenleftbigg∂ux ∂x+∂uy ∂y/parenrightbigg +2µ∂uy ∂y=(λ+2µ)/parenleftbigg∂ux ∂x+∂uy ∂y/parenrightbigg −2µ∂ux ∂x =(λ+2µ)/parenleftbigg∂2φ ∂x2+∂2φ ∂y2/parenrightbigg −2µ/parenleftbigg∂2φ ∂x2+∂2Hz ∂x∂y/parenrightbigg (4.18) τzz=λ 2(λ+µ)(τxx+τyy)=ν(τxx+τyy)=λ/parenleftbigg∂2φ ∂x2+∂2φ ∂y2/parenrightbigg (4.19) τxy=µ/parenleftbigg∂uy ∂x+∂ux ∂y/parenrightbigg =µ/parenleftbigg 2∂2φ ∂x∂y−∂2Hz ∂x2+∂2Hz ∂y2/parenrightbigg (4.20) 3.5. Reflection of elastic waves from a plane boundary 17 τyz=µ∂uz ∂y=µ/parenleftbigg −∂2Hx ∂y2+∂2Hy ∂y∂x/parenrightbigg (4.21) τxz=µ∂uz ∂x=µ/parenleftbigg −∂2Hx ∂∂x∂y+∂2Hy ∂x2/parenrightbigg (4.22) Different physical situations arise for different boundary conditions. We shall con- sider first the half plane problem bounded by the plane y=0 . 5 Reflection of elastic waves from a plane boundary Consider the half space y> 0,−∞<x< ∞. Several types of boundary conditions can be prescribed on the plane boundary : (i) dynamic: the stress components only (the traction condition); (ii) kinematic: the displacement components only, or (iii). a combination of stress components and displacement components. Most difficult are (iv) the mixed conditions in which stresses are given over part of the boundary and displacements over the other. We consider the simplest case where the plane y= 0 is completely free of external stresses, τyy=τxy=0, (5.23) and τyz= 0 (5.24) It is clear that (5.23) affects the P and SV waves only, while (5.24) affects the SH wave only. Therefore we have two uncoupled problems each of which can be treated separately. 5.1 P and SV waves Consider the case where only PandSVwaves are present, then Hx=Hy= 0. Let all waves have wavenumber vectors inclined in the positve xdirection: φ=f(y)eiξx−iωt,H z=hz(y)eiζx−iωt(5.25) It follows from (4.15) and (4.16) that d2f dy2+α2f=0,d2hz dy2+β2hz=0, (5.26) 3.5. Reflection of elastic waves from a plane boundary 18 with α=/radicalBigg ω2 c2 L−ξ2=/radicalBig k2 L−ξ2,β =/radicalBigg ω2 c2 T−ζ2=/radicalBig k2 T−ζ2 (5.27) We first take the square roots to be real; the general solution to (5.26) are sinusoids, hence, φ=APei(ξx−αy−ωt)+BPei(ξx+αy−ωt),H z=ASei(ζx−βy−ωt)+BSei(ζx+βy−ωt)(5.28) On the right-hand sides the first terms are the incident waves and the second are the reflected waves. If the incident amplitudes AP,ASand are given, what are the properties of the reflected waves BP,BS? The wave number components can be written in the polar form: (ξ,α)=kL(sinθL,cosθL),(ζ,β)=kT(sinθT,cosθT) (5.29) where (kL,kT) are the wavenumbers, the ( θL,θT) the directions of the P wave and SV wave, respectively. In terms of these we rewrite (5.28) φ=APeikL(sinθLx−cosθLy−ωt)+BPeikL(sinθLx+cos θLy−ωt)(5.30) Hz=ASeikT(sinθTx−cosθLy−ωt)+BSeikT(sinθTx+cos θTy−ωt)(5.31) In order to satisfy (5.23) (τ yy=τxy= 0) ony= 0 for allx, we must insist: kLsinθL=kTsinθT,(ξ=ζ) (5.32) This is in the form of Snell’s law: sinθL cL=sinθT cT(5.33) implying sinθL sinθT=cL cT=/radicalBigg λ+2µ µ=kT kL≡κ (5.34) When (5.23) are applied on y= 0 the exponential factors cancel, and we get two algebraic conditions for the two unknown amplitudes of the reflected waves ( BP,BS): k2 L(2 sin2θL−κ2)(AP+BP)−k2 Tsin 2θT(AS−BS) = 0 (5.35) k2 Lsin 2θL(AP−BP)−k2 TcosθT(AS+BS)=0. (5.36) 3.5. Reflection of elastic waves from a plane boundary 19 Using (5.34), we get 2 sin2θL−κ2=κ2(2 sin2θT−1) =−κ2cos 2θT The two equations can be solved and the solution expressed in matrix form:   BP BS  = SPPSSP SPSSSS   AP AS  (5.37) where S= SPPSSP SPSSSS  (5.38) denotes the scattering matrrix. Thus SPSrepresents the reflected S-wave due to incident Pwave of unit amplitude, etc. It is straightforward to verify that SPP=sin 2θLsin 2θT−κ2cos22θT sin 2θLsin 2θT+κ2cos22θT(5.39) SSP=−2κ2sin 2θTcos 2θT sin 2θLsin 2θT+κ2cos22θT(5.40) SPS=2 sin 2θLcos 2θT sin 2θLsin 2θT+κ2cos22θT(5.41) SSS=sin 2θLsin 2θT−κ2cos22θT sin 2θLsin 2θT+κ2cos22θT(5.42) In view of (5.33) and κ=cL cT=/radicalbigg 2−2ν 1−2ν(5.43) The scattering matrix is a function of Poisson’s ratio and the angle of incidence. (i) P- wave Incidence : In this case θLis the incidence angle. Consider the special case when the only incident wave is a P wave. Then AP/negationslash= 0 andAS= 0 and only SPPandSSPare relevant. . Note first that θL>θ Tin general . For normal incidence, θL= 0, hence θT= 0. We find SPP=−1,S PS= 0 (5.44) there is no SV wave. The refelcted wave is a P wave. On the other hand if sin 2θLsin 2θT−κ2cos22θT= 0 (5.45) 3.5. Reflection of elastic waves from a plane boundary 20 thenSPP= 0, hence BP= 0 butBS/negationslash= 0; onlySVwave is reflected. This is the case of mode conversion, whereby an incident P waves changes to a SV wave after reflection. The amplitude of the reflected SV wave is BS AP=SPS=tan 2θT κ2(5.46) (ii) SV wave Incidence : Let AP= 0 butAS/negationslash= 0. In this case θTis the incidence angle. Then only SSPandSSSare relevant. For normal incidence, θL=θT=0 , SSS=−1, andSSP= 0; no P wave is reflected. Mode conversion ( BP/negationslash=0,BS= 0) also happens when (5.45) is satisfied. Since θL>θ T, there is a critical incidence angle θT beyond which the P wave cannot be reflected back into the solid and propagates only along the x axis. At the critical angle sinθL=1,(θL=π/2),or sinθT=1/κ (5.47) by Snell’s law. Thus for ν=1/3,κ= 2 and the critical incidence angle is θT=3 0◦. Beyond the critical angle of incidence, the Pwaves decay exponentially away from the free surface. The amplitude of the SV wave is linear in ywhich is unphysical, suggesting the limitation of unbounded space assumption. 5.2 SH wave Because of (5.2) ∂Hx ∂x+∂Hy ∂y=0 we can introduce a stream function ψso that Hx=−∂ψ ∂y,H y=∂ψ ∂x(5.48) where ∇2ψ=1 c2 T∂2ψ ∂t2(5.49) Clearly the out-of-plane dispacement is uz=−∂Hx ∂y+∂Hy ∂x=∂2ψ ∂x2+∂2ψ ∂y2=∇2ψ (5.50) 3.5. Reflection of elastic waves from a plane boundary 21 Figure 3: Amplitude ratios for incident P waves for various Possion’s ratios. From Graff: Waves in Elastic Solids . Symbols should be converted according to : A1→AP,A2→ BP,B1→AS,B2→BS. 3.5. Reflection of elastic waves from a plane boundary 22 Figure 4: Reflected wave amplitude ratios for incident SV waves for various Possion’s ratios. From Graff: Waves in Elastic Solids . Symbols should be converted according to :A1→AP,A2→BP,B1→AS,B2→BS. 3.6. Scattering of SH waves by a circular cavity 23 and τyz=µ∂ ∂y∇2ψ=µ c2 T∂ ∂y∂2ψ ∂t2(5.51) The zero-stress boundary condition implies ∂ψ ∂y= 0 (5.52) Thus the problem for ψis analogous to one for sound waves reflected by a solid plane. Again for monochromatic incident waves, the solution is easily shown to be ψ=/parenleftbig Ae−iβy−Aeiβy/parenrightbig eiαx−iωt(5.53) where α2+β2=k2 T (5.54) We remark that when the boundary is any cylindrical surface with axis parallel to thezaxis, the the stress-free condition reads τzn=0,onB. (5.55) wherenis the unit outward normal to B. Since in the pure SH wave problem τzn=µ∂uz ∂n=∂ ∂n∇2ψ=µ c2 T∂ ∂n∂2ψ ∂t2 Condition (5.55) implies ∂ψ ∂n=0,onB. (5.56) Thus the analogy to acoustic scattering by a hard object is true irrespective of the geomntry of the scatterer. 6 Scattering of monochromatic SH waves by a cir- cular cavity 6.1 Solution in polar coordinates We consider the scattering of two-dimensional SH waves of single frequency. The time- dependent potential can be wirtten as ψ(x,y,t )=/Rfractur/bracketleftbig φ(x,y)e−iωt/bracketrightbig (6.1) 3.6. Scattering of SH waves by a circular cavity 24 where the potential φis governed by the Helmholtz equation ∇2φ+k2φ=∂2φ ∂x2+∂2φ2 ∂y2+k2φ=0,k =ω cT(6.2) To be specific consider the scatterer to be a finite cavity of some general geometry. On the stress-free boundary Bthe shear stress vanishes, τzn=−µω2 c2 T/Rfractur/parenleftbigg∂φ ∂ne−iωt/parenrightbigg = 0 (6.3) hence ∂φ ∂n= 0 (6.4) Let the incident waves be a plane wave φI=Aeik·x(6.5) and the angle of incidence is θowith respect to the positive xaxis. In polar coordinates we write k=k(cosθo,sinθo),x=r(cosθ,sinθ) (6.6) φI=Aexp [ikr(cosθocosθ+ sinθosinθ)] =Aeikrcos(θ−θo)(6.7) It can be shown (see Appendix A) that the plane wave can be expanded in Fourier-Bessel series : eikrcos(θ−θo)=∞/summationdisplay n=0/epsilon1ninJn(kr) cosn(θ−θo) (6.8) where/epsilon1nis the Jacobi symbol: /epsilon10=0,/epsilon1 n=2,n=1,2,3,... (6.9) Each term in the series (6.8) is called a partial wave. Let the total wave be the sum of the incident and scattered waves φ=φI+φS (6.10) then the scattered waves must satisfy the radiation condition at infinity, i.e., it can only radiate energy outward from the scatterer. The boundary condition on the cavity surface is ∂φ ∂r=0,r =a (6.11) 3.6. Scattering of SH waves by a circular cavity 25 In polar coordinates the governing equation reads 1 r∂ ∂r/parenleftbigg r∂φ ∂r/parenrightbigg +1 r2∂2φ ∂θ2+k2φ= 0 (6.12) SinceφIsatisfies the preceding equation, so does φS. By the method of separation of variables, φS(r,θ)=R(r)Θ(θ) we find r2R/prime/prime+rR/prime+(k2r2−n2)R=0,and Θ/prime/prime+n2Θ=0 wheren=0,1,2,... are eigenvalues in order that Θ is periodic in θwith period 2 π. For each eigenvalue nthe possible solutions are Θn= (sinnθ,cosnθ), Rn=/parenleftbig H(1) n(kr),H(2) n(kr)/parenrightbig , whereH(1) n(kr),H(2) n(kr) are Hankel functions of the first and second kind, related to the Bessel and Weber functions by H(1) n(kr)=Jn(kr)+iYn(kr),H(2) n(kr)=Jn(kr)−iYn(kr) (6.13) The most general solution to the Helmholtz equation is φS=A∞/summationdisplay n=0(Ansinnθ+Bncosnθ)/bracketleftbig CnH(1) n(kr)+DnH(2) n(kr)/bracketrightbig , (6.14) For large radius the asymptotic form of the Hankel functions are H(1) n∼/radicalbigg 2 πkrei(kr−π 4−nπ 2),H(2) n∼/radicalbigg 2 πkre−i(kr−π 4−nπ 2)(6.15) In conjunction with the time factor exp( −iωt),H(1) ngives an outgoing wave while H(2) n gives an incoming wave. To satisfy the radiation condition, we must discard all terms involvingH(2) n. From here on we shall abbreviate H(1) nsimply byHn. The scattered wave is now φS=A∞/summationdisplay n=0(Ansinnθ+Bncosnθ)Hn(kr) (6.16) 3.6. Scattering of SH waves by a circular cavity 26 The expansion coefficients ( An,Bn) must be chosen to satisfy the boundary condition on the cavity surface1Once they are determined, the wave is found everywhere. In particular in the far field, we can use the asymptotic formula to get φS∼A∞/summationdisplay n=0(Ansinnθ+Bncosnθ)e−inπ/ 2/radicalbigg 2 πkreikr−iπ/4(6.17) Let us define the dimensionless directivity factor A(θ)=∞/summationdisplay n=0(Ansinnθ+Bncosnθ)e−inπ/ 2(6.18) which indicates the angular variation of the far-field amplitude, then φS∼AA(θ)/radicalbigg 2 πkreikr−iπ/4(6.19) This expression exhibits clearly the asymptotic behaviour of φSas an outgoing wave. By differentiation, we readily see that lim kr→∞√ r/parenleftbigg∂φS ∂r−φS/parenrightbigg = 0 (6.20) which is one way of stating the radiation condition for two dimensional SH waves. At any radius rthe total rate of energy outflux by the scattered wave is r/integraldisplay2π 0dθ τrz∂uz ∂t=µr/integraldisplay2π 0dθ /Rfractur/bracketleftbigg −µk2∂φ ∂re−iωt/bracketrightbigg /Rfractur[iωk2φe−iωt] =−µωk4r 2/integraldisplay2π 0dθ/Rfractur/bracketleftbigg iφ∗∂φ ∂r/bracketrightbigg =−µωk4r 2/Ifractur/integraldisplay2π 0dθ/bracketleftbigg φ∗∂φ ∂r/bracketrightbigg (6.21) where overline indicates time averaging over a wave period 2 π/ω. We remark that in the analogous case of plane acoustics where the sound pressure and radial fluid velocity are respectively, p=−ρo∂φ ∂t,andur=∂φ ∂r(6.22) the energy scattering rate is r/integraldisplay∞ 0dθ pur=ωρor 2/Rfractur/integraldisplay Cdθ/parenleftbigg −iφ∗∂φ ∂r/parenrightbigg =−ωρor 2/Ifractur/integraldisplay Cdθ/parenleftbigg φ∗∂φ ∂r/parenrightbigg (6.23) 1In one of the numerical solution techniques, one divides the physical region by a circle enclosing the cavity. Between the cavity and the circle, finite elements are used. Outside the circle, (6.16) is used. By constructing a suitable variational principle, finite element computation yields the nodal coefficients as well as the expansion coefficients. See (Chen & Mei , 1974). 3.6. Scattering of SH waves by a circular cavity 27 6.2 A general theorem on scattering For the same scatterer and the same frequency ω, different angles of incidence θjdefine different scattering problems φj. In particular at infinty, we have φj∼Aj/braceleftBigg eikrcos(θ−θj)+Aj(θ)/radicalbigg 2 πkreikr−iπ/4/bracerightBigg (6.24) Let us apply Green’s formula to φ1andφ2over a closed area bounded by a closed contourC,/integraldisplay/integraldisplay S/parenleftbig φ2∇2φ1−φ1∇2φ2/parenrightbig dA=/integraldisplay B/parenleftbigg φ2∂φ1 ∂n−φ1∂φ2 ∂n/parenrightbigg ds+/integraldisplay Cds/parenleftbigg φ2∂φ1 ∂n−φ1∂φ1 ∂n/parenrightbigg ds where nrefers to the unit normal vector pointing out of S. The surface integral vanishes on account of the Helmholtz equation, while the line integral along the cavity surface vanishes by virture of the boundary condition, hence/integraldisplay Cds/parenleftbigg φ2∂φ1 ∂n−φ1∂φ2 ∂n/parenrightbigg ds= 0 (6.25) By similar reasioning, we get/integraldisplay Cds/parenleftbigg φ2∂φ∗ 1 ∂n−φ∗ 1∂φ2 ∂n/parenrightbigg ds= 0 (6.26) whereφ∗ 1denotes the complex conjugate of φ1. Let us choose φ1=φ2=φin (6.26), and get /integraldisplay Cds/parenleftbigg φ∂φ∗ ∂n−φ∗∂φ ∂n/parenrightbigg ds=2/Ifractur/parenleftbigg/integraldisplay Cdsφ∂φ∗ ∂n/parenrightbigg = 0 (6.27) Physically, across any circle the net rate of energy flux vanishes, i.e., the scattered power must be balanced by the incident power. Making use of (6.24) we get 0=/Ifractur/integraldisplay2π 0rdθ/bracketleftBigg eikrcos(θ−θo)+/radicalbigg 2 πkrAo(θ)eikr−iπ/4/bracketrightBigg ·/bracketleftBigg −ikcos(θ−θo)e−ikrcos(θ−θo)−ik/radicalbigg 2 πkrA∗ o(θ)e−ikr+iπ/4/bracketrightBigg =/Ifractur/integraldisplay2π 0rdθ/braceleftbigg −ikcos(θ−θo)+2 πkr(−ik)|Ao|2 +eikr[cosθ−θo)−1]+iπ/4(−ik)/radicalbigg 2 πkrA∗ o +e−ikr[cosθ−θo)−1]−iπ/4(−ik) cos(θ−θo)/radicalbigg 2 πkrAo/bracerightBigg 3.6. Scattering of SH waves by a circular cavity 28 The first term in the integrand gives no contribution to the integral above because of periodicity. Since /Ifractur(if)=/Ifractur(if∗), we get 0=−2 π/integraldisplay2π 0|Ao(θ)|2dθ +/Ifractur/integraldisplay2π 0rdθ/braceleftBigg Ao(−ik)/radicalbigg 2 πkr[1 + cos(θ−θo)]eiπ/4eikr(1−cos(θ−θo))/bracerightBigg =−2 π/integraldisplay2π 0|Ao(θ)|2dθ −/Rfractur/braceleftBigg e−iπ/4/bracketleftBigg Ao(k)r/radicalbigg 2 πkr/integraldisplay2π 0dθ[1 + cos(θ−θo)]eikr(1−cos(θ−θo))/bracketrightBigg/bracerightBigg For largekrthe remaining integral can be found approximately by the method of sta- tionary phase (see Appendix B), with the result /integraldisplay2π 0dθ[1 + cos(θ−θo)]eikr(1−cos(θ−θo))∼/radicalbigg 2π kreiπ/4(6.28) We get finally/integraldisplay2π 0|A|2dθ=−2/RfracturA(θo) (6.29) Thus the total scattered energy in all directions is related to the amplitude of the scattered wave in the forward direction. In atomic physics, where this theorem was originated (by Niels Bohr), measurement of the scattering amplitude in all directions is not easy. This theorem suggests an econmical alternative. Homework For the same scatterer, clonsider two scattering problems φ1andφ2. Show that A1(θ2)=A2(θ1) (6.30) For general elastic waves, see Mei (1978) : Extensions of some identities in elastody- namics with rigid inclusions. J . Acoust. Soc. Am. 64(5), 1514-1522. 6.3 Scattering by a circular cavity continued Without loss of generality we can take θo= 0. On the surface of the cylindrical cavity r=a, we impose ∂φI ∂r+∂φS ∂r=0,r=a 3.6. Scattering of SH waves by a circular cavity 29 It follows that An= 0 and /epsilon1ninAJ/prime n(ka)+BnkH/prime n(ka)=0,n =0,1,2,3,...n where primes denote differentiation with respect to the argument. Hence Bn=−A/epsilon1ninJ/prime n(ka) H/prime n(ka) The sum of incident and scattered waves is φ=A∞/summationdisplay n=0enin/bracketleftbigg Jn(kr)−J/prime n(ka) H/primen(ka)Hn(kr)/bracketrightbigg cosnθ (6.31) and ψ=Ae−iωt∞/summationdisplay n=0enin/bracketleftbigg Jn(kr)−J/prime n(ka) H/primen(ka)Hn(kr)/bracketrightbigg cosnθ (6.32) The limit of long waves can be approximatedly analyzed by using the expansions for Bessel functions for small argument Jn(x)∼xn 2nn!,Y n(x)∼2 πlogx, Y n(x)∼2n(n−1)! πxn(6.33) Then the scattered wave has the potential φS A∼−H0(kr)J/prime 0(ka) H/prime 0(ka)−2iH1(kr)J/prime 1(ka) H/prime 1(ka)cosθ+O(ka)3 =π 2(ka)2/parenleftbigg −i 2H0(kr)−H1(kr) cosθ/parenrightbigg +O(ka)3(6.34) The termH0(kr) coresponds to a oscillating source which sends istropic waves in all directions. The second term is a dipole sending scattered waves mostly in forward and backward directions. For large kr, the angular variation is a lot more complex. The far field pattern for various kais shown in fig 4. On the cavity surface surface, the displacement is proportional to ψ(a,θ)o rφ(a,θ). The angular variation is plotted for several kain figure 5. 3.6. Scattering of SH waves by a circular cavity 30 Figure 5: Angular distribution of scattered energy in the far field in cylindrical scattering Figure 6: Polar distribution of φ(a,θ) on a circular cylinder. 3.7. Diffraction of SH waves by a long crack 31 7 Diffraction of SH wave by a long crack References Morse & Ingard, Theoretical Acoustics Series expansions. Born & Wolf, Principle of Optics Fourier Transform and the method of steepest descent. B. Noble. The Wiener-Hopf Technique . If the obstacle is large, there is always a shadow behind where the incident wave cannot penetrate deeply. The phenomenon of scattering by large obstacles is usually refered to as diffraction. Diffraction of plane incident SH waves by a long crack is identical to that of a hard screen in acoustics. The exact solution was due to A. Sommerfeld. We shall apply the boundary layer idea and give the approximate solution valid far away from the tip kr/greatermuch1b yt h e parabolic approximation , due to V. Fock. Refering to figure () let us make a crude division of the entire field into the illuminated zone I , dominated by the incident wave alone, the reflection zone II dominated the sumof the incident and the reflected wave, and the shadow zone III where there is no wave. The boundaries of these zones are the rays touching the crack tip. According to this crude picture the solution is φ=  A oexp(ikcosθx+iksinθy),I Ao[exp(ikcosθx+iksinθy)+e x p (ikcosθx−iksinθy)],I I 0, III(7.1) Clearly (7.1) is inadquate because the potential cannot be discontinuous across the boundaries. A remedy to provide smooth transitions is needed. Consider the shadow boundary Ox/prime. Let us introduce a new cartesian coordinate system so that x/primeaxis is along, while the y/primeaxis is normal to, the shadow boundary. The relations between ( x,y) and (x/prime,y/prime) are x/prime=xcosθ+ysinθ, y/prime=ycosθ−xsinθ (7.2) Thus the incident wave is simply φI=Aoeikx/prime(7.3) 3.7. Diffraction of SH waves by a long crack 32 □□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□ □□□□□□□□□□□□□□□□□□□□□□ □□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□ □□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□ □□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□□ Figure 7: Wave zones near a long crack Following the chain rule of differentiation, ∂φ ∂x=∂φ ∂x/prime∂x/prime ∂x+∂φ ∂y/prime∂y/prime ∂x= cosθ∂φ ∂x/prime−sinθ∂φ ∂y/prime ∂φ ∂y=∂φ ∂x/prime∂x/prime ∂y+∂φ ∂y/prime∂y/prime ∂y= sinθ∂φ ∂x/prime+ cosθ∂φ ∂y/prime we can show straightforwardly that ∂2φ ∂x2+∂2φ ∂y2=∂2φ ∂x/prime2+∂2φ ∂y/prime2 so that the Helmholtz equation is unchanged in form in the x/prime,y/primesystem. We try to fit a boundary layer along the x’ axis and expect the potential to be almost like a plane wave φ(x,/prime,y/prime)=A(x/prime,y/prime)eikx/prime(7.4) , but the amplitude is slowly modulated in both x/primeandy/primedirections. Substituting (7.4 into the Helmholtz equation, we get eikx/prime/braceleftbigg∂2A ∂x/prime2+2ik∂A ∂x/prime−k2A+∂2A ∂y/prime2+k2A/bracerightbigg = 0 (7.5) Expecting that the characteristic scale LxofAalongx/primeis much longer than a wavelength, kLx/greatermuch1, we have 2ik∂A ∂x/prime/greatermuch∂2A ∂x/prime2 3.7. Diffraction of SH waves by a long crack 33 Hence we get as the first approximation the Schr´ odinger equation2 2ik∂A ∂x/prime+∂2A ∂y/prime2≈0 (7.7) In this transition zone where the remaining terms are of comparable importance, hence the length scales must be related by k x/prime∼1 y/prime2,implyingky/prime∼√ kx/prime Thus the transition zone is the interior of a parabola. Equation (7.7) is of the parabolic type. The boundary conditions are A(x,∞) = 0 (7.8) A(x,−∞)=Ao (7.9) The initial condition is A(0,y/prime)=  0,y/prime>0, A0,y/prime<0(7.10) he initial-boundary value for Ahas no intrinsic length scales except x/prime,y/primethemselves. Therefore the condition kLx/greatermuch1 meanskx/prime/greatermuch1 i.e., far away from the tip. This problem is somwhat analogous to the problem of one-dimensional heat diffusion across a boundary. A convenient way of solution is the method of similarity. Assume the solution A=Aof(γ) (7.11) where γ=−ky/prime √ πkx/prime(7.12) is the similarity variable. We find upon subsitution that fsatisfies the ordinary differ- ential equation f/prime/prime−iπγf/prime= 0 (7.13) 2In one-dimensional quantum mechanics the wave function in a potential-free field is governed by the Schr¨ odinger equation ih∂ψ ∂t+1 2M∂2ψ ∂x2= 0 (7.6) 3.7. Diffraction of SH waves by a long crack 34 subject to the boundary conditions that f→0,γ→− ∞ ;f→1,γ→∞. (7.14) Rewriting (7.13) as f/prime/prime f/prime=iπγ we get logf/prime=iπγ/2 + constant . One more integration gives f=C/integraldisplayγ −∞exp/parenleftbiggiπu2 2/parenrightbigg du Since/integraldisplay∞ 0exp/parenleftbiggiπu2 2/parenrightbigg du=eiπ/4 √ 2 we get C=e−iπ/4 √ 2 and f=A Ao=e−iπ/4 √ 2/integraldisplayγ −∞exp/parenleftbiggiπu2 2/parenrightbigg du=e−iπ/4 √ 2/braceleftbiggeiπ/4 √ 2+/integraldisplayγ 0exp/parenleftbiggiπu2 2/parenrightbigg du/bracerightbigg (7.15) Defining the cosine and sine Fresnel integrals by C(γ)=/integraldisplayγ 0cos/parenleftbiggπv2 2/parenrightbigg dv, S (γ)=/integraldisplayγ 0sin/parenleftbiggπv2 2/parenrightbigg dv (7.16) we can then write e−iπ/4 √ 2/braceleftbigg/bracketleftbigg1 2+C(γ)/bracketrightbigg +i/bracketleftbigg1 2+S(γ)/bracketrightbigg/bracerightbigg (7.17) In the complex plane the plot of C(γ)+iS(γ) vs.γis the famous Cornu’s spiral, shown in figure ( ??). The wave intensity is given by |A|2 A2 o=1 2/braceleftBigg/bracketleftbigg1 2+C(γ)/bracketrightbigg2 +/bracketleftbigg1 2+S(γ)/bracketrightbigg2/bracerightBigg (7.18) SinceC,S→0a sγ→− ∞ , the wave intensity diminshes to zero gradually into the shadow. However, C,S→1/2a sγ→∞ in an oscillatory manner. The wave 3.7. Diffraction of SH waves by a long crack 35 Figure 8: Amplitude contours of Hz. From Born and Wolf Optics according to the exact theory. Figure 9: Phase contours of Hz. From Born and Wolf Optics according to the exact theory. 3.7. Diffraction of SH waves by a long crack 36 Figure 10: Diffraction of a normally incident E-polarized plane wave Figure 11: Cornu’s spiral, a plot of the Fresnel integrals 3.7. Diffraction of SH waves by a long crack 37 intensity oscillates while approaching to unity asymptotically. In optics this shows up as alternately light and dark diffraction bands. In more complex propagation problems, the parabolic approximation can simplify the numerical task in that an elliptic boundary value problem involving an infinite domainis reduced to an initial boundary value problem. One can use Crank-Nicholson scheme to march in ”time”, i.e., x /prime. Homework Find by the parabolic approximation the transition solution along the edge of the reflection zone. 3.8. Rayleigh surface waves 38 8 Rayleigh surface waves Refs. Graff, Achenbach, Fung In a homogeneous elastic half plane, in addition to P, SV and SH waves, another wave which is trapped along the surface of a half plane can also be present. Because most of the action is near the surface, this surface wave is of special importance to seismic effects on the ground surface. Let us start from the governing equations again ∂2φ ∂x2+∂2φ ∂y2=1 c2 L∂2φ ∂t2, (8.1) ∂2Hz ∂x2+∂2Hz ∂y2=1 c2 T∂2Hz ∂t2(8.2) We now seek waves propagating along the x direction φ=/Rfractur/parenleftbig f(y)eiξx−iωt/parenrightbig ,H z=/Rfractur/parenleftbig h(y)eiξx−iωt/parenrightbig (8.3) Thenf(y),h(y) must satisfy d2f dy2+/parenleftbig ω2/c2 L−ξ2/parenrightbig f=0,d2h dy2+/parenleftbig ω2/c2 T−ξ2/parenrightbig h=0, (8.4) To have surface waves we insist that α=/radicalBig ξ2−ω2/c2 L, β=/radicalBig ξ2−ω2/c2 T (8.5) be real and postive. Keeping only the solutions which are bounded for y∼∞ ,w eg e t φ=Ae− αyei(ξx−ωt),H z=Be− βyei(ξx−ωt). (8.6) The expressions for the displacements and stresses can be found straightforwardly. ux=/parenleftBig iξAe− αy− βBe− βy/parenrightBig ei(ξx−ωt), (8.7) uy=−/parenleftBig αAe− αy+iξBe− βy/parenrightBig ei(ξx−ωt), (8.8) τxx=µ/braceleftBig/parenleftBig β2−ξ2−2 α2/parenrightBig Ae− αy−2i βξBe− βy/bracerightBig ei(ξx−ωt), (8.9) τyy=µ/braceleftBig/parenleftBig β2+ξ2/parenrightBig Ae− αy+2i βξBe− βy/bracerightBig ei(ξx−ωt), (8.10) τxy=µ/braceleftBig −2i αξAe− αy+/parenleftBig ξ2+ β2/parenrightBig Be− βy/bracerightBig ei(ξx−ωt)(8.11) 3.8. Rayleigh surface waves 39 On the free surface the traction-free conditions τyy=τxy= 0 require that /parenleftBig β2+ξ2/parenrightBig A+2i βξB =0, (8.12) −2i αξA +/parenleftBig β2+ξ2/parenrightBig B=0. (8.13) For nontrivial solutions of A,B the coefficient determinant must vanish, /parenleftBig β2+ξ2/parenrightBig2 −4 α βξ2=0, (8.14) or/bracketleftbigg 2ξ2−ω2 c2 T/bracketrightbigg2 −4ξ2/radicalBigg ξ2−ω2 c2 L/radicalBigg ξ2−ω2 c2 T= 0 (8.15) which is the dispersion relation between frequency ωand wavenumber ξ. From either (8.12) or (8.13) we get the amplitude ratio: A B=−2i βξ β2+ξ2= β2+ξ2 2i αξ, (8.16) In terms of the wave velocity c=ω/ξ, (8.15) becomes /parenleftbigg 2−c2 c2 T/parenrightbigg2 =4/parenleftbigg 1−c2 c2 L/parenrightbigg1 2/parenleftbigg 1−c2 c2 T/parenrightbigg1 2 . (8.17) or, upon squaring both sides, finally c2 c2 T/braceleftBigg/parenleftbiggc cT/parenrightbigg6 −8/parenleftbiggc cT/parenrightbigg4 +/parenleftbigg 24−16 κ2/parenrightbigg/parenleftbiggc cT/parenrightbigg2 −16/parenleftbigg 1−1 κ2/parenrightbigg/bracerightBigg =0. (8.18) where k=cL cT=/radicalBigg λ+2µ µ=/radicalbigg 2−2ν 1−2ν The first solution c=ω= 0 is at best a static problem. In fact α= β=ξandA=−iB, so thatux=uy≡0 which is of no interest. We need only consider the cubic equation for c2. Note that the roots of the cubic equation depend only on Poisson’s ratio, through κ2= 2(1−ν)/(1−2ν). There can be three real roots for corω, or one real root and two complex-conjugate roots. We rule out the latter because the complex roots imply either temporal damping or instability; neither of which is a propagating wave. When all three roots are real we must pick the one so that both αand βare real. We shall denote the speed of Rayleigh wave by cR. 3.8. Rayleigh surface waves 40 Figure 12: The velocity of Rayleigh surface waves cR. From Fung Foundations of Solid Mechanics . Forc= 0, the factor in curley brackets is {.}=−16/parenleftbigg 1−c2 T c2 L/parenrightbigg <0 Forc=cTthe same factor is equal to unity and hence positive. There must be a solution forcsuch that 0 <c<c T. Furthermore, we cannot have roots in the range c/cT>1. If so, β2=ξ2/parenleftbigg 1−c2 c2 T/parenrightbigg <0 which is not a surface wave. Thus the surface wave, if it exists, is slower than the shear wave. Numerical studies for the entire range of Poisson’s ratio (0 <ν< 0.5) have shown that there are one real and two complex conjugate roots if ν>0.263...and three real roots ifν<0.263.... But there is only one real root that gives the surface wave velocity cR. A graph of cRfor all values of Poisson’s ratio, due to Knopoff , is shown in Fig. 8. A curve-fitted expression for the Rayleigh wave velocity is cR/cT=( 0·87 + 1 ·12ν)/(1 +ν). (8.19) For rocks,λ=µandν=1 4, the roots are (c/cT)2=4,2+2/√ 3,2−2/√ 3. (8.20) 3.9. Moving load on the ground surface 41 The only acceptable root for Rayleigh wave speed cRis (cR/cT)2=( 2−2/√ 3)1 2=0·9194 (8.21) or cR=0.9588cT. (8.22) The particle displacement of a particle on the free surface is, from (8.7) and (8.8) ux=iA/parenleftBigg ξ− β2+ξ2 2ξ/parenrightBigg ei(ξx−ωt)(8.23) uy=A/parenleftBigg − α+ β2+ξ2 2ξ/parenrightBigg ei(ξx−ωt)(8.24) Note that a=A/bracketleftBigg ξ− β2+ξ2 2ξ/bracketrightBigg =A/bracketleftbigg ξ+k2 T 2ξ/bracketrightbigg >0 b=A/bracketleftBigg − α+ β2+ξ2 2ξ/bracketrightBigg =A/bracketleftbigg( α− β)2+k2 L 2 β/bracketrightbigg >0 hence ux=asin(ωt−ξx),u y=bcos(ωt−ξx) and u2 x a2+u2 y b2= 1 (8.25) The particle trajectory is an ellipse. In complex form we have ux a+iuy b= exp {i(ωt−ξx−π/2)} (8.26) Hence astincreases,a particle at ( x,0) traces the ellipse in the counter-clockwise direc- tion. See figure (8). 9 Elastic waves due to a load traveling on the ground surface Refs: Fung: Foundations of Solid Mechanics Cole and Huth: (1956, Elastic half space ; J Appl Mech 25, 433-436.) 3.9. Moving load on the ground surface 42 Figure 13: Displacement of particles on the ground surface in Rayleigh surface wave From Fung Foundations of Solid Mechanics . Mei, Si & Chen , (1985, Poro-elastic half space, Wave Motion , 7, 129-141.). In this section the yaxis is positive if pointing upwards. Let the traction on the ground surface be : τyy=−P(x+Ut).τ xy=0,ony= 0 (9.1) Let us make a (Galilean ) transformation to a coodinate system moving to the left at the speed of U, so that the load appears stationary, Then, by the chain rule, derivatives are changed accorindg ∂ ∂x→∂ ∂x,∂ ∂y→∂ ∂y,∂ ∂t→∂ ∂t+U∂ ∂x(9.2) In the moving coordinates, the wave equations are changed to ∂2Φ ∂x2+∂2Φ ∂y2=1 c2 L/parenleftbigg∂ ∂t+U∂ ∂x/parenrightbigg2 Φ, (9.3) ∂2H ∂x2+∂2H ∂y2=1 c2 T/parenleftbigg∂ ∂t+U∂ ∂x/parenrightbigg2 H (9.4) where we have abbreviated Hzsimply byH. In the steady state limit they become ∂2Φ ∂x2+∂2Φ ∂y2=U2 c2 L∂2Φ ∂x2(9.5) 3.9. Moving load on the ground surface 43 ∂2H ∂x2+∂2H ∂y2=U2 c2 T∂2H ∂x2(9.6) Introducing the Mach numbers: M1=U cL,M 2=U cT(9.7) then (9.5 ) amd (9.6) become /parenleftbig 1−M2 1/parenrightbig∂2Φ ∂x2+∂2Φ ∂y2= 0 (9.8) /parenleftbig 1−M2 2/parenrightbig∂2H ∂x2+∂2H ∂y2= 0 (9.9) BecausecL>c T, we must have M2>M 1. The stress components can be derived straigthtforwardly in terms of the potentials τxx=(λ+2µ)∇2Φ−2µ/parenleftbigg∂2Φ ∂y2−∂2H ∂x∂y/parenrightbigg =(λ+2µ)U2 c2 L∂2Φ ∂x2−2µ(M2 1−1)∂2Φ ∂x2+2µ∂2H ∂x∂y Using ther fact that λ+2µ)/µ=c2 L/c2 T, we further get τxx=µ/bracketleftbigg/parenleftbig M2 2−2M2 1+2/parenrightbig∂2Φ ∂x2+2∂2H ∂x∂y/bracketrightbigg (9.10) Similarly we find τyy=µ/bracketleftbigg (M2 2−2)∂2Φ ∂x2−2∂2H ∂x∂y/bracketrightbigg (9.11) τxy=µ/bracketleftbigg 2∂2Φ ∂x∂y+(M2 2−2)∂2H ∂x2/bracketrightbigg (9.12) We now examine the special pressure distribution, as shown in figure (9): p(x,0) =PoP(x)≡P0  e−α1x,x > 0, eα2x,x < 0(9.13) We jshall define a length scale L=/integraldisplay∞ −∞P(x)dx=1 α1+1 α2(9.14) Thus the traction boundary condtions (9.3) on the ground surface become (M2 2−2)∂2Φ ∂x2−2∂2H ∂x∂y=−Po µP(x),y= 0; (9.15) 2∂2Φ ∂x∂y+(M2 2−2)∂2H ∂x2=0,y=0. (9.16) Three cases will be distinguished: 3.9. Moving load on the ground surface 44 Figure 14: A moving pressure distibution on an elastic half space. Shown in a moving coordinate system, the pressure appears stationary. 3.9. Moving load on the ground surface 45 9.1 Supersonic: M2>M 1>1 Let us apply the exponential Fourier transform defined by f(λ)=/integraldisplay∞ −∞F(x)e−iλxdx, F(x)=1 2π/integraldisplay∞ −∞f(λ)eiλxdλ. (9.17) From the governing wave equations we get d2φ ∂y2+¯β2 1λ2φ=0,d2h ∂y2+¯β2 2λ2h=0 ;y<0. (9.18) where β2 j=M2 j−1,j =1,2. (9.19) The general solutions of the Fourier Transforms are φ=A(λ)ei¯β1y+B(λ)e−i¯β1y, h=C(λ)ei¯β2y+D(λ)ei¯β2y so that Φ(x,y)=1 2π/integraldisplay∞ −∞/bracketleftBig A(λ)eiλ(x+¯β1y)+B(λ)eiλ(x−¯β1y)/bracketrightBig dλ. (9.20) H(x,y)=1 2π/integraldisplay∞ −∞/bracketleftBig C(λ)eiλ(x+¯β2y)+D(λ)eiλ(x−¯β2y)/bracketrightBig dλ. (9.21) In order that waves below the ground surface trail behind the surface load, we discard the second term in each integral. Thus φ=A(λ)ei¯β1y,h =C(λ)ei¯β2y(9.22) Now the boundary conditions require 2iλdφ dy−λ2¯β2 2h= 0 (9.23) and −λ2¯β2 2φ−2iλdh dy=iPo µ/parenleftbigg1 λ−iα1−1 λ+iα2/parenrightbigg (9.24) 3.9. Moving load on the ground surface 46 Use has been made of the result /integraldisplay∞ −∞e−iλxP(x)dx=/integraldisplay0 −∞e−iλxeα2xdx+/integraldisplay∞ 0e−iλxe−α1xdx =1 i/parenleftbigg1 λ−iα1−1 λ+iα2/parenrightbigg (9.25) It follows that −2λ2¯β1A−¯β2 2λ2C= 0 (9.26) and −β2 2λ2A+2λ2¯β2C=iPo µ/parenleftbigg1 λ−iα1−1 λ+iα2/parenrightbigg (9.27) The last two equations can be solved to give A=iPo µ/parenleftbigg1 λ−iα1+1 λ+iα2/parenrightbigg−¯β2 2 λ2(¯β4 2+4¯β1¯β2)(9.28) C==iPo µ/parenleftbigg1 λ−iα1+1 λ+iα2/parenrightbigg2¯β1 λ2(¯β4 2+4¯β1¯β2)(9.29) The inverse transforms of φandhare:  Φ −H =iP0 2πµ k1 k2 /integraldisplay∞ −∞/parenleftbigg1 λ−iα1−1 λ+iα2/parenrightbigg · expiλ x+  β1 β2 y  dλ λ2(9.30) where k1=−(M2 2−2) (M2 1−2)2+4 β1 β2, k2=−2 β1 (M2 2−2)2+4 β1 β2. (9.31) Using (9.10), ( 9.11) and ( 9.12), we get the stress components. For example τxy Po=2¯β1k1 2πi/integraldisplay∞ −∞/parenleftbigg1 λ−iα1−1 λ+iα2/parenrightbigg eiλξ1dλ (9.32) −(M2 2−2)k2 2πi/integraldisplay∞ −∞/parenleftbigg1 λ−iα1−1 λ+iα2/parenrightbigg eiλξ2dλ (9.33) where ξ1=x+ β1y, ξ 2=x+ β2y (9.34) 3.9. Moving load on the ground surface 47 In view of (9.25), the inverse transform is immediate, τxy Po=2¯β1k1P(ξ1)−(M2 2−2)k2P(ξ2) (9.35) As a check, the shear stress on the ground surface y=0i s τxy Po=/bracketleftbig 2¯β1k1−(M2 2−2)k2/bracketrightbig P(x) = 0 (9.36) in view of (9.31). It can be shown that ˜τxx=τ0 xx P0=/parenleftbig M2 2−M2 1+2/parenrightbig k1P(ξ1)−2 β2k2P(ξ2), ˜τyy= (to be worked out) (9.37) Note that the disturbances in the half space indeed trail behind the surface pressure. The front of the P wave forms the first Mach wedge, followed by the front of the SV wave. Disturbances are concentrated essentially along the chracteristics x+ β1|y|= constant and x+ β2|y|=constant. Homework : Verify the above results by the method of characteristics. 9.2 Subsonic case, 1>M 2>M 1 In this case (9.8) and (9.9) are elliptic. Let β2 1=1−M2 1,β 2=1−M2 2, k3=−M2 2+2 (M2 2−2)2−4β1β2, k4=−2β1 (M2 2−2)2−4β1β2. (9.38) The formal solutions for Φ and Hare Φ=iP0k3 2πµ/integraldisplay∞ −∞/parenleftbigg1 λ−iα1−1 λ+iα2/parenrightbigg ×e|λ|β1yeiλxdλ λ2, −H=P0k4 2πµ/integraldisplay∞ −∞/parenleftbigg1 λ−iα1−1 λ+iα2/parenrightbigg ×e|λ|β2yeiλxsgnλ λ2dλ. (9.39) 3.9. Moving load on the ground surface 48 Figure 15: Stress variations in the ground under supersonic load on the surface. From Mei, Si & Chen, 1985. (In this article, the ground is poroelastic, pstands fo the pore pressure.) 3.9. Moving load on the ground surface 49 By using (9.10), (9.11) and (9.12), the stress components can be expressed as Fourier integrals, which can be evaluated in terms of the exponential integral defined by E1(z)=/integraldisplay∞ xe−τ τdτ =−γ−lnz−∞/summationdisplay n=1(−1)nzn nn!. (9.40) Let z1=x+iβ1y, z 2=x+iβ2y (9.41) and G(z)=e−α1zE1(−α1z)−eα2zE1(α2z). (9.42) Then the stress components are ˜τxx=τo xx P0=−/parenleftbig M2 2−2M2 1+2/parenrightbigk3 π/IfracturG(z1)−2β2k4 π/IfracturG(z2), ˜τyy=τo yy P0=−/parenleftbig M2 2−2/parenrightbigk3 π/IfracturG(z1)+2β2k4 π/IfracturG(z2) ˜τxy=τo xy P0=2β1k3 π/Rfractur[G(z2)−G(z1)] (9.43) Note that lim y↑0−E1(−α1z1)=  E1(α|x|)i fx<0, −Ei(αx)−iπifx>0, lim y↑0−E1(αz1)=  −Ei(α|x|)+πiifx<0, E1(αx)i f x>0(9.44) where Ei(x)=−PV/integraldisplay∞ xe−τ τdτ (9.45) with the integral being a principal value. From the definitions k3andk4become infinite when their denominators vanish. This occurs when the external load travels at the speed of the Rayleigh surface waveand indicates resonance. This is possible because Rayleigh wave speed is less that both c LandcT. The unbounded resonance need not be a threat in practise because the model of steady two-dimensional line load is an idealization not usually realized. 3.9. Moving load on the ground surface 50 Figure 16: Stress variations in the ground under subsonic load on the surface. From Mei, Si & Chen, 1985. 3.9. Moving load on the ground surface 51 9.3 Transonic case, M1>1>M 2 The scalar potentials are Φ=P0 2πµ/integraldisplay∞ −∞dλeiλxA(λ)e|λ|β1y, −H=P0 2πµ/integraldisplay∞ −∞dλeiλxB(λ)eiλ β2y(9.46) where A(λ)= −/parenleftbig M2 2−2/parenrightbig/parenleftbigg k5+i|λ| λk6/parenrightbiggi λ2/parenleftbigg1 λ−iα1−1 λ+iα2/parenrightbigg , B(α)=2iβ1/parenleftbigg|λ| λk5+ik6/parenrightbiggi λ2/parenleftbigg1 λ−iα1−1 λ+iα2/parenrightbigg , (9.47) k5=−(M2 2−2)2 (M2 2−2)4+1 6β2 1 β2 2, k6=−4β1 β2 (M2 2−2)4+1 6β2 1 β2 2. In terms of z1=x+iβ1yandξ2=x+ β|y|all the integrals in (3.18) can again be evaluated. The results involve the following functions: H(ξ)=eα2ξE1(α2ξ)+e−α1ξEi(α1ξ), H∗(ξ)=eα1|ξ|E1(α1|ξ|)+e−α2|ξ|Ei(α2|ξ|) ≡H(|ξ|). (9.48) The stresses are ˜τxx=τ0 xx P0=−/parenleftbig M2 2−2M2 1+2/parenrightbig/parenleftbig M2 2−2/parenrightbig1 π/Ifractur{(k5−ik6)G(z1)} −4β1 β2 π  k5H(ξ2)+πk6e−α1ξ2,ξ 2>0, −k5H∗(ξ2)+πk6e−α2|ξ2|,ξ 2<0,(9.49) ˜τyy=τ0 yy P0=−/parenleftbig M2 2−2/parenrightbig21 π/Ifractur{(k5−ik6)G(z1)} +4β1 β2 π  k5H(ξ2)+πk6e−α1ξ2,ξ 2>0, −k5H∗(ξ2)+πk6e−α2|ξ2|,ξ 2<0(9.50) 3.9. Moving load on the ground surface 52 Figure 17: Stress variations in the ground under transonic load on the surface.From Mei, Si & Chen, 1985. ˜τxy=τ0 xy P0=−/parenleftbig M2 2−2/parenrightbig2β1 π/Rfractur{(k5−ik6)G(z1)} −/parenleftbig M2 2−2/parenrightbig2β1 π  k5H(ξ2)+πk6e−α1ξ2,ξ 2>0, −k5H∗(ξ2)+πk6e−α2|ξ2|,ξ 2<0(9.51) Note Computations are for the following inputs: ν=1/2,µ=1 08N/m2,α 1=0.005m, α 2=0.1m, L = 210m A Partial wave expansion A useful result in wave theory is the expansion of the plane wave in a Fourier series of the polar angle θ. In polar coordinates the spatial factor of a plane wave of unit 3.9. Moving load on the ground surface 53 amplitude is eikx=eikrcosθ. Consider the following product of exponential functions ezt/2e−z/2t=/bracketleftBigg∞/summationdisplay n=01 n!/parenleftbiggzt 2/parenrightbiggn/bracketrightBigg/bracketleftBigg∞/summationdisplay n=01 n!/parenleftbigg−z 2t/parenrightbiggn/bracketrightBigg ∞/summationdisplay −∞tn/bracketleftbigg(z/2)n n!−(z/2)n+2 1!(n+ 1)!+(z/2)n+4 2!(n+ 2)!+···+(−1)r(z/2)n+2r r!(n+r)!+···/bracketrightbigg . The coefficient of tnis nothing but Jn(z), hence exp/bracketleftbiggz 2/parenleftbigg t−1 t/parenrightbigg/bracketrightbigg =∞/summationdisplay −∞tnJn(z). Now we set t=ieiθz=kr. The plane wave then becomes eikx=∞/summationdisplay N=−∞ein(θ+π/2Jn(z). Using the fact that J−n=(−1)nJn, we finally get eikx=eikrcosθ=∞/summationdisplay n=0/epsilon1ninJn(kr) cosnθ, (A.1) where/epsilon1nis the Jacobi symbol. The above result may be viewed as the Fourier expan- sion of the plane wave with Bessel functions being the expansion coefficients. In wave propagation theories, each term in the series represents a distinct angular variation and is called a partial wave . Using the orthogonality of cos nθ, we may evaluate the Fourier coefficient Jn(kr)=2 /epsilon1ninπ/integraldisplayπ 0eikrcosθcosnθdθ, (A.2) which is one of a host of integral representations of Bessel functions. 3.9. Moving load on the ground surface 54 B Approximate evaluation of an integral Consider the integral /integraldisplay2π 0dθ[1 + cos(θ−θo)]eikr(1−cos(θ−θo)) For largekrthe stationary phase points are found from ∂ ∂θ[1−cos(θ−θo)] = sin(θ−θo)=0 orθ=θo,θo+πwithin the range [0 ,2π]. Near the first stationary point the integrand is dominated by 2A(θo)eikt(θ−θo)2/2. When the limits are approximated by ( −∞,∞), the inegral can be evaluated to give A(θo)/integraldisplay∞ −∞eikrθ2/2dθ=/radicalbigg 2π kreiπ/4A(θo) Near the second stationary point the integral vanishes since 1+cos( θ−θo)= =1 −1=0 . Hence the result (6.28) follows.