Theory of Elasti Waves Gerhard Müller
Editors Mi hael Weber  GFZ & Universität Potsdam Georg Rümpker  Universität Fra...
96 downloads
720 Views
3MB Size
Report
This content was uploaded by our users and we assume good faith they have the permission to share this book. If you own the copyright to this book and it is wrongfully on our website, we offer a simple DMCA procedure to remove your content from our site. Start by pressing the button below!
Report copyright / DMCA form
Theory of Elasti Waves Gerhard Müller
Editors Mi hael Weber  GFZ & Universität Potsdam Georg Rümpker  Universität Frankfurt Dirk Gajewski  Universität Hamburg Germany tew_2007.ps + tew_2007.pdf
2
Prefa e When Gerhard Müller hose to leave us on 9 July 2002 be ause of his illness, we lost a tea her and olleague. Part of his lega y is several le ture notes whi h he had worked on for more then 20 years. These notes have be ome the ba kbone of tea hing seismi s and seismology at basi ally all German universities. When asked some years before his death if he had onsidered to translate "Theorie Elastis her Wellen" into English and publish it as a book, his answer was "I plan to do it when I am retired". We hope that our eort would have found his approval. We would like to thank R. and I. Coman (Universität Hamburg) for preparing a rst, German draft in LATEX of this s ript, A. Siebert (GFZ Potsdam) for her help in preparing the gures and our students for pointing out errors and asking questions. We would like to thank A. Priestley for proofreading the s ript and turning Deuts hlish into English and K. Priestley for his many omments. We thank the GFZ Potsdam and the Dublin Institute for Advan ed Studies for their support during a sabbati al of MW in Dublin, where most of this book was prepared. We would also like to thank the GFZ for ontinuing support in the preparation of this book.
M. Weber
G. Rümpker
D. Gajewski
Potsdam, Frankfurt, Hamburg January 2007
This le an be downloaded from
http://gfzpotsdam.de/mhw/tew/
tew _2007.ps(64MB) + tew _2007.pdf (3.5MB) 3
4
Prefa e of the German Le ture Notes This s ript is the revised and extended version of a manus ript whi h was used for several years in a 1 to 2semester le ture on the theory of elasti waves at the universities of Karlsruhe and Frankfurt. The aim of this manus ript is to give students with some ba kground in mathemati s and theoreti al physi s the basi knowledge of the theory of elasti waves, whi h is ne essary for the study of spe ial literature in monographs and s ienti journals. Sin e this is an introdu tory text, theory and methods are explained with simple models to keep the omputational omplexity and the formulae as simple as possible. This is why often liquid media instead of solid media are onsidered, and only horizontally polarised waves (SHwaves) are dis ussed, when shear waves in layered, solid media are onsidered. A third example is that the normal mode theory for point sour es is derived for an ideal wave guide with free or rigid boundaries.
These simpli ations o
asionally hide the dire t onne tion to
seismology. In my opinion, there is no other approa h if one aims at presenting theory and methods in detail and introdu ing at least some aspe t from the wide eld of seismology. After working through this s ript students should, I hope, be better prepared to read the advan ed text books of Pilant (1979), Aki and Ri hards (1980, 2000), BenMenahem and Sing (1981), Dahlen and Tromp (1998), Kennett (2002) and Chapman (2004), whi h treat models as realisti ally as possible. This manus ript has its emphasis in the wave seismi treatment of elasti body and surfa e waves in layered media. The understanding of the dynami properties of these two wave types, i.e., their amplitudes, frequen ies and impulse forms, are a basi prerequisiterequisite for the study of the stru ture of the Earth, may it be in the rust, the mantle or the ore, and for the study of pro esses in the earthquake sour e. Ray seismi s in inhomogeneous media and their relation with wave seismi s are dis ussed in more detail than in earlier versions of the s ript, but seismologi ally interesting topi s like eigenmodes of the Earth and extended sour es of elasti waves are still not treated, sin e they would ex eed the s ope of an introdu tory le ture. At several pla es of the manus ript, exer ises are in luded, the solution of these is an important part in understanding the material. One of the appendi es tries to over in ompa t form the basi s of the Lapla e and Fourier transform and of the delta fun tion, so that these topi s an be used in the main part of the s ript. I would like to thank Ingrid Hörn hen for the often tedious writing and orre ting of this manus ript.
Gerhard Müller
Contents 1 Literature
9
2 Foundations of elasti ity theory
11
2.1
Analysis of strain
. . . . . . . . . . . . . . . . . . . . . . . . . .
12
2.2
Analysis of stress . . . . . . . . . . . . . . . . . . . . . . . . . . .
21
2.3
Equilibrium onditions . . . . . . . . . . . . . . . . . . . . . . . .
23
2.4
Stressstrain relations
. . . . . . . . . . . . . . . . . . . . . . . .
26
2.5
Equation of motion, boundary and initial ... . . . . . . . . . . . .
31
2.6
Displa ement potentials and wave types
34
. . . . . . . . . . . . . .
3 Body waves
39
3.1
Plane body waves . . . . . . . . . . . . . . . . . . . . . . . . . . .
39
3.2
The initial value problem for plane waves
40
3.3
Simple boundary value problems for plane waves
3.4
Spheri al waves from explosion point sour es
3.5
Spheri al waves from single for e and dipole ...
3.6
. . . . . . . . . . . . . . . . . . . . . .
43
. . . . . . . . . . .
45
. . . . . . . . . .
49
3.5.1
Single for e point sour e . . . . . . . . . . . . . . . . . . .
49
3.5.2
Dipole point sour es . . . . . . . . . . . . . . . . . . . . .
56
Ree tion and refra tion of plane waves ... . . . . . . . . . . . . .
62
3.6.1
Plane waves with arbitrary propagation dire tion . . . . .
62
3.6.2
Basi equations . . . . . . . . . . . . . . . . . . . . . . . .
63
3.6.3
Ree tion and refra tion of SHwaves
. . . . . . . . . . .
66
3.6.4
Ree tion of Pwaves at a free surfa e . . . . . . . . . . .
76
5
6
CONTENTS
3.6.5 3.7
3.8
3.9
3.10
Ree tion and refra tion oe ients for layered media . .
Ree tivity method: Ree tion of ...
84
. . . . . . . . . . . . . . . .
92
3.7.1
Theory
. . . . . . . . . . . . . . . . . . . . . . . . . . . .
92
3.7.2
Ree tion and head waves . . . . . . . . . . . . . . . . . .
95
3.7.3
Complete seismograms . . . . . . . . . . . . . . . . . . . .
97
Exa t or generalised ray theory  GRT . . . . . . . . . . . . . . .
98
3.8.1
In ident ylindri al wave . . . . . . . . . . . . . . . . . . .
99
3.8.2
Wavefront approximation for
3.8.3
Ree tion and refra tion of the ylindri al wave . . . . . . 103
3.8.4
Dis ussion of ree ted wave types
UR
. . . . . . . . . . . . . . 102
. . . . . . . . . . . . . 110
Ray seismi s in ontinuous inhomogeneous media . . . . . . . . . 113 3.9.1
Fermat's prin iple and the ray equation
3.9.2
High frequen y approximation of the equation of motion . 119
3.9.3
Eikonal equation and seismi rays
3.9.4
Amplitudes in ray seismi approximation
. . . . . . . . . . 114
. . . . . . . . . . . . . 121 . . . . . . . . . 123
WKBJ method . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126 3.10.1 Harmoni ex itation and ree tion oe ient . . . . . . . 126 3.10.2 Impulsive ex itation and WKBJseismograms . . . . . . . 131
4 Surfa e waves 4.1
135
Free surfa e waves in layered media . . . . . . . . . . . . . . . . . 135 4.1.1
Basi equations . . . . . . . . . . . . . . . . . . . . . . . . 135
4.1.2
Rayleigh waves at the surfa e of an homogeneous halfspa e138
4.1.3
Love waves at the surfa e of a layered halfspa e
4.1.4
Determination of the phase velo ity of surfa e waves from
. . . . . 142
observations . . . . . . . . . . . . . . . . . . . . . . . . . . 150 4.1.5
The group velo ity . . . . . . . . . . . . . . . . . . . . . . 153
4.1.6
Des ription of surfa e waves by onstru tive interferen e of body waves . . . . . . . . . . . . . . . . . . . . . . . . . 157
4.2
Surfa e waves from point sour es . . . . . . . . . . . . . . . . . . 160 4.2.1
Ideal wave guide for harmoni ex itation . . . . . . . . . . 160
4.2.2
The modal seismogram of the ideal wave guide
4.2.3
Computation of modal seismograms with the method of stationary phase
4.2.4
. . . . . . 167
. . . . . . . . . . . . . . . . . . . . . . . 170
Ray representation of the eld in an ideal wave guide
. . 173
7
CONTENTS
A Lapla e transform and delta fun tion A.1
A.2
179
Introdu tion to the Lapla e transform
. . . . . . . . . . . . . . . 179
A.1.1
Literature . . . . . . . . . . . . . . . . . . . . . . . . . . . 179
A.1.2
Denition of the Lapla e transform . . . . . . . . . . . . . 179
A.1.3
Assumptions on F(t) . . . . . . . . . . . . . . . . . . . . . 180
A.1.4
Examples
A.1.5
Properties of the Lapla e transform
A.1.6
Ba ktransform . . . . . . . . . . . . . . . . . . . . . . . . 183
A.1.7
Relation with the Fourier transform
. . . . . . . . . . . . . . . . . . . . . . . . . . . 180 . . . . . . . . . . . . 181
. . . . . . . . . . . . 184
Appli ation of the Lapla e transform . . . . . . . . . . . . . . . . 184 A.2.1
Linear ordinary dierential equations with onstant oef ients . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 184
A.2.2 A.3
Partial dierential equations
The delta fun tion
δ(t)
. . . . . . . . . . . . . . . . 192
. . . . . . . . . . . . . . . . . . . . . . . . 194
δ(t)
. . . . . . . . . . . . . . . . . . . . . 194
A.3.1
Introdu tion of
A.3.2
Properties of
A.3.3
Appli ation of
A.3.4
Duhamel's law and linear systems
A.3.5
Pra ti al approa h for the onsideration of nonzero initial
δ(t)
. . . . . . . . . . . . . . . . . . . . . . . 197
δ(t)
. . . . . . . . . . . . . . . . . . . . . . 200 . . . . . . . . . . . . . 202
values of the perturbation fun tion F(t) of a linear problem205
B Hilbert transform
209
B.1
The Hilbert transform pair
. . . . . . . . . . . . . . . . . . . . . 209
B.2
The Hilbert transform as a lter
. . . . . . . . . . . . . . . . . . 210
C Bessel fun tions
215
D The Sommerfeld integral
219
E The omputation of modal seismograms
221
E.1
Numeri al al ulations . . . . . . . . . . . . . . . . . . . . . . . . 221
E.2
Method of stationary phase
E.3
Airy phases . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 224
Index
226
. . . . . . . . . . . . . . . . . . . . . 222
8
CONTENTS
Chapter 1
Literature The following list ontains only books, that treat the propagation of elasti waves and a few text books on ontinuum me hani s. Arti les in journals are mentioned if ne essary. Their number is kept to a minimum.
A henba h, J.D.: Wave propagation in elasti solid, NorthHolland Publ. Comp., Amsterdam, 1973 Aki, K. and P.G. Ri hards: Quantitative seismology  theory and methods (2 volumes), Freeman and Co., San Fran is o, 1980 and 2002 BenMenahem, A. and S.J. Singh: Seismi waves and sour es, Springer, Heidelberg, 1981 Bleistein, N.:
Mathemati al methods for wave phenomena, A ademi Press,
New York, 1984 Brekhovskikh, L.M.: Waves in layered media, A ademi Press, New York, 1960 Brekhovskikh, L., and Gon harov, V.: Me hani s of ontinua and wave dynami s, SpringerVerlag, Berlin, 1985 Budden, K.G.: The waveguide mode theory of wave propagation, Logos Press, London, 1961 Bullen, K.E., and Bolt, B.A.:
An introdu tion to the theory of seismology,
Cambridge University Press, Cambridge, 1985 Cagniard, L.: Ree tion and refra tion of progressive waves, M GrawHill Book Comp., New York, 1962
ervéný, V., I.A. Molotov and I. P² en ík:
Ray method in seismology, Univerzita
Karlova, Prague, 1977 Chapman, Ch.: Fundamentals of seismi wave propagation, Cambridge University Press, 2004 9
10
CHAPTER 1.
LITERATURE
Dahlen, F. A. and J. Tromp: Theoreti al global seismology, Prin eton University, Prin eton, New Jersey, 1998 Ewing, M., W.S. Jardetzky and F. Press:
Elasti waves in layered media,
M GrawHill Book Comp., New York, 1957 Fung, Y.C.: Foundations of solid me hani s, Prenti eHall, Englewood Clis, N.Y., 1965 Grant, F.S. and G.F. West: Interpretation theory in applied geophysi s, M GrawHill Book Comp., New York, 1965 Hudson, J.A.: The ex itation and propagation of elasti waves, Cambridge University Press, Cambridge, 1980 Kennett, B.L.N.:
The seismi wave eld (2 volumes), Cambridge University
Press, Cambridge, 2002 Landau, L.D. and E.M. Lifs hitz: Elastizitätstheorie, Akademie Verlag, Berlin, 1977 Love, A.E.H.: A treatise on the mathemati al theory of elasti ity, 4th edition, Dover Publi ations, New York, 1944 Pilant, W.L.: Elasti waves in the Earth, Elsevier, Amsterdam, 1979 Riley, K.F., M.P. Hobson and S.J. Ben e: Mathemati al methods for physi s and engineering, A omprehensive guide, Cambridge University Press, Cambridge, 2nd edition, 2002 Sommerfeld, A.: Me hanik der deformierbaren Medien, Akad. Verlagsgesells haft, Leipzig, 1964 Tolstoy, I.: Wave propagation, M GrawHill Book Comp., New York, 1973 Tolstoy, I. and C.S. Clay: O ean a ousti stheory and experiment in underwater sound, M GrawHill Book Comp., New York, 1966 White, J.E.: Seismi wavesradiation, transmission and attenuation, M GrawHill Book Comp., New York, 1965
Chapter 2
Foundations of elasti ity theory Comments In this hapter symboli and index notation is used, i.e., a ve tor (symboli
− → f ) is also written as fi ( omponents f1 , f2 , . . . , fn ), the lo ation ve tor → (symboli − x ) as xi ( omponents x1, x2, x3 ), and a matrix (symboli a) as aij
notation (i
= line index = 1, 2, . . . , m, j = row index= 1, 2, . . . , n). aij with the ve tor fj is the ve tor gi =
n X
aij fj
The produ t of matrix
(i = 1, 2, . . . , m).
j=1
summation onvention = SC )
A short notation for this is (
gi = aij fj . In the following text, if a produ t on the right o
urs in whi h there is a repeated index, this index takes all values from 1, 2, ..., n (usually n = 3) and all produ ts are summed. If the symboli notation is simpler, e.g., for the ross produ t of two ve tors or for divergen e or rotation, the symboli notation is used.
11
12
CHAPTER 2.
FOUNDATIONS OF ELASTICITY THEORY
2.1 Analysis of strain 2.1.1 Components of the displa ement ve tor Consider a body that is deformed by an external for e. Before deformation, the
xi and the innitesimal lose point Q has the xi + yi . The omponents of yi are assumed to be independent variables; this is why dxi was not used. After deformation, P has been displa ed by the displa ement ve tor ui to P', and Q has been displa ed to Q' by the point
P
has the lo ation ve tor
lo ation ve tor
ve tor (expansion up to linear terms)
zi = ui + dui = ui +
∂ui yj ∂xj
(SC).
Fig. 2.1: Neighbourhood of P and Q before and after deformation.
Ve tor
P
zi
des ribes (for variable
Q
in the neighbourhood of
P ) the hanges near
due to the deformation. In general, these hanges onsist of: a translation,
a rotation of the
whole
region around an axis through
P
and the a tual de
formation, whi h hanges the length of lines (rotation and deformation will be dis ussed later in more detail)
zi = ui + dui =
ui translation
ǫij
=
ǫij ξij
= =
+
ǫij yj deformation
+
ξij yj rotation
1 ∂ui 1 ∂ui ∂uj ∂uj , ξij := + − 2 ∂xj ∂xi 2 ∂xj ∂xj ǫji −ξji (⇒ ξ11 = ξ22 = ξ33 = 0) .
(2.1) (2.2) (2.3)
2.1.
The matri es and
13
ANALYSIS OF STRAIN
ξij
ǫij
and
ξij
are
tensors of 2nd degree. ǫij is symmetri due to (2.2) ǫij is alled deformation tensor and ξij
is antisymmetri due to (2.3).
is alled
rotation tensor.
2.1.2 Tensors of 2nd degree A tensor of 2nd degree, transforms ve tor
yi
tij ,
transforms a ve tor into another ve tor (e.g.,
into the deformation part of
dui ;
ǫij
another example of this is
the inertial tensor transforms the ve tor of the angular velo ity into the rotation impulse ve tor (rotation of a rigid body)). If the oordinate system is rotated, the tensor omponents have to be transformed as follows to yield the original ve tor
t,kl amn
= =
aik ajl tij cos γmn
(SC twi e) (see sket h).
(2.4)
t,kl = Tensor omponent in the rotated oordinate system (dashed line in sket h).
Fig. 2.2: Coordinate system of tensors of 2nd degree.
For a ertain orientation of the rotated system, the nondiagonal elements
t′12 , t′13 , t′21 ,
... vanish. These oordinate axis are alled
main axes
of the tensor,
and the tensor is in diagonal form. In the diagonal form, many physi al relations be ome simpler. Certain ombinations of tensor omponents are independent of the oordinate system of the tensor. These are the three
T3
invariants (T1 , T2 ,
are the diagonal elements of the tensor in diagonal form). More on tensors
an be found in, e.g., Riley, Hobson and Ben e.
14
CHAPTER 2.
I1 = I2 = I3 =
t11 t12 t13 t21 t22 t23 t31 t32 t33 t11 + t22 + t33
FOUNDATIONS OF ELASTICITY THEORY
= T1 T2 T3
(determinant)
= T1 + T2 + T3
(tra e)
t11 t22 + t22 t33 + t33 t11 − t12 t21 − t23 t32 − t31 t13
= T1 T2 + T2 T3 + T3 T1 .
2.1.3 Rotation omponent of displa ement The rotation omponent of displa ement follows from
0 −ξ12 −ξ13
ξ12 y2 + ξ13 y3 ξ13 y1 − → → ξ23 y2 = −ξ12 y1 + ξ23 y3 = ξ × − y −ξ13 y1 − ξ23 y2 0 y3
ξ12 0 −ξ23
with
1 − → → u. ξ = (−ξ23 , ξ13 , −ξ12 ) = ∇ × − 2
− → − ξ ×→ y des ribes an innitesimal rotation of the region of P around − → an axis through P with the dire tion of ξ . The rotation angle has the abso− → → − → lute value ξ and is independent of − y (show). A prerequisite is that ξ is Ve tor
innitesimal. A su ient ondition for this is that
∂ui ∂xj ≪ 1
for all
i
and
j.
(2.5)
2.1.4 Deformation omponent of displa ement After separating out the rotation term, only the deformation term is of interest sin e it des ribes the for es whi h a t in a body. The deformation is des ribed
ompletely by the six omponents
ǫij
whi h are, in general, dierent.
These
dimensionless omponents will now be interpreted physi ally. The starting point is
dui = ǫij yj ,
i.e., we assume no rotation.
a) During this transformation, a line remains a line, a plane remains a plane, a sphere be omes an ellipsoid and parallel lines remain parallel.
Deformation omponents ǫ11 , ǫ22 , ǫ33 Coordinate origin at P and spe ial sele tion of Q : y1 6= 0, y2 = y3 = 0.
b)
2.1.
15
ANALYSIS OF STRAIN
2 Q’
du 2
y1
P
Q
du 1
1
Fig. 2.3: Sket h for deformation omponents.
du1
=
ǫ11 y1
du2 du3
= =
ǫ21 y1 ǫ31 y1 = 0
(assumption :
ǫ
31
not the relative hange
du1 y1 is the relative hange in length in dire tion 1 ( in length of , see also d). Stret hing o
urs if ǫ11 >
ǫ11 =
ǫ11 < 0.
PQ
Similarly,
ǫ22
and
ǫ33
0
and shortening if
are the relative length hanges in dire tion 2 and
3.
)
= 0).
Shear omponents ǫ12 , ǫ13 , ǫ23
Fig. 2.4: Sket h for shear omponents.
Q1 → Q′1 : du2 Q2 → Q′2 : du1
= =
ǫ21 y1 = ǫ12 y1 ǫ12 y2
16
CHAPTER 2.
FOUNDATIONS OF ELASTICITY THEORY
du2 y1 du1 tan β ≃ β ≃ y2
tan α ≃ α ≃
This means, right angle at
ǫ12
P
=
ǫ12
=
ǫ12 .
is the angle around whi h the 1 or 2axis is rotated. is redu ed by
after deformation (sin e
2ǫ12 .
ǫ13 , ǫ23
or
The
If the parallelogram is not in the 12 plane
ǫ33
is nonzero), these statements hold for
the verti al proje tion in this plane. Similar results hold for d)
ǫ13
ǫ23 .
and
Length hanges of distan e P Q
P Q.
Fig. 2.5: Sket h for length hanges of distan e
PQ =
P ′ Q′
=
=
l0 =
l= (
(
(
3 X
3 X i=1
3 X
yi2
)1/2 2
(yi + ǫij yj )
i=1
yi2
+ 2ǫij yi yj +
)1/2
3 X i=1
i=1
= 2
(ǫij yj )
)1/2
.
The 1st, 2nd and 3rd term require SC on e, twi e and three times, respe tively. The 3rd term ontains only squares of
innitesimal strain theory
ǫij
and an, within the framework of
treated here, be negle ted relative to the 2nd term
(the prerequisite for this is (25))
2.1.
17
ANALYSIS OF STRAIN
12 2 1 l = l0 1 + 2 ǫij yi yj = l0 + ǫij yi yj . l0 l0 Relative length hanges
yi yj l − l0 = ǫij 2 = ǫij ni nj (SC l0 l0 ni = Approa hes for
yi = unit l0
twi e; quadrati form in
ve tor in dire tion of
nite strain theory
nk )
yi .
exist (see, e.g., Bullen and Bolt). Su h a
theory has to be developed from the very start. Then, for example, the simple separation of the rotation term in the displa ement ve tor, whi h is possible for innitesimal deformations, is no longer possible. The deformation tensor
ǫij
also be omes more ompli ated. e)
Volume hanges ( ubi dilatation)
V ontaining point P surfa e S. After deformation, for whi h we assume without loss of generality that P remains in its position, volume V is hanged by ∆V .
We onsider a nite (not innitesimal) volume with
Fig. 2.6: Sket h for volume hanges.
∆V =
Z
un df.
S
Transformation of this surfa e integral with Gauss' law gives
∆V =
Z
V
→ ∇·− u dV,
18
CHAPTER 2.
FOUNDATIONS OF ELASTICITY THEORY
and the relative volume hange an be written as
∆V 1 = V V Into the limit
V →0
V
(shrinking to point
lim
V →0 This limit is alled
Z
→ ∇·− u dV.
(2.6)
P ), this be omes
∆V = Θ. V
ubi dilatation.
From (2.6) with (2.1), it follows that
∂u1 ∂u2 ∂u3 → + + = ǫ11 +ǫ22 +ǫ33 Θ = ∇·− u := ∂x1 ∂x2 ∂x3 For
Θ>0
the volume in reases, for
Θ x1
x2 − x1 x2 − x1 x2 x1 =F t− = F1 t − . u(x2 , t) = F t − − c c c c
t at distan e x2 the same ee ts o
ur as at distan e t − (x2 − x1 )/c. This orresponds to a wave whi h has x2 in the time (x2 − x1 )/c. The propagation velo ity is,
This means that for time
x1
at the earlier time
travelled from therefore,
.
x1
The
to
wavefronts
of this wave, i.e., the surfa es between perturbed
and unperturbed regions, are the planes
x=
onst. Therefore, these are
plane
waves. If G(x) in (3.2) or G(t) in (3.3) are not zero, two plane waves propagate in opposite dire tions. In the ase of
u = ux , we have a longitudinal wave (polarisation in u = uy , we have a transverse wave
of propagation); in ase of
perpendi ular to the dire tion of propagation). Harmoni waves an be represented as
the dire tion (polarisation
h x i u(x, t) = A exp iω(t − ) = A exp [i(ωt − kx)] c
A=
Amplitude (real or omplex), ω = angular frequen y ν = ω/2π = freT = 1/ν = period, k = ω/c= wavenumber and Λ = 2π/k = wave length. Between c, Λ and ν the wellknown relation c = Λν holds. The use of the omwith
quen y,
plex exponential fun tion in the des ription of plane harmoni waves is more
onvenient than the use of the real sine and osine fun tions. In the following, only the exponential fun tion will be used.
3.2 The initial value problem for plane waves We look for the solution of the onedimensional wave equation (3.1) whi h satises the initial onditions
3.2.
THE INITIAL VALUE PROBLEM FOR PLANE WAVES
u(x, 0) = f (x)
41
for displa ement and
∂u (x, 0) = g(x) ∂t
for parti le velo ity.
This is an initial value problem of a linear
ordinary
dierential equation, e.g.,
the problem to determine the movement of a pendulum, if initial displa ement
t=0, it also holds that
and initial velo ity are given. We start from (3.2). For
F (x) + G(x)
= f (x)
(3.4)
′
= g(x).
(3.5)
′
−cF (x) + cG (x) Integrating (3.5) with respe t to
x, gives
F (x) − G(x) = −
1 c
Z
x
g(ξ)dξ.
(3.6)
−∞
From the addition of (3.4) and (3.6), it follows that
1 F (x) = 2
Z 1 x f (x) − g(ξ)dξ , c −∞
and from the subtra tion of these two equations that
G(x) =
1 2
Z 1 x f (x) + g(ξ)dξ . c −∞
From this, it follows that
1 1 u(x, t) = {f (x − ct) + f (x + ct)} + 2 2c
Z
x+ct
g(ξ)dξ. x−ct
This solution satises the wave equation and the initial onditions ( he k). We will dis uss
two spe ial ases.
3.2.1 Case 1 g(x) = 0,
i.e., the initial velo ity is zero. Then
u(x, t) =
1 {f (x − ct) + f (x + ct)} . 2
42
CHAPTER 3.
BODY WAVES
Two snap shots (Fig. 3.1) for t=0 and for t>0, illustrate this result.
u(x,0) f(x)
t=0
x
ct
ct u(x,t)
f(xct)
f(x+ct)
t>0
x
Fig. 3.1: Snap shots of two plane waves.
Two plane waves propagate from the point of ex itation in both dire tions with the velo ity
t=0.
.
A pra ti al example is a stret hed rope with the form
f (x)
for
3.2.2 Case 2 f (x) = 0, i.e., the initial displa ement is zero. Furthermore, we assume g(x) = V0 δ(x). δ(x) is Dira 's delta fun tion (see appendix A.3). g(x) orresponds to an impulse at x = 0. V0 has the dimension of velo ity times length. Then
u(x, t) =
The sket h (Fig.
V0 2c
Z
x+ct
δ(ξ)dξ.
x−ct
3.2) shows the value of the integrand and the integration
interval for a xed point in time
t>0
and for a lo ation
x > 0.
3.3.
SIMPLE BOUNDARY VALUE PROBLEMS FOR PLANE WAVES
43
Fig. 3.2: Value of the integrand and the integration interval of u(x,t).
Only when the integration interval in ludes the point
ξ = 0,
does the integral
be ome nonzero. Then it always has the value of 1
u(x, t) =
H(t)
is the
V0 2c H(t
− xc )
for x>0
V0 2c H(t
+ xc )
for x r1 .
We start from the equation of motion (2.24) with
− → f = 0.
This is how the
problem is solved in appendix A (appendix A.2.2) using the Lapla e transform.
46
CHAPTER 3.
BODY WAVES
Here, more simply, the displa ement potential from se tion 2.6 will be used. In this ase, the shear potential is zero, sin e a radially symmetri radial ve tor
an be derived solely from the ompressional potential
Φ(r, t) =
Z
r
U (r′ , t)dr′ .
In spheri al oordinates (r, ϑ, λ), it holds that
∇Φ = For
Φ,
1 ∂Φ ∂Φ 1 ∂Φ ∂r , r ∂ϑ , r sin ϑ ∂λ
the wave equation with
∇2 Φ =
ϕ=0
= (U (r, t), 0, 0).
an by written a
ording to (2.31) as
∂ 2 Φ 2 ∂Φ 1 ∂2Φ 1 ∂ 2 (rΦ) + = 2 2 = 2 2 ∂r r ∂r r ∂r α ∂t ∂ 2 (rΦ) 1 ∂ 2 (rΦ) = 2 . 2 ∂r α ∂t2
(3.7)
In the ase of radial symmetry, the wave equation an be redu ed to the form of a onedimensional wave equation for Cartesian oordinates for the fun tion
rΦ, 1 ∂2u ∂2u = , ∂x2 α2 ∂t2 . The most general solution for (3.7) is ,therefore,
Φ(r, t) =
r r o 1n F (t − ) + G(t + ) . r α α
This des ribes the superposition of two ompressional waves, one propagating outward from the point sour e and the other propagating inwards towards the point sour e. In realisti problems, the se ond term is always zero and
Φ(r, t) = Fun tion potential.
F (t)
is often alled the
1 r . F t− r α
ex itation fun tion
The wavefronts are the spheres
(3.8)
or redu ed displa ement
r = onst.
The potential as a
fun tion of time has the same form everywhere, and the amplitudes de rease with distan e as
1/r.
The radial displa ement of the spheri al wave onsists of
two ontributions with dierent dependen e on
r
3.4.
SPHERICAL WAVES FROM EXPLOSION POINT SOURCES
U (r, t) =
∂Φ 1 ′ r r 1 − . = − 2F t− F t− ∂r r α rα α
47
(3.9)
These two terms, therefore, hange their form with in reasing distan e. Generally, this holds for the displa ement of waves from a point sour e. The rst term in (3.9) is alled
neareld term sin e it dominates for su iently small r. fareld term and des ribes with su ient a
ura y the
The se ond term is the
displa ement for distan es from the point sour e whi h are larger then several F (t) = eiωt ). That means there the displa ement
wave lengths (show this for redu es proportional to
1/r.
From the boundary ondition
−
r = r1 ,
it follows that
1 ′ 1 r1 r1 − 2F t− = U1 (t). F t− r1 α α r1 α
We hoose the origin time so that
U1 (t) only begins to be nonzero for t = r1 /α.
It, therefore, appears as if the wave starts at time t = 0 at the point sour e r (r = 0). If U 1 t − 1 = U1 (t), it holds that U 1 (t) is already nonzero for t > 0. α r With τ = t − α1 , it follows that
1 1 ′ F (τ ) + 2 F (τ ) = −U 1 (τ ). r1 α r1
(3.10)
The solution of (3.10) an be found with the Lapla e transform (see se tion A.2.1.1 of appendix A). Sin e the riterion (A.16) of appendix A is satised for all physi ally realisti displa ements
U 1 (τ ),
1 r1
F (+0) of F (τ ) is zero. Therefore, F (τ ) ←→ f (s) and U 1 (τ ) ←→ u1 (s) gives
the initial value
transformation of (3.10)with
s 1 + α r1
f (s) = −u1 (s)
f (s) = −r1 α With
s+
α r1
−1
1 u1 (s) s + rα1
.
(3.11)
α
←→ e− r1 τ
(see appendix A, se tion A.1.4), and using on
volution (see appendix A, equation A.7), the inverse transformation of (3.11) reads as
F (τ ) = −r1 α From this, it follows
Z
0
τ
U 1 (ϑ)e
− rα (τ −ϑ) 1
dϑ.
48
CHAPTER 3.
F ′ (τ ) = −r1 αU 1 (τ ) + α2 The radial displa ement for
U (r, t) = with the
r1 r
r > r1
Z
τ
BODY WAVES
α
U 1 (ϑ)e− r1 (τ −ϑ) dϑ.
0
then an be written using (3.9) as
Z τ 1 1 − α (τ −ϑ) − dϑ U 1 (τ ) + α U 1 (ϑ)e r1 r r1 0
retarded time τ = t − αr .
(3.12)
This solves the boundary problem.
Appli ations 1.
U1 (t) = U 0 δ t −
r1 α
The dimension of
U0
, i.e., U 1 (t) = U 0 δ(t). is time times length. Equation (3.12) is valid also in
this ase (see appendix A)
1 r1 1 − rα τ 1 U (r, t) = U 0 δ(τ ) + α e − H(τ ) . r r r1
τ=0 (t=r/α)
τ
Fig. 3.3: U(r,t) as a fun tion of time.
2.
U1 (t) = U0 H t −
r1 α
The dimension of
U0
U (r, t)
= = =
, i.e. U 1 (t) = U0 H(t). is length. From (3.12), it follows
h r α iϑ=τ α r1 1 1 1 r ϑ H(τ ) e − r1 τ U0 H(τ ) + α − e 1 r r r1 α ϑ=0 n r o r1 1 − rα τ U0 H(τ ) 1 + −1 1−e 1 r r n r1 r1 − rα τ o r1 U0 H(τ ) + 1− e 1 . r r r
3.5.
SPHERICAL WAVES FROM SINGLE FORCE AND DIPOLE ...
49
r1 r U0
r12 U0 r2 τ=0
τ
Fig. 3.4: U(r,t) as a fun tion of time.
Exer ise 3.4 r1 a ts in a spheri al avity with radius r1 . α What is the dierential equation for the ex itation fun tion F (t) (analogue to Pressure
P (t)
P (t) = 0 for t
1
3.5.
SPHERICAL WAVES FROM SINGLE FORCE AND DIPOLE ...
53
it follows that
Z r z r′ ′ dr′ r K t − 4πρα2 r3 0 α Z αr z K(t − τ )τ dτ. 4πρr3 0
Φ(x, y, z, t) = = The wave equations for
Ψx
and
Ψy
are solved in a similar fashion. Therefore,
it is possible to write the potentials of the single for e point sour e as
Φ(x, y, z, t)
=
Ψx (x, y, z, t) = Ψy (x, y, z, t)
=
Ψz (x, y, z, t)
=
with r2
=
r α
r R y β − 4πρr3 0 K(t − τ )τ dτ r R β x K(t − τ )τ dτ 4πρr 3 0 0 2 2 2 x +y +z . z 4πρr 3
R
0
K(t − τ )τ dτ
Before we derive the displa ements, we hange to spheri al oordinates relative to the single for e point sour e
z
K(t)
P υ r
x
Fig. 3.7: Spheri al oordinates
λ
(r, ϑ, λ).
y
(3.16)
(r, ϑ, λ)
54
CHAPTER 3.
x
= r sin ϑ cos λ
y
= r sin ϑ sin λ
z
= r cos ϑ.
In spheri al oordinates, the shear potential has no and for the
λ omponent
BODY WAVES
r and ϑ omponent (show),
it holds that
Ψλ
=
−Ψx sin λ + Ψy cos λ.
y
Ψy
Ψx λ x Fig. 3.8: xyplane of Fig. 3.7.
This gives
Φ(r, ϑ, t) Ψλ (r, ϑ, t)
=
cos ϑ 4πρr 2
=
sin ϑ 4πρr 2
This equation does not depend on
− → − → u = ∇Φ + ∇ × Ψ
λ.
R
r α
0
R
The
r β
0
K(t − τ )τ dτ
K(t − τ )τ dτ.
(3.17)
displa ement ve tor
an be written in spheri al oordinates (show) as
ur
=
∂Φ ∂r
uϑ
=
1 ∂Φ r ∂ϑ
uλ
=
0.
+
1 ∂ r sin ϑ ∂ϑ (sin ϑΨλ )
−
1 ∂ r ∂r (rΨλ )
(3.18)
3.5.
55
SPHERICAL WAVES FROM SINGLE FORCE AND DIPOLE ...
This shows that the
P wave following from Φ is not purely longitudinal, but it (in uϑ ). Similarly, the S wave following from
ontains a transverse omponent
Ψλ is not purely transverse sin e ur ontains a shear omponent. The rst term in uϑ and the se ond in ur are near eld terms ( ompare exer ise 3.5). Here we ompute only the fareld terms of ur and uϑ (only dierentiation of the integrals in(3.17))
ur uϑ
≃
cos ϑ 4πρα2 r K
≃
− sin ϑ 4πρβ 2 r K
t−
t−
r α r β
(longitudinal P − wave) (transversal
The fareld displa ements have, therefore, the form of the for e with
1/r.
and the
(3.19)
S − wave).
K(t) de reasing
The single for e point sour e has dire tionally dependent radiation,
fareld radiation hara teristi s
are shown in Fig. 3.9.
Fig. 3.9: Far eld radiation hara teristi s of single for e point sour e.
The radiation hara teristi s (P and Swaves) are ea h two ir les. Those for 2 2 the waves have a radius whi h is α /β larger then those of the Pwaves. If the
S radiation angle ϑ is varied for xed r, the displa ements ur are proportional
to the distan e
OP 2 .
OP 1 ,
and the displa ements
The sign of the displa ement
radiation ir le to the se ond. from that shown in Fig.
ur
uϑ
are proportional to the distan e
hanges in transition from the rst
P
The full 3D radiation hara teristi s follows
3.9 by rotation around the dire tion of the for e.
S wave is radiated in P wave is radiated (but
Within the framework of the fareld equations (3.19), no the dire tion of the for e, and perpendi ular to it, no
ompare exer ise 3.5).
The pra ti al use of the single for e point sour e, a ting perpendi ular on the free surfa e, is that it is a good model for the ee t of a by
vibroseis
and often also for
explosions
drop weight, ex itation
detonated lose to the surfa e.
A
omplete solution requires the onsideration of the ee ts of the free surfa e, but that is signi antly more ompli ated. Furthermore, the dieren es to the
56
CHAPTER 3.
fullspa e models for all
P waves
and for
S waves,
BODY WAVES
for radiation angles smaller
then 30 to 40 degrees, are small, respe tively.
Exer ise 3.5 Compute the omplete displa ement (3.18) using (3.17) and examine in parti 0 ular, the dire tions ϑ = 0 and ϑ = 90 . Whi h polarisation does the displa ement ve tor have, and at whi h times are arrivals to be expe ted? Compute for
K(t) = K0 H(t)
the stati displa ement
(t > r/β).
3.5.2 Dipole point sour es A for e dipole an be onstru ted from two opposing single for es whi h are
dipole with moment for whi h the line onne ting the for es is perpendi ular to the dire tion of the for e. The onne ting line for a "dipole without moment " points in the a ting on two neighbouring points. Fig. 3.10 shows, on the left, a
dire tion of the for e.
Fig. 3.10: Single ouple and double ouple onstru ted from single for es.
Two dipoles with moment for whi h the sum over the moments is zero (right in Fig. 3.10), are a good model for many earthquake sour es, i.e., in the ase where the spatial radiation of earthquake waves of su iently large wave length is similar to that of a double ouple model. The a tual pro esses a ting in the earthquake sour e are naturally not four single for es. Usually, the ro k breaks along a surfa e if the shear strength is ex eeded by the a
umulation of shear
shear rupture ).
stress (
Another possibility is that the shear stress ex eeds the
stati fri tion on a preexisting rupture surfa e. Sour e models from single for es and dipoles are only
equivalent point sour es.
In the following, we derive the fareld displa ement of the single ouple model and give the results for the double ouple model.
ouple (with
x0 6= 0)
in Fig.
We start from the single
3.11 and ompute rst from (3.19) the
P wave
3.5.
SPHERICAL WAVES FROM SINGLE FORCE AND DIPOLE ...
displa ement of for e
K(t)
with
x0 )2 + y 2 + z 2 ux uy = uz
Cartesian
z 4πρα2 r 2 K
t−
z
omponents
r
α
57
cos ϑ = z/r, r2 = (x −
(x − x0 )/r y/r · z/r.
P(x,y,z) y r
K(t) υ
λ
x0 ε
x K(t)
Fig. 3.11: Single ouple model.
The displa ements a tion shifted by
u′x , u′y , u′z
ǫ,
of for e
−K(t) with the two neighbouring points of ux , uy , uz
an be determined using the Taylor expansion of
at the sour e oordinate
x0
and trun ating after the linear term. This leads, for
example, to
u′x
∂ux = − ux + ǫ . ∂x0
The single ouple displa ement then follows by superposition
u′′x = ux + u′x = −
∂ux ǫ. ∂x0
To obtain the fareld displa ement requires only the dierentiation of the for e
K(t − r/α) with respe t to r, and additional dierentiation ∂r/∂x0 = −(x − x0 )/r. The other terms with x0 ontribute only to the near eld, the amplitude of whi h de reases faster then 1/r. This leads to
u′′x = −
z r −1 −(x − x0 ) x − x0 ′ t − K ǫ . 4πρα2 r2 α α r r
58
The
CHAPTER 3.
y
and
z displa ement u′′x u′′y u′′z
As expe ted, the
are treated similarly. Therefore,
z(x−x0 ) ′ = − 4πρα t− 3 r3 K
P displa ement
r
α
(x − x0 )/r y/r ǫ· z/r.
K(t),
(3.20)
of the single ouple is also longitudinal.
The for e dipole is dened stri tly by the limit taneous in rease of
BODY WAVES
so that
ǫ → 0,
ombined with a simul
lim K(t) ǫ = M (t)
ǫ→0 remains nite (but nonzero).
M (t)
is alled
moment fun tion
of the dipole
with the dimensions of a rotational moment. From (3.20) with
P wave
z/r = cos ϑ
(x − x0 )/r = sin ϑ cos λ, rdire tion is
and
it follows that the
displa ement of the single ouple in
ur = −
r cos ϑ sin ϑ cos λ ′ . t − M 4πρα3 r α
In on luding, we now assume that
uϑ =
x0 = 0.
For the
S wave,
it follows similarly
sin ϑ sin ϑ cos λ ′ r . t − M 4πρβ 3 r β
As for the single for e, the azimuthal omponent is zero. The following shows the results for the single ouple and the radiation in the
z P y r
υ
λ x
Fig. 3.12: Single ouple in the xzplane.
x − z plane (y = 0)
3.5.
SPHERICAL WAVES FROM SINGLE FORCE AND DIPOLE ...
ur
2ϑ cos λ ′ t− = − sin8πρα 3r M
uϑ
=
uλ
sin2 ϑ cos λ ′ 4πρβ 3 r M
t−
r β
59
r α
(3.21)
= 0.
z
S +  P  +
x
Fig. 3.13: Far eld displa ement of a single ouple.
0 The ratio of the maximum S radiation (for ϑ = 90 ) to the maximum P √ 0 radiation (for ϑ = 45 ) is about 10, if α ≈ β 3. The radiation hara teristi s in planes other then
y = 0
cos λ. z=0
P.
Plane
x=0
is a
is one only for
follow from the one shown by multipli ation with
nodal plane
for
P
as well as for
The fareld displa ements for a double ouple in the
z P y r
λ x
Fig. 3.14: Double ouple in the xzplane.
the plane
x−z plane are (see exer ise
3.6)
υ
S radiation;
60
CHAPTER 3.
ur
=
uϑ
=
uλ
=
2ϑ cos λ ′ − sin4πρα t− 3r M
r α
2ϑ cos λ M′ t − − cos4πρβ 3r
r β
cos ϑ sin λ ′ 4πρβ 3 r M
t − βr .
one
The moment fun tion in (3.22) is that of
ouple. The radiation hara teristi s in the
BODY WAVES
(3.22)
of the two dipoles of the double
x − z plane
are shown in Fig. 3.15.
z
+
 P

+
S x
Fig. 3.15: Far eld displa ement of a double ouple. The
P radiation
of the double ouple has the same form as that of a single
S to P is x = 0 and
ouple but is twi e as large; the ratio of the maximum radiation of now about 5 (for
z = 0.
The
S wave
An (innitesimal) pla ement in
z dire tion,
√ λ ≈ β 3). P nodal
planes are the planes with
has no nodal planes, but only nodal dire tions (whi h?).
shear rupture,
xdire tion
either in the plane
or in the plane
and perpendi ular to it. If, by using the
P
z = 0
with relative dis
with relative displa ement in
radiates waves as a double ouple, i.e., (3.22) holds. A shear rup
ture or earthquake, therefore, radiates no of the
x = 0
P waves in the dire tion of its rupture
distributions of the signs of rst motion
wave, the two nodal planes have been determined, the two possible rup
ture surfa es are found. The determination of the
P nodal plane of earthquakes
3.5.
SPHERICAL WAVES FROM SINGLE FORCE AND DIPOLE ...
fault plane solution ),
(
61
is an important aid in the study of sour e pro esses as
well as the study of larges ale te toni s of a sour e region. Often the de ision between the two options for the rupture surfa e an be made based on geologi al arguments.
The moment fun tion of an earthquake with a smooth rupture, is, to a good approximation, a step fun tion with nonvanishing rise time
M0 ,
the
moment
T
and nal value
of the earthquake (see Fig. 3.16). The fareld displa ements
are then, a
ording to (3.22), onesided impulses.
M M0
0
T
t
T
t
M’
0
Fig. 3.16: Moment fun tion and fareld displa ement of a smooth rupture.
Propagation ee ts in layered media, e.g., the Earth's rust, an hange the impulse form. In reality, the displa ements look very often dierent, relative to the one shown here, due to ompli ated rupture pro esses.
Exer ise 3.6 Derive the double ouple displa ement
ur
in (3.22) from the orresponding single
ouple displa ement in (3.21). Use equation (3.20) in Cartesian oordinates.
62
CHAPTER 3.
BODY WAVES
3.6 Ree tion and refra tion of plane waves at plane interfa es 3.6.1 Plane waves with arbitrary propagation dire tion In se tions 3.1 to 3.3 plane waves travelling in the dire tion of a oordinate axis were used. In the following, we need plane waves with an arbitrary dire tion of propagation. They an be des ribed by the following potentials
Φ
=
− → Ψ =
− → !# − → x k A exp iω t − α " − → !# − → x k − → B exp iω t − n. β "
Their variation with time is also harmoni .
(3.23)
(3.24)
This assumption is su ient for
− → → n are onstant unit ve tors, B are onstant, k and − − → x is the lo ation ve tor, ω is the angular frequen y and i the imaginary unit. − → Φ and the omponents of Ψ satisfy the wave equation (please onrm) most on lusions.
∇2 Φ =
A
and
1 ∂ 2Φ , α2 ∂t2
∇2 Ψj =
1 ∂ 2 Ψj β 2 ∂t2
(Cartesian
oordinates).
Sin e, a
ording to (3.23) and (3.24), the movement at all times and lo ations is nonzero, the wavefronts an no longer be dened as surfa es separating undisturbeddisturbed from disturbed regions. fronts as
We, therefore, onsider wave
− → → x k /c) with c = α or c = β . surfa es of onstant phase ω(t − −
These
surfa es are dened by
d dt
− → − → x k t− c
!
= 0.
− → k,
whi h also gives the dire tion of propagation. The wavefronts move parallel with respe t to themselves with the phase → k multiplied by the wavenumber ω/α or ω/β , is alled the velo ity . Ve tor − wavenumber ve tor. The polarisation dire tion of the ompressional part They are perpendi ular to ve tor
" − → !# − → iω x k − → ∇ Φ = − A exp iω t − k α α
(3.25)
3.6.
63
REFLECTION AND REFRACTION OF PLANE WAVES ...
is longitudinal (parallel to
− → k)
and that of the shear omponent
" − → !# − → iω x k − → − → − ∇ × Ψ = − B exp iω t − k ×→ n β β
(3.26)
− → → → → (rot(f · − n) = f · ∇ × − n −− n × ∇ f ) is transversal (perpendi ular to k ). − → From (3.26), it follows, that for Ψ , without loss of generality, the additional − → − → − →
ondition of orthogonality of k and n an be introdu ed. (Separation of n in − →
omponents parallel and perpendi ular to k ).
3.6.2 Basi equations We onsider a ombination of two halfspa es whi h are separated by a plane at z = 0. The ombination is arbitrary (solidsolid, solidva uum, liquidliquid, ...). We use Cartesian oordinates as shown in Fig. 3.17.
Fig. 3.17: Two halfspa es in Cartesian oordinates.
The
y axis
points out of the plane. The displa ement ve tor is
− → u = (u, v, w), and its omponents are
independent of y , i.e., we treat a plane problem in whi h
on all planes parallel to the
x − z plane,
the same onditions hold. The most
simple way to study elasti waves, under these onditions, is to derive but not
v,
from potentials. Writing
− → − → u = ∇Φ + ∇ × Ψ by omponents,
u = v
=
w
=
∂Φ ∂Ψ2 − ∂x ∂z ∂Ψ3 ∂Ψ1 − ∂z ∂x ∂Φ ∂Ψ2 + , ∂z ∂x
u
and
w
64
CHAPTER 3.
it is obvious that for o
ur in
u and w.
BODY WAVES
v two potentials Ψ1 and Ψ3 are required, and these do not v , it is better to use dire tly the equation of motion (2.21)
For
without body for es, whi h under these onditions be omes a wave equation
∇2 v = The
The
1 ∂2v . β 2 ∂t2
basi equations, therefore, are, if Ψ instead of Ψ2
boundary onditions
∇2 Φ
=
∇2 Ψ
=
∇2 v
=
∇2
=
1 ∂2Φ α2 ∂t2 1 ∂2Ψ β 2 ∂t2 1 ∂2v β 2 ∂t2 ∂2 ∂2 + 2 2 ∂x ∂z
u =
∂Φ ∂x
−
∂Ψ ∂z
w
∂Φ ∂z
+
∂Ψ ∂x .
=
is used
on the surfa e
z=0
(3.27)
(3.28)
between the halfspa es requires
ontinuity of the stress omponents
pzz
=
pzx
=
pzy
=
∂w ∂w → λ∇ · − u + 2µ = λ∇2 Φ + 2µ ∂z ∂z ∂w ∂u µ + ∂x ∂z ∂v µ , ∂z
or
λ ∂2Φ α2 ∂t2
+ 2µ
pzz
=
pzx
2 ∂ Φ + = µ 2 ∂x∂z
pzy
= µ ∂v ∂z .
2 ∂ Ψ − ∂z2
∂2Φ ∂z 2
∂2 Ψ ∂x2
+
∂2Ψ ∂x∂z
(3.29)
Whi h of the displa ement omponents is ontinuous depends on the spe ial
ombination of the halfspa es.
3.6.
65
REFLECTION AND REFRACTION OF PLANE WAVES ...
v with Φ and Ψ exists via the boundary onditions and, u and w, it follows, that the S waves, the displa ement of whi h is only horizontal (in y dire tion: SHwaves ), propagate independently from the P waves, following from Φ, and the S waves, following from Ψ, that also have a verti al omponent (in z dire tion: SVwaves ). If a SHwave impinges on an Sin e no onne tion of
therefore, with
interfa e, only ree ted and refra ted If, on the other hand, a and refra ted
SHwaves
P − (SV −)wave
o
ur, but no
P or SVwaves.
intera ts with an interfa e, ree ted
SV −(P −)waves o
ur, but no SH waves o
ur.
These statements
hold, in general, only for the ase of an interfa e between two solid halfspa es.
SH nor SVwaves propagate; in a va uum a rigid halfspa e,
In liquids, neither
or no waves propagate at all. Correspondingly, the situation is even more simple if su h halfspa es are involved. The
de oupling
of
PSV
SHwaves
and
z = onst
the simple ase of an interfa e
holds for plane problems not only in
between two homogeneous halfspa es,
but also in the more ompli ated ase of an inhomogeneous medium, as long as density, wave velo ity, and module are only fun tions of
x
z.
and
One on
sequen e of this de oupling is that in the following, ree tion and refra tion of
P  and SVwaves an
SHwaves. Swave of arbitrary polarisation in its
be treated independently from that of the
Furthermore, it is possible to disse t an
SV and SH omponent and to study
their respe tive ree tion and refra tion
independently from ea h other. In ea h ase, we assume for the in ident plane wave a potential the form of (3.23) or (3.24), respe tively, (in the se ond ase
y  omponent Ψ).
In ase of a
SHwave, we assume that v
an equation in the form of (3.23) with is part of the dire tion ve tor
− → k
β
instead of
α.
Fig. 3.18: In ident plane wave and angle of in iden e
− → k = (sin ϕ, 0, cos ϕ).
− → Ψ
Φ
or
− → Ψ
in
has only the
an be des ribed by
The angle of in iden e
ϕ
ϕ.
(3.30)
66
CHAPTER 3.
For the ree ted and refra ted wave an ansatz is made with
A and B, respe tively, and dierent
dire tion ve tors
− → k.
BODY WAVES
dierent amplitudes
The relation between
the new dire tion ve tors and (3.30) is via Snell's law. The relation between the displa ement amplitudes of the ree ted and the refra ted wave with the in ident wave, is alled
ree tion oe ient
and
refra tion oe ient,
respe 
tively, and it depends on the angle of in iden e and the material properties in the halfspa es. for
Rpp , Rps , Bpp , Bps , Rss , Rsp , Bss , Bsp will be the oe ients and bss those for the SHwaves. The rst index indi ates
PSVwaves, rss
the type of in ident wave, the se ond the ree ted and refra ted wave type, respe tively. We dis uss, in the following, only relatively simple ases, for whi h illustrate the main ee ts to be studied.
3.6.3 Ree tion and refra tion of SHwaves Ree tion and refra tion oe ients The displa ement
v0
of the in ident
SHwave in y dire tion is
cos ϕ sin ϕ x− z . v0 = C0 exp iω t − β1 β1
v0
(3.31)
v1 ϕ
ϕ1
z=0 ϕ2
µ1,ρ1 ,β1 µ2,ρ2,β2
x
v2 z Fig. 3.19: In ident, ree ted and dira ted SHwaves at a plane interfa e.
The
ansatz for the ree ted and refra ted SHwave as plane waves with ree tion ϕ1 and the refra tion angle ϕ2 , respe tively, and the same frequen y as
angle
the in ident wave is
ree tion
:
v1
=
refra tion
:
v2
=
cos ϕ1 sin ϕ1 x+ z C1 exp iω t − β1 β1 sin ϕ2 cos ϕ2 C2 exp iω t − x− z . β2 β2
(3.32)
(3.33)
3.6.
The unknowns are the angles and the refra tion oe ient The
67
REFLECTION AND REFRACTION OF PLANE WAVES ...
ϕ1 and ϕ2 , the bss = C2 /C0 .
ree tion oe ient
rss = C1 /C0
boundary onditions require at z = 0 the ontinuity of displa ement (that is
a reasonable requirement) and ontinuity for the normal and tangential stresses. This leads to
v0 + v1 = ∂ µ1 ∂z (v0 + v1 ) = The stress omponents
pzz
and
pzx
v2 2 µ2 ∂v ∂z
for
z = 0.
(3.34)
are zero everywhere, sin e no
P and/or SV
wave o
ur. Insert (3.31), (3.32) and (3.33) into (3.34). From the rst boundary
ondition this leads to
sin ϕ sin ϕ1 sin ϕ2 C0 exp iω t − x +C1 exp iω t − x = C2 exp iω t − x . β1 β1 β2 (3.35)
We plan to nd solutions
C1
and
and
v2
C2
v1
and
v2
of the problem, for whi h the amplitudes
are independent of lo ation, sin e only then an we be sure that
are solutions of the orresponding wave equation.
C1
and
C2
v1
be ome
only independent of lo ation if in (3.35)
sin ϕ1 sin ϕ2 sin ϕ = = , β1 β1 β2
(3.36)
sin e only then the exponential term an be an elled. Equation (3.36) is the wellknown
Snell's Law
angle of in iden e
ϕ
whi h states that the ree tion angle
and that for the refra tion angle
ϕ2
ϕ1
is equal to the
is
β2 sin ϕ2 . = sin ϕ β1 With (3.35), this leads to
C2 − C1 = C0 . The se ond boundary ondition in (3.34) gives
cos ϕ cos ϕ2 cos ϕ1 µ1 iω − C0 + C1 = −µ2 iω C2 . β1 β1 β2 With
ϕ1 = ϕ
and
µ1,2 /β1,2 = ρ1,2 β1,2 ,
it follows that
(3.37)
68
CHAPTER 3.
BODY WAVES
ρ1 β1 cos ϕ(C1 − C0 ) = −ρ2 β2 cos ϕ2 C2 or
ρ2 β2 cos ϕ2 C2 + C1 = C0 . ρ1 β1 cos ϕ
(3.38)
From (3.37) and (3.38) follow the ree tion and refra tion oe ients
C1 C0 C2 = C0
rss =
=
bss
=
ρ1 β1 cos ϕ − ρ2 β2 cos ϕ2 ρ1 β1 cos ϕ + ρ2 β2 cos ϕ2 2ρ1 β1 cos ϕ . ρ1 β1 cos ϕ + ρ2 β2 cos ϕ2
(3.39)
(3.40)
With (3.36), this leads to
12 β22 2 cos ϕ2 = (1 − sin ϕ2 ) = 1 − 2 sin ϕ . β1 1 2
2
For
perpendi ular
in iden e (ϕ
rss =
In this ase,
rss
and
two halfspa es. For
bss
(3.41)
= 0)
ρ 1 β1 − ρ 2 β2 ρ 1 β1 + ρ 2 β2
and
bss =
2ρ1 β1 . ρ 1 β1 + ρ 2 β 2
impedan es ρ1 β1 and ρ2 β2 of the = π/2), rss = −1 and bss = 0. The
depend only on the
grazing
in iden e (ϕ
absolute value of the amplitude of the ree ted wave is never larger then that of the in ident wave; that of the refra ted wave an be larger if for If
ρ 2 β2 < ρ 1 β1
(e.g.,
ϕ = 0).
rss
is negative, this means that in one point of the interfa e the displa ement
−y dire tion, if the displa ement ve tor +y dire tion. For impulsive ex itation (see also
ve tor of the ree ted wave points in of the in ident wave points in
later), this means that the dire tion of rst motion of the in ident and the ree ted wave are opposite. The following gure shows
β1 /β2 > 1
and
ρ1 = ρ2 .
rss 
as a fun tion of
ϕ
for dierent velo ity ratios
3.6.
REFLECTION AND REFRACTION OF PLANE WAVES ...
Fig. 3.20:
rss 
as a fun tion of
ϕ
69
for dierent velo ity ratios.
Total ree tion If
β2 < β1 , as in Fig. 3.20, cos ϕ2 is real for all angles of in ident ϕ, the same rss and bss . Total ree tion,, i.e., rss  = 1, is then only possible for
is true for
grazing in iden e. If
β2 > β1 , cos ϕ2
is only real as long as
ϕ =< ϕ∗ = arcsin
β1 . β2
ϕ∗ is the riti al angle (or limiting angle of total ree tion ). A
ording to (3.41), ϕ = ϕ∗ is onne ted to the ase with grazing propagation of the wave in the se ond halfspa e (ϕ2 = π/2). ϕ > ϕ∗ , cos ϕ2 be omes imaginary, or, to be more exa t, negative imaginary ω and positive imaginary for negative ω , sin e only then v2 for z → +∞ remains limited. rss and bss be ome omplex. v1 and v2 still solve
If
for positive
the wave equations and satisfy the boundary onditions, even when posing the
ansatz (3.32) and (3.33) have not expli itly been hosen in omplex form. The ree tion oe ient an then be written as
= exp −2i arctan ab
rss
=
a−ib a+ib
a
=
ρ1 β1 cos ϕ
b
=
−ρ2 β2
β22 β12
sin2 ϕ − 1
12
ω ω
.
(3.42)
70
CHAPTER 3.
BODY WAVES
It has the absolute value 1 and a phase that depends on the angle of in iden e.
Its sign hanges with the sign of the frequen y.
Values of
3.21 for a few ombinations.
Fig. 3.21:
rss 
as a fun tion of
ϕ
β1 / sin ϕ.
are shown in Fig.
for dierent velo ity ratios.
The refra ted wave propagates for lo ity
rss 
ϕ > ϕ∗
parallel to the interfa e with the ve
Its amplitude is not only ontrolled by
by the exponential term whi h depends on
z.
bss , but is also ontrolled
The amplitude of the refra ted
wave de ays, therefore, exponentially with in reasing distan e from the interfa e
inhomogeneous
(
or
boundary layer wave ). "
ω v2 = bss C0 exp − β2
It follows that (please he k)
β22 sin2 ϕ − 1 β12
21 # sin ϕ x . z exp iω t − β1
Other ases The treatment of the ree tion of plane
liquids
P waves
at an interfa e between two
gives similar results to the one dis ussed above (see also exer ise 3.9).
If the interfa e between two
solid halfspa es
is onsidered, the omputational
eort is signi antly larger, sin e now ree ted and refra ted
SV waves have to
be in luded. We, therefore, skip the details. The absolute value of the ree tion
oe ient
Rpp
is shown in Fig. 3.22.
3.6.
 R pp 
 R pp 
α1 sin ϕ∗ = α −
β 2 < α 1 < α2
0
π/2
ϕ∗
0
ϕ
ϕ = 0, Rpp =
Rpp 
for
SV wave SV wave
1 sin ϕ∗∗ = α β−2
α 1 < β 2 < α2
2
0
ϕ∗
0
Fig. 3.22: Absolute value of the ree tion oe ient
For
71
REFLECTION AND REFRACTION OF PLANE WAVES ...
ϕ∗∗
π/2
ϕ
Rpp .
ρ2 α2 −ρ1 α1 ρ2 α2 +ρ1 α1 ( ompare also exer ise 3.9).
ϕ∗ < ϕ < π/2
is smaller then 1 for two reasons. First, the ree ted
also arries energy; se ond, for the ase on the left of Fig.
all
3.22, a
propagates in the lower halfspa e for ϕ, and, similarly, for the ∗∗ ∗∗
ase on the right of Fig. 3.22, for ϕ < ϕ . ϕ < ϕ is the whi h exists only for
ϕ∗∗ = arcsin
For angles
ϕ
larger then
total ree tion o
urs. transported in the
se ond riti al angle
α1 < β2 < α2
ϕ∗∗ ,
α1 α1 > ϕ∗ = arcsin . β2 α2
the se ond energy loss is no longer possible, and
The ree ted energy is then, to a smaller part, also
SV wave.
Some numeri al results for ree tion and refra tion oe ients for a
PSV ase
are given in Fig. 3.23 (model of the rustmantle boundary (Moho) with α1 = 6.5km/sec, β1 = 3.6km/sec, ρ1 = 2.8g/cm3, α2 = 8.2km/sec, β2 = 4.5km/sec, ρ2 3.3g/cm3).
=
72
CHAPTER 3.
BODY WAVES
Fig. 3.23: Absolute value of ree tion and refra tion oe ients and
Rpp , Rps , Bpp
Rss .
Transition to impulsive ex itation The transition from the harmoni ase, treated up to now, to the impulse ase,
an be done with the (3.31), the
SHwave
Fourier transform
v0 = F
( ompare appendix A.1.7). Instead of
sin ϕ cos ϕ t− x− z β1 β1
3.6.
73
REFLECTION AND REFRACTION OF PLANE WAVES ...
may impinge on the interfa e, and we disse t integral in
partial vibration
F (t) F (ω)
Z +∞ 1 F (ω)eiωt dω 2π −∞ Z +∞ = F (t)e−iωt dt
F (t)
with the aid of the Fourier
=
F (t).
Fourier transform of
−∞
We then study, as before, ree tion and refra tion of the
partial waves
1 sin ϕ cos ϕ dv0 = F (ω) exp iω t − x− z dω 2π β1 β1 and then sum the ree ted partial waves to derive the ree ted
v1 =
1 2π
Z
+∞
−∞
SH wave
cos ϕ sin ϕ x+ z dω. rss F (ω) exp iω t − β1 β1
(3.43)
As long as the ree tion oe ient rss is frequen y independent (whi h is the β2 < β1 or for ϕ < ϕ∗ with β2 > β1 ), it an be moved before the
ase for
integral, thus, yielding
v1 = rss F
sin ϕ cos ϕ t− x+ z . β1 β1
The ree ted impulse has, in this ase, the same form as the in ident impulse. The amplitude ratio of the two impulses is equal to the ree tion oe ient.
rss , β2 > β1 . Then
a
ording to (3.42), be omes dependent from One then has to pro eed dierently.
ω
We disse t
for
ϕ > ϕ∗
rss
into real and
with
imaginary parts
rss
=
R(ϕ) = I(ϕ)
=
R(ϕ) + iI(ϕ) a2 − b 2 a2 + b 2 2ab − 2 a + b2
ω ω
(b
for
ω > 0).
A
ording to (3.43), it holds that
1 v1 = R(ϕ)F (τ ) + I(ϕ) 2π
Z
+∞
−∞
iω F (ω)eiωτ dω ω
(3.44)
74
CHAPTER 3.
BODY WAVES
with
τ =t−
sin ϕ cos ϕ x+ z. β1 β1
I(ϕ) in (3.44) (iω/ ω) · F (ω).
The fun tion, with whi h
is multiplied, an be written, here via
its Fourier transform, as
This is a simple
lter
of fun tion
F (τ ) ω in +900
( ompare general omments on lters in appendix A.3.4). Ea h frequen y
F (τ ) keeps its amplitude, but its phase is hanged. The phase hange is 0 for ω > 0 and −90 for ω < 0. This orresponds to a Hilbert transform and is shown in appendix B. The fun tion with whi h I(ϕ) in (3.44) is multiplied is, therefore, the Hilbert transform FH (τ ) of F (τ )
1 FH (τ ) = P π
P
Z
+∞
−∞
1 F (t) dt = t−τ π
Z
+∞
−∞
ln t F ′ (τ − t)dt.
indi ates the main value (without the singularity at
form of
FH (τ )
t = τ ),
(3.45)
and the se ond
follows from the rst by partial integration. Thus,
v1 = R(ϕ)F (τ ) + I(ϕ)FH (τ ).
(3.46)
Due to the se ond term in (3.46), the form of the ree ted wave is from that of the in ident wave. Fig. 3.24 shows the 0 waves and angle of in iden e ϕ from 0 to 90 .
SH
For pre riti al angles of in iden e
ϕ < ϕ∗ = 480 ,
results
dierent
of the ree tion of
the ree tion has the form of
the in ident wave with positive and negative signs. Beyond the riti al angle, 0 in the range of total ree tion, impulse deformations o
ur until at ϕ = 90 the in ident wave form appears again, but with opposite sign ( orresponding to a 0 0 ree tion oe ient rss = −1). The phase shift of rss at ϕ = 55 is about ±90 , with the onsequen e that
R(ϕ) ≈ 0.
The ree tion impulse for this angle of
in iden e is, therefore, lose to the Hilbert transform
FH (τ )
of
F (τ )
(the exa t
Hilbert transform is an impulse that is symmetri with respe t to its minimum).
3.6.
REFLECTION AND REFRACTION OF PLANE WAVES ...
Fig. 3.24: Ree tion of
SHwaves for dierent angle of in iden e ϕ.
75
76
CHAPTER 3.
BODY WAVES
Exer ise 3.7: Whi h sign does the ree tion oe ient
rss
have in Fig. 3.20 and Fig. 3.21, in
the regions where it is real?
Exer ise 3.8: Determine the angle of in iden e for whi h
rss
Brewster angle),
is zero (
give the onditions under whi h this a tually happens ( ompare 3.24).
ϕ ≈ 400
and
in Fig.
Exer ise 3.9 Compute the ree tion and refra tion oe ients for a plane surfa e between two liquids and for a plane harmoni longitudinal wave under angle of in iden e
ϕ impinges. Give, qualitatively, the trend of the oe ients for ρ1 = ρ2 with α1 > α2 and α1 < α2 . Hint: Use an ansatz for the displa ement potential in the form of (3.31) to (3.33) and express the boundary onditions via potentials as dis ussed in se tion 3.6.2.
3.6.4 Ree tion of Pwaves at a free surfa e Ree tion oe ients Pwaves from a free surfa e is of pra ti al imP waves from earthquakes and explosions propagate
The study of the ree tion of portan e for seismology.
through the Earth and impinge at the seismi station from below. tal and verti al displa ement are modied by the free surfa e. ree ted
P
and
S waves
Horizon
Furthermore,
are ree ted downwards and re orded at larger dis
tan es, sometimes with large amplitudes. It is, therefore, useful and ne essary to know the ree tion oe ient of the Earth's surfa e. For the moment, we negle t the layered nature of the rust in our model, thus, only giving a rst approximation to reality. Based on the omments given at the end of se tion 3.6.2, we sele t the following ansatz for the potentials
in ident
P − wave Φ0
ree ted
sin ϕ cos ϕ A0 exp iω t − x− z α α
(3.47)
=
sin ϕ1 cos ϕ1 A1 exp iω t − x+ z α α
(3.48)
=
sin ϕ′1 cos ϕ′1 B1 exp iω t − x+ z . β β
(3.49)
P − wave Φ1
ree ted
=
SV − wave Ψ1
3.6.
REFLECTION AND REFRACTION OF PLANE WAVES ...
Fig. 3.25: In ident
77
Pwave and ree ted P and Swave .
The boundary onditions at z = 0 require vanishing normal and tangential stress pzz = pzx = 0. No boundary onditions for the displa ement exist. With (3.29) − → and Φ = Φ0 + Φ1 (Ψ = Ψ1 = y − omponent of Ψ ), it follows that
2µ ∂ 2 ∂ 2 Ψ1 1 ∂2 (Φ + Φ ) + (Φ + Φ ) + 0 1 0 1 α2 ∂t2 λ ∂z 2 ∂x∂z 2 ∂ ∂ 2 Ψ1 ∂ 2 Ψ1 2 (Φ0 + Φ1 ) + − ∂x∂z ∂x2 ∂z 2
=
0
z=0
=
0
z = 0. (3.51)
(3.50)
As in the last se tion, Snell's law follows from the boundary onditions
From this, it follows that
sin ϕ1 sin ϕ′1 sin ϕ = = . α α β β · sin ϕ < ϕ. ϕ1 = ϕ and ϕ′1 = arcsin α
(3.52)
With (3.47), (3.48), (3.49) and
µ β2 µ ρβ 2 = = = λ λ + 2µ − 2µ ρα2 − 2ρβ 2 α2 − 2β 2 (3.50) leads to
1 (A0 + A1 ) (iω)2 α2 " 2 # 2β 2 iω iω iω ′ ′ + 2 (A0 + A1 ) cos ϕ + B1 − sin ϕ1 cos ϕ1 = 0. α − 2β 2 α β β Then
78
CHAPTER 3.
BODY WAVES
2α2 β 2 sin ϕ′1 cos ϕ′1 cos2 ϕ A0 + A1 + 2 − B1 (A0 + A1 ) = 0. α − 2β 2 α2 β2 with
1+
2β 2 cos2 ϕ = α2 − 2β 2 = = =
and
γ=
α2 β2
>2
2 2β 2 α 2 − 1 + cos ϕ α2 − 2β 2 2β 2 2 α 2β 2 2 − sin ϕ α2 − 2β 2 2β 2 2 α β2 2 − 2 sin ϕ α2 − 2β 2 β 2
γ − 2 sin2 γ γ−2
, it follows that
γ − 2 sin2 ϕ 2γ sin ϕ′1 cos ϕ′1 (A0 + A1 ) − B1 = 0. γ−2 γ−2 From this
(γ − 2 sin2 ϕ)
1 B1 A1 − 2 sin ϕ(γ − sin2 ϕ) 2 = 2 sin2 ϕ − γ. A0 A0
(3.53)
Equation (3.51) then gives
iω iω iω iω − cos ϕ + 2A1 − sin ϕ cos ϕ 2A0 − sin ϕ α α α α 2 2 iω iω +B1 − sin ϕ′1 − B1 cos ϕ′1 =0 β β or
sin2 ϕ′1 − cos2 ϕ′1 2 sin ϕ cos ϕ (A0 − A1 ) + B1 = 0. 2 α β2 Equation (3.52) then gives
2 sin ϕ cos ϕ
B1 A1 + γ − 2 sin2 ϕ = 2 sin ϕ cos ϕ. A0 A0
(3.54)
3.6.
79
REFLECTION AND REFRACTION OF PLANE WAVES ...
From (3.53) and (3.54), it follows that the amplitude ratios are
1 2 4 sin2 ϕ cos ϕ γ − sin2 ϕ 2 − γ − 2 sin2 ϕ A1 = 2 1 A0 4 sin2 ϕ cos ϕ γ − sin2 ϕ 2 + γ − 2 sin2 ϕ
(3.55)
4 sin ϕ cos ϕ γ − 2 sin2 ϕ B1 = 1 2 . A0 4 sin2 ϕ cos ϕ γ − sin2 ϕ 2 + γ − 2 sin2 ϕ
To derive displa ement amplitudes (that is how the oe ients
(3.56)
Rpp
and
Rps
in
se tion 3.6.2 were dened) from the ratios of potential amplitudes given here,
P wave is PP ree tion
we use (3.25) and (3.26). The displa ement amplitude of the in ident
iω − iω α A0 ; that of the ree ted P wave is − α A1 .
This then gives the
oe ient (see also (3.55))
Rpp =
A1 . A0
(3.57)
Equation (3.26) gives the displa ement amplitude of the ree ted
− iω β B1 .
Thus, the
PS ree tion oe ient is (see also (3.56)) Rps =
Rpp
and
Rps are real
and
is always positive. For and only a
α B1 . β A0
and
ϕ=
π 2,
Rpp = −1
Fig. 3.26: Ree tion and refra tion oe ients of
ϕ.
(3.58)
frequen y independent for all angles of in ident ϕ. Rps
ϕ=0
P wave is ree ted.
of in iden e
SV wave as
and
Rps = 0,
respe tively,
Pwaves for dierent
angles
80
CHAPTER 3.
The meaning of the
negative signs
BODY WAVES
in the ree tion oe ients be omes lear, if
the displa ement ve tor of the in ident and ree ted waves are represented via (3.25) and (3.26)
− → u0
=
− → u1
=
− → u ′1
=
sin ϕ cos ϕ iω − → k0 ∇ Φ0 = − A0 exp iω t − x− z α α α iω sin ϕ cos ϕ − → ∇ Φ1 = − A0 Rpp exp iω t − x+ z k1 α α α sin ϕ′1 cos ϕ′1 iω − →′ − − → x+ z k1×→ n. ∇ × Ψ = − A0 Rps exp iω t − α β β
Fig. 3.27: Polarity of ree ted
Rpp < 0,
P and SVwaves.
therefore, means that if the displa ement of the ree ted
P wave in a
− → = 0) points in the dire tion of − k 1 , the in ident wave − → points in the dire tion of k 0 . For Rps < 0, the displa ement of the ree ted SV − →′ → wave would, for su h an in ident wave, be pointing in the dire tion of − k 1 × − n. point on the interfa e (z
These onne tions be ome more obvious if we go from the harmoni ase to the
impulsive ase
( ompare se tion 3.6.3, transition to impulse ex itation).
For the problem studied, the ree tion oe ients are frequen y independent. Therefore, the ree ted waves have always the same form as the in ident wave
− → u0 − → u1
sin ϕ cos ϕ − → = F t− x− z k0 α α cos ϕ sin ϕ − → x+ z k1 = Rpp F t − α α
(3.59)
(3.60)
3.6.
REFLECTION AND REFRACTION OF PLANE WAVES ...
In the ase
− → u ′1
=
ϕ′1
=
sin ϕ′1 cos ϕ′1 − →′ → Rps F t − x+ z k1×− n β β β sin ϕ . arcsin α
(3.61)
Rpp < 0, if the rst motion of the in ident P wave is dire ted towards z = 0, this also holds for the ree ted SV wave and the ree ted
the interfa e
P wave;
81
otherwise, the rst motion of the ree ted
the interfa e. Fig. 3.28 shows the ase for
P wave
points away from
Rpp < 0.
Fig. 3.28: Denition of the rst motion of ree ted
P and SVwaves.
Displa ements at the surfa e Finally, we ompute the resulting displa ement at the free surfa e (z=0) in whi h the three waves (3.59), (3.60) and (3.61) superimpose.
Horizontal displa ement u = u = fu (ϕ)
and
=
(positive in
Rps cos ϕ′1 ] F
[(1 + Rpp ) sin ϕ + sin ϕ fu (ϕ)F t − x α
=
sin ϕ x t− α
1 4γ sin ϕ cos ϕ γ − sin2 ϕ 2 2 1 4 sin2 ϕ cos ϕ γ − sin2 ϕ 2 + γ − 2 sin2 ϕ
Verti al displa ement
w
x dire tion):
(positive in
z dire tion):
[(1 − Rpp ) cos ϕ +
Rps sin ϕ′1 ] F
sin ϕ t− x α
(3.62)
82
CHAPTER 3.
w fw (ϕ)
BODY WAVES
sin ϕ = fw (ϕ)F t − x α =
2γ cos ϕ γ − 2 sin2 ϕ 2 . 1 4 sin2 ϕ cos ϕ γ − sin2 ϕ 2 + γ − 2 sin2 ϕ
The ampli ation fa tors (or transfer fun tions of the surfa e) respe tively, are given in Fig. 3.29 for the ase
(3.63)
fu (ϕ) and fw (ϕ),
γ = 3.
Fig. 3.29: Transfer fun tions of the free surfa e.
linearly polarised wave with the apparent velo ity α/ sin ϕ propapolarisation angle ǫ, (see Fig. 3.30), is not identi al to the angle of in iden e ϕ. ǫ is also alled the apparent angle of in iden e. Therefore, a
gates at the surfa e. The
ǫ = arctan
u w
12 2 sin ϕ γ − sin ϕ . = arctan γ − 2 sin2 ϕ 2
3.6.
REFLECTION AND REFRACTION OF PLANE WAVES ...
Fig. 3.30: Polarisation angle
ǫ
and angle of in iden e
ε
ϕ.
ε=ϕ
π 2
ε(ϕ)
π 2
O Fig. 3.31: Qualitative relationship between
ǫ
83
and
ϕ
.
1
π
= arctan
ǫ′ (0)
=
2
ǫ
ϕ
2 1
γ2
2 (γ − 1) 2 γ−2
!
.
In ident SV wave SV wave, instead of the P wave onsidered up until now, impinges on the P wave is ree ted for angles of in iden e ϕ > ϕ∗ = arcsin αβ , but only an SV wave (Rss  = 1) is ree ted. This follows from onsiderations ∗ similar to that for an in ident P wave. For ϕ < ϕ , the displa ement at the ∗ free surfa e is linearly polarised, but for ϕ > ϕ , it is polarised ellipti ally. If a
free surfa e, no
84
CHAPTER 3.
BODY WAVES
SV
This property is observed: waves from earthquakes for distan es smaller 0 then about 40 are ellipti ally polarised, but are linearly polarised for larger distan es.
Fig. 3.32: Polarisation of
SV waves from earthquakes .
3.6.5 Ree tion and refra tion oe ients for layered media Matrix formalism In the last two se tions, we studied the ree tion and refra tion of plane waves at
one interfa e.
The ree tion and refra tion oe ients depend, then, mainly on
the properties of the halfspa es and the angle of in iden e. Only if the riti al angle is ex eeded, a weak frequen y dependen e o
urs: the sign of the phase (the oe ients be ome omplex) is ontrolled by the sign of the frequen y of the in ident wave ( ompare se tion 3.6.3). The frequen y dependen e be omes mu h more pronoun ed when the ree tion and refra tion of plane waves in
layered media is onsidered (two or more interfa es). Then, interferen e phenomena o
ur and for spe ial frequen ies (or wave lengths) onstru tive or destru tive interferen es o
ur. a (subparallel)
generally,
Here, we will study the ree tion and refra tion of
liquid layers between two liquid halfspa es.
P waves
from a pa ket of
The orresponding problem for
SH 
waves in solid media an be solved similarly. There is a lose similarity between
P waves
SH waves in layered solid media. The PSVwaves in solid media (possibly with interspersed liquid layers)
in layered liquid media and
treatment of
is, in prin iple, the same, but the derivation is signi antly more ompli ated. In all these approa hes, a
matrix formalism
is used, whi h is espe ially ee tive
for implementing on omputers. We hoose the annotation of the liquidlayered medium as given in Fig. 3.33.
3.6.
85
REFLECTION AND REFRACTION OF PLANE WAVES ...
Fig. 3.33: Liquidlayered medium with n layers. The displa ement potential
Φj
in the
j th
layer
(j = 1, 2, . . . , n)
satises the
wave equation
∂ 2 Φj 1 ∂ 2 Φj ∂ 2 Φj + = 2 . 2 2 ∂x ∂z αj ∂t2 Solutions of this equation, whi h an be interpreted as harmoni plane waves, have the form
exp [i (ωt ± kj x ± lj z)] with
kj2 + lj2 = ω/α2j ,
where
kj
is the
horizontal, lj
respe tively.
We assume positive frequen ies
wavenumbers
kj .
ω
the
verti al wavenumber,
and nonnegative horizontal
Then, we an disregard the sign + of
kj x,
sin e it orre
sponds to waves whi h propagate in xdire tion. This is not possible for our sele tion of the in ident wave, (see Fig.
3.33).
The two signs of
lj z
z
have to
be kept, sin e in all layers (ex ept the nth) waves propagate in +  and in
z dire tion.
We then ome to the
Φj
=
z1 Bn
potential ansatz
(3.64)
=
Aj exp [i (ωt − kj x − lj (z − zj ))] +Bj exp i ωt − kj′ x + lj′ (z − zj )
z2 = 0,
=
0.
(3.66)
kj2 + lj2 = kj′2 + lj′2 =
ω2 α2j
(3.65)
86
CHAPTER 3.
BODY WAVES
In (3.64), we have assumed, for the moment, that the wavenumbers of the waves propagating in +z and zdire tion are dierent. Furthermore, we have repla ed
z
by
z − zj .
This does not hange the meaning of the terms but simplies the
omputations. The part
A1 exp [i (ωt − k1 x − l1 z)] of Φ1
( ompare, e.g., (3.47)). This means that of in iden e
ϕ
will be interpreted as in ident
k1
P wave
and l1 are onne ted with the angle
as
k1
=
ω α1
l1
=
ω α1
sin ϕ cos ϕ.
B1 exp [i (ωt − k1′ x + l1′ z)] of Φ1 is ered halfspa e z > 0. We want to ompute the refra tion oe ient Bpp (again dened The part
amplitudes )
Rpp
=
B1 A1
Bpp
=
α1 αn
·
The boundary onditions for the interfa es
(3.67)
the wave ree ted from the laythe ree tion oe ient as the ratio of the
An A1 .
Rpp
and
displa ement
(3.68)
z = z 2 , z3 , . . . , zn
require ontinuity pzz = λ∇2 Φ =
of the verti al displa ement ∂Φ/∂z and of the normal stress ρ∂ 2 Φ/∂t2 . For z = zj , this gives
∂Φj ∂z
∂Φj−1 ∂z
=
and
ρj
∂ 2 Φj ∂t2
=
ρj−1
From the rst relation, it follows that (the phase term
∂ 2 Φj−1 ∂t2 .
eiωt
is negle ted in the
following sin e it an els out),
−lj Aj exp [−ikj x] + lj′ Bj exp −ikj′ x The se ond relation gives
ρj Aj exp [−ikj x] + ρj Bj exp −ikj′ x Both equations hold for
= −lj−1 Aj−1 exp [i (−kj−1 x − lj−1 dj−1 )] ′ ′ ′ + lj−1 Bj−1 exp i −kj−1 x + lj−1 dj−1 .
= +
ρj−1 Aj−1 exp [i (−kj−1 x − lj−1 dj−1 )] ′ ′ ρj−1 Bj−1 exp i −kj−1 x + lj−1 dj−1 .
j = 2, 3, . . . , n, and dj−1 = zj −zj−1 (d1 = 0).
x
As before,
we require that the exponential terms depending on must an el, leading to ′ kj = kj′ = kj−1 = kj−1 . This, then, gives (with (3.67))
3.6.
87
REFLECTION AND REFRACTION OF PLANE WAVES ...
′ kn′ = kn = kn−1 = kn−1 = . . . = k1′ = k1 =
This is an alternative form of
lj′
Snell's law.
With (3.65), this leads to
! 12
! 21 α2j 2 . 1 − 2 sin ϕ α1
ω2 − k12 α2j
= lj =
ω sin ϕ. α1
ω = αj
(3.69)
If sin ϕ > α1 /αj , lj is imaginary (and even negative imaginary ), only then for j = n is the amplitude of the potentials limited for z → ∞. This leads to the following system of equations, whi h onne ts Aj and Bj with Aj−1 and Bj−1 , respe tively
Aj − Bj
=
Aj + Bj
=
lj−1 Aj−1 e−ilj−1 dj−1 − Bj−1 eilj−1 dj−1 lj ρj−1 Aj−1 e−ilj−1 dj−1 + Bj−1 eilj−1 dj−1 . ρj
In matrix form, this an be written as (please he k)
Aj Bj
where
mj
e−ilj−1 dj−1 lj−1 ρj + lj ρj−1 = −l 2lj ρj j−1 ρj + lj ρj−1 Aj−1 · Bj−1 Aj−1 = mj · Bj−1 is the
(−lj−1 ρj + lj ρj−1 )e2ilj−1 dj−1 (lj−1 ρj + lj ρj−1 )e2ilj−1 dj−1
(3.70)
layer matrix.
Repeated appli ation of (3.70) gives
An Bn
= = =
On omputers, the produ t
M
mn · mn−1 · . . . · m3 · m2 A1 M B1 A1 M11 M12 . B1 M21 M22 of the layer matri es
qui kly and e iently. First, the angular frequen y
A1 B1
mn to m2 an be determined ω and the angle of in iden e
88
ϕ
CHAPTER 3.
are given; then, the
lj 's
BODY WAVES
are determined with (3.69), and the matri es are
multiplied. This gives the elements of
An = M11 A1 + M12 B1
M.
and
From
Bn = M21 A1 + M22 B1
with (3.66), it follows that
M21 B1 =− A1 M22 The ree tion oe ient
Rpp
and
M12 M21 An = M11 − . A1 M22 Bpp
of the layered
M12 M21 M11 − . M22
(3.71)
and the refra tion oe ient
medium, therefore, an be written a
ording to (3.68) as
Rpp = −
M21 M22
and
α1 αn
Bpp =
Two homogeneous halfspa es In this very simple ase, it follows (with
M = m2 =
1 2l2 ρ2
d1 = 0)
l 1 ρ2 + l 2 ρ1 −l1 ρ2 + l2 ρ1
that
−l1 ρ2 + l2 ρ1 l 1 ρ2 + l 2 ρ1
and, therefore, a
ording to (3.71)
With
Rpp
=
Bpp
=
l1 =
ω α1
−l2 ρ1 + l1 ρ2 l 2 ρ1 + l 1 ρ2 α1 (l1 ρ2 + l2 ρ1 )2 − (l2 ρ1 − l1 ρ2 )2 2l1 ρ1 α1 . = α2 2l2 ρ2 (l2 ρ1 + l1 ρ2 ) α2 l2 ρ1 + l1 ρ2
cos ϕ
and l2
=
ω α2
refra tion), it follows that
Rpp
=
Bpp
=
( ompare with exer ise 3.9). For
1−
α22 α21
12 sin2 ϕ =
ω α2
cos ϕ2 (ϕ2 =angle
ρ2 α2 cos ϕ − ρ1 α1 cos ϕ2 ρ2 α2 cos ϕ + ρ1 α1 cos ϕ2 2ρ1 α1 cos ϕ ρ2 α2 cos ϕ + ρ1 α1 cos ϕ2 ϕ = 0 (→ ϕ2 = 0),
it follows that
of
3.6.
89
REFLECTION AND REFRACTION OF PLANE WAVES ...
2ρ1 α1 . ρ2 α2 + ρ1 α1
(3.72)
These are equations that also hold for an interfa e between two
solid halfspa es.
Rpp =
ρ2 α2 − ρ1 α1 ρ2 α2 + ρ1 α1
and
Bpp =
Lamella in fullspa e We limit our study here to
verti al ree tions
from a lamella.
α1 ,ρ1
z 2= 0
x
α2,ρ2
d 2= d z
α3=α1 ,ρ3= ρ1
Fig. 3.34: Lamella of thi kness d. In this ase,
n = 3, l1 = l3 = ω/α1
and l2
= ω/α2 .
Then with (3.70) and
d1 = 0,
d2 = d
m2
=
α2 2ρ2
=
α2 2α1
−iω αd
α1 e m3 = 2α2 with
γ=
ρ1 α1 ρ2 α2 and
γ′ = d
e−iω α2 M = m3 m2 = 4γ We now ompute the
ρ2 α2 ρ1 α1
2
=
ρ2 α1 − αρ21
+ +
1+γ −1 + γ
1 + γ′ −1 + γ ′
ρ1 α2 ρ1 α2
− αρ21 + αρ12 ρ2 ρ1 α1 + α2 −1 + γ 1+γ
d
(−1 + γ ′ )e2iω α2 2iω d (1 + γ ′ )e α2
d
(1 + γ)2 − (1 − γ)2 e2iω α2 2iω d 1 − γ 2 + (γ 2 − 1)e α2
ree tion oe ient Rpp
d
γ 2 − 1 + (1 − γ 2 )e2iω α2 2iω d −(1 − γ)2 + (1 + γ)2 e α2
(a
ording to (3.71))
d
=
M21 1 − e−2iω α2 (1 − γ 2 )(1 − e2iω α2 ) − = = R0 d 2iω −2iω αd M22 2 (1 − γ)2 − (1 + γ)2 e α2 1 − R02 e
=
1−γ ρ2 α2 − ρ1 α1 . = 1+γ ρ2 α2 + ρ1 α1
with
R0
!
1 γ . This leads to
d
Rpp
(3.73)
!
.
90
R0 ,
CHAPTER 3.
BODY WAVES
a
ording to (3.72), is the ree tion oe ient of the interfa e
relatively small ree tion oe ients the Earth
(R0  < 0.2),
R0 ,
z = 0.
For
whi h are typi al for dis ontinuities in
one an write as a good approximation
−2iω αd 2 Rpp = R0 1 − e .
Dis ussion of Rpp
(3.74)
Rpp , in the form of (3.73) or (3.74), is zero for angular frequen ies ω , for whi h 2ω αd2 is an even multiple of π . With the frequen y ν and the wave length Λ in the lamella (α2 = νΛ), the ondition for destru tive interferen e is 1 3 d = , 1, , 2, . . . Λ 2 2
(3.75)
The lamella has to have a thi kness of a multiple of the half wavelength so
that in ree tion destru tive interferen e refra tions show onstru tive interferen e. A
ording to (3.74), multiple of
π.
Rpp
is maximum
Then
o
urs with
Rpp = 0.
(Rpp  = 2 R0 )
if
2ω αd2
In this ase
is an uneven
d 1 3 5 7 = , , , ,... Λ 4 4 4 4 In this ase, the waves interfere
for refra tion.
The periodi ity of
Rpp ,
(3.76)
onstru tively for ree tion
and
destru tively
visible in (3.75) and (3.76), holds generally so
α2 π Rpp ω + n = Rpp (ω), d
n = 1, 2, 3, . . .
impulP wave has the verti al
To on lude, we dis uss how the ree tion from a lamella looks for an
sive ex itation.
w0 = F t − displa ement w1 of
displa ement verti al
We assume that the verti ally in ident
w1 = with
Rpp (ω)
z α1
and that
F (ω)
is the spe trum of
F (t).
The
the ree ted wave is then ( ompare se tion 3.6.3)
1 2π
Z
+∞
Rpp (ω)F (ω)e
iω t+ αz
1
−∞
dω
(3.77)
from (3.73). In pra tise, integral (3.77) is omputed numeri ally,
sin e fast numeri al methods for
Fourier analysis
the spe trum from the time fun tion (
F (ω)
from
exist and omputation of
F (t))
and
Fourier synthesis,
i.e., omputation of the time fun tion from its spe trum (w1 from its spe trum
3.6.
REFLECTION AND REFRACTION OF PLANE WAVES ...
Rpp (ω)F (ω)eiωz/α1 ).
form (FFT).
Su h numeri al methods are known as
91
Fast Fourier trans
Insight into the pro esses o
urring during ree tion, the topi of this hapter, 2
an be a hieved as follows: we expand (3.73) (whi h due to R0 < 1 always
onverges)and get
Rpp (ω)
∞ X n 2d 2d −iω α −iω α 2 2 R02 e = R0 1 − e n=0
−iω 4d 2d −iω α 2 3 α2 2 − R e = R0 − R0 1 − 0 1 − R0 e −iω 6d α2 − ... −R05 1 − R02 e R02
(3.78)
Substitution of this into (3.77) and taking the inverse transform of ea h element gives
w1
z z 2d 2 = R0 F t + (3.79) − R0 1 − R0 F t + − α1 α1 α2 z z 4d 6d − R03 1 − R02 F t + − R05 1 − R02 F t + − ... − − α1 α2 α1 α2
The rst term is the
ree tion from the interfa e z = 0.
expe ted, is the ree tion oe ient
R0
Its amplitude, as
of this interfa e.
The se ond term
des ribes a wave whi h is delayed by twi e the travel time through the lamella, thus, orresponding to the
ree tion from the interfa e z = d.
the expe ted size; the ree tion oe ient of this interfa e is
Its amplitude has
−R0 .
The produ t
z = 0 for waves travelling in +z− 2 and −z−dire tion, 2ρ1 α1 /(ρ2 α2 + ρ1 α1 ) and 2ρ2 α2 /(ρ1 α1 + ρ2 α2 ), is 1 − R0 . In
of the ree tion oe ients of the interfa e
multiple ree tions within the lamella (with three and ve ree tions, respe tively). The terms in (3.79) orrespond to the rays shown in Fig. 3.35.
the same way, the third and fourth term of (3.79) an be interpreted as
z=0 z =d Fig. 3.35: Ree ted and multiple ree ted rays in a lamella.
Equation (3.79) is a de omposition of the ree ted wave eld in (innite many) ray ontributions. It is fully equivalent to (3.77).
92
CHAPTER 3.
The approximation (3.74) for
BODY WAVES
Rpp (ω) orresponds to the trun ation of the exn = 0 and, therefore, the limitation on the the interfa es z = 0 and z = d, respe tively (and
pansion in (3.78) after the term for two
primary ree tions
negle ting
R02
from
relative to 1).
Exer ise 3.10 Show that for the refra tion oe ient in (3.71), it holds that
Bpp =
α1 detM αn M22
with
detM =
l 1 ρ1 . l n ρn
Apply this formula in the lamella, in ases in whi h (3.75) and (3.76) hold.
Exer ise 3.11 The
P velo ity
α2 > α1 .
of the lamella is larger then that of the surrounding medium:
Does then total ree tion o
ur? Dis uss this qualitatively.
3.7 Ree tivity method: Ree tion of spheri al waves from layered media 3.7.1 Theory The results of se tion 3.6.5 an, with relative ease, be extended to the ex itation by spheri al waves. For simpli ation of representation, we again assume that we deal, at the moment, only with
P waves in liquids.
Fig. 3.36: Explosive point sour e over liquid, layered medium.
3.7.
93
REFLECTIVITY METHOD: REFLECTION OF ...
The spheri al waves are ex ited by an explosion point sour e lo ated at hight
h
above the layered medium.
The displa ement potential of this sour e for
harmoni ex itation is ( ompare se tion 3.4)
Φ1e = with
R2 = r2 + (z + h)2 .
z axis,
Φj
in the
t− αR
1
(3.80)
Be ause of the symmetry under rotation around the
r j th layer is
ylindri al oordinates
potential
1 iω e R
and
z
are used.
The wave equation for the
1 ∂Φj 1 ∂ 2 Φj ∂ 2 Φj ∂ 2 Φj + = . + ∂r2 r ∂r ∂z 2 α2j ∂t2
(3.81)
Elementary solutions of this equation are (please he k)
J0 (kr) exp [i (ωt ± lj (z − zj ))]
with
( ompare se tion 3.6.5 for notation).
ω2 k 2 + lj2 = 2 , lj = αj J0 (kr)is
the
ω2 − k2 α2j
Bessel fun tion
! 12
(3.82)
of rst kind
and zeroth order ( ompare appendix C).
−ikj x i(ωt±lj (z−zj )) Equation (3.82) is an analogue to the solutions e ·e of the wave 2 2 2 2 2 2 2 equation ∂ Φj /∂x + ∂ Φj /∂z = (1/αj )∂ Φj /∂t dis ussed in the last hapter. In (3.82), the index
k
j
of the horizontal wavenumber
k
has been dropped, sin e
is a parameter over whi h one an integrate (furthermore, it was shown in
se tion 3.6.5 that all
kj 's
are identi al).
With (3.82), the fun tions
Z
∞
f (k)J0 (kr)ei(ωt±lj (z−zj )) dk
(3.83)
0
are also solutions of (3.81) if the integral onverges.
Thus, we ome to the
potential ansatz
Φj =
Z
0
∞
o n J0 (kr) Aj (k)ei(ωt−lj (z−zj )) + Bj (k)ei(ωt+lj (z−zj )) dk.
(3.84)
Note the lose relation of (3.84) to (3.64). Whether this ansatz a tually has a
rstly, if Φ1e se ondly, if Φj
solution, depends
from (3.80) an be represented in the integral
form (3.83) and
in (3.84) satises the boundary onditions for
z = z2 , z3 , . . . zn .
94
CHAPTER 3.
BODY WAVES
The rst requirement is satised sin e the following integral representation is
Sommerfeld integral, ompare appendix D)
valid (
1 iω e R
t− αR
1
=
Z
∞
J0 (kr)
0
k i(ωt−l1 z+h e dk. il1
(3.85)
We, therefore, an interpret the rst part
Z
∞
J0 (kr)A1 ei(ωt−l1 z) dk
(3.86)
0
of
Φ1
(with
z1 = 0)
as the in ident wave (see also se tion 3.6.5)
Φ1r ).
part is the ree ted wave
Φ1e
(the se ond
We have to ompare (3.85) and (3.86) for
lo ations in whi h the spheri al wave passes on in iden e at the interfa e
−h < z ≤ 0. A1 (k) = (k/il1 )e−il1 h . i.e., for
The
z + h = z + h
In this ase,
boundary onditions
for the interfa es an be taken from se tion 3.6.5. The
potentials (3.84) are dierentiated under the integral. from the boundary onditions is only satised
The identity following
for all r,
if the integrands are
identi al. This leads to the same system of equations for in se tion 3.6.5, i.e., (3.70). wavenumbers
lj
Aj (k)
ϕ).
k
and
and
Bj (k)
as
In ontrast to the previous se tion, the verti al
have to be onsidered now as fun tions of
angle of in iden e
z = 0,
and the omparison gives
ϕ
k
(and not of the
are onne ted via
k=
ω sin ϕ. α1
(3.87)
Rpp = B1 /A1 = −M21 /M22 has ω and angle of in iden e introdu ed via (3.87): Rpp = Rpp (ω, k).
Following se tion 3.6.5, the ree tion oe ient
been omputed as a fun tion of the angular frequen y
ϕ;;
then the dependen e on
The se ond part of
Φ1r
Φ1 ,
= =
k
an be
the ree ted wave, an then be written as
Z
∞
Z0 ∞ 0
J0 (kr)A1 (k)Rpp (ω, k)ei(ωt+l1 z) dk k J0 (kr)Rpp (ω, k)ei(ωt+l1 (z−h)) dk. il1
The orresponding verti al displa ement is
∂Φ1r = eiωt w1r (r, z, ω, t) = ∂z
Z
∞
kJ0 (kr)Rpp (ω, k)eil1 (z−h) dk
0
and the horizontal displa ement (with
J0′ (x) = −J1 (x))
(3.88)
3.7.
95
REFLECTIVITY METHOD: REFLECTION OF ...
∂Φ1r = eiωt u1r (r, z, ω, t) = ∂r
∞
Z
0
−k 2 J1 (kr)Rpp (ω, k)eil1 (z−h) dk. il1
(3.89)
The integrals in (3.88) and (3.89) are best omputed numeri ally, espe ially,
solid
in the ase of many layers. For the ree tion oe ient and
w1r
and
u1r
Rpp (ω, k)
ompressional part of the ree tion from the shear part, similar results hold, whi h only
des ribe only the
layered halfspa e
z ≥ 0.
For the
now ontain the ree tion oe ients The transition to
media, (3.88) and (3.89) also hold, but
is more ompli ated than for liquid media
Rps (ω, k).
impulse ex itation Φ1e =
1 F R
R t− α1
instead of (3.80) is relatively simple (see se tion 3.6.3). If of
F (t),
F (ω)
is the spe trum
it holds that
Φ1e =
1 2πR
Z
+∞
F (ω)e
−∞
iω t− αR
1
dω.
The orresponding displa ements of the ree ted wave are
W1r (r, z, t) U1r (r, z, t) with
w1r
=
from (3.88) and
u1r
1 2π
R +∞ −∞
F (ω)
w1r (r, z, ω, t) u1r (r, z, ω, t)
dω
(3.90)
from (3.89). The integrals in (3.88) and (3.89),
Fourier transforms of the displa ement. The numeri al omputation of (3.88), (3.89) and (3.90) is alled the Ree tivity method ; it is a pra ti al approa h for the omputation of theoreti al seismograms of body waves. With it, the amplitudes of body waves from explosions and earthquakes an be studied, thus, progressing beyond the more lassi al travel time interpretation. multiplied by
F (ω),
are, therefore, the
3.7.2 Ree tion and head waves An example for theoreti al seismograms is given in Fig. 3.37 (from K. Fu hs: The ree tion of spheri al waves from transition zones with arbitrary depthdepended elasti moduli and density. Journ. of Physi s of the Earth, vol. 16, Spe ial Issue, S. 2741, 1968). It is the result for a simple model of the rust, assumed to be homogeneous.
The point sour e and the re eivers are at the
Earth's surfa e; the inuen e of whi h has been negle ted here. The transition
96
CHAPTER 3.
Mohorovi£i¢zone,
of the rust to the upper mantle ( order
dis ontinuity,
short
BODY WAVES
Moho )
is a rst
i.e., the wave velo ities and the density hange abruptly
(for a dis ontinuity of 2nd order these parameters would still be ontinuous, but their derivative with depth would have a jump).
Fig. 3.37: Syntheti seismogram for ree tion and refra tion from a 1st order dis ontinuity (from K. Fu hs, 1968, Journ. of Physi s of the Earth).
ree tion from the Moho. For distan es from the
riti al point r∗ = 74.91 km, orresponding to the riti al angle of in iden e ϕ∗ , the rst onset is the head wave with the apparent velo ity The dominant wave is the sour e beyond the
of 8.2 km/se . Its amplitude de ays rapidly with in reasing distan e, and its ∗ form is the time integral of the ree tion for r < r . For pre riti al distan e
r,
the form of the ree tion is pra ti ally identi al to that of the in ident wave.
At the riti al point, it begins to hange its form. This was already dis ussed in se tion 3.6.3 in terms of the properties of the ree tion oe ient for plane waves (this holds for
P
and
SHwaves, respe tively).
impulse form is roughly opposite to that for sin e the ree tion oe ient to
−1
Rpp
For large distan es, the
r < r∗ .
This is also expe ted,
for the angle of in iden e
ϕ = π/2
is equal
(for liquids, this follows from the formulae given in se tion 3.5.6). The
amplitude behaviour of the ree tion is relatively similar to the trend of the absolute value
Rpp 
of the ree tion oe ient, if
Rpp 
is divided by the path
3.7.
97
REFLECTIVITY METHOD: REFLECTION OF ...
length and if one onsiders the verti al omponent (see, e.g., Fig.
3.22 and
orresponding equations). The main dis repan ies are near the riti al point. A
ording to
Rpp ,
the ree ted wave should have its maximum dire tly at the
riti al point, whereas in reality, it is shifted to larger distan es. This shift is larger, the lower the frequen y of the in ident wave.
Fig. 3.38: Ree tion amplitude versus oset as a fun tion of frequen y.
The onsideration of this shift is important when determining the riti al point from observed ree tions, e.g., in ree tion seismi s.
3.7.3 Complete seismograms Fig. 3.39 shows the
SH seismograms
potential of the ree tivity method.
This shows omplete
for a prole at the surfa e of a realisti Earth model.
The
sour e is a horizontal single for e at the Earth's surfa e, a ting perpendi ular to the prole. The dominant period is 20 se . The most pronoun ed phases are the dispersive Love waves (for surfa e waves, see hapter 4), whose amplitudes are mostly lipped. (mantle wave
S
and
The propagation paths of the largest body wave phases
SS, ore ree tion S S
and dira tion at the ore
Sdif f
are
also sket hed.) A detailed des ription of the ree tivity method is given in G. Müller: The ree tivity method: A tutorial, Journ. eophys,., vol. 58, 153174, 1985.
98
CHAPTER 3.
Fig. 3.39: Complete
SH seismograms for a
BODY WAVES
prole at the surfa e of a realisti
Earth model.
3.8 Exa t or generalised ray theory  GRT We ontinue the hapter on elasti body waves with the treatment of ree tion and refra tion of
ylindri al waves
radiated from a
line sour e
and ree ted
and refra ted at a plane interfa e whi h is parallel to the line sour e.
This
problem is more simple and less pra ti al than the ase onsidered in se tion 3.7 of a point sour e over a layered medium. On the other hand, we will learn a totally dierent way of treating wave propagation whi h leads to relatively simple
analyti al
(and not only numeri ally solvable) results. This is the main
aim of this se tion. This method, originally developed by Cagniard, de Hoop and Garvin (see, e.g., W.W. Garvin: Exa t transient solution of the buried line sour e problem, Pro . Roy. So . London, Ser. A, vol 234, pg. 528541, 1956),
an also be applied for layered media and be modied for point sour es.
In
that form it is, similar to the ree tivity method, usable for the omputation of theoreti al bodywave seismograms in the interpretation of observations.
3.8.
99
EXACT OR GENERALISED RAY THEORY  GRT
Again, we limit ourselves to treat the problem of a
liquid
model (see Fig. 3.40),
sin e we an then study the main ideas with a minimum of omputation.
Fig. 3.40: Explosive line sour e in a liquid medium.
We work with the displa ement potentials
(Φ1e =
in ident, Φ1r
=
ree ted
P − wave)
Φ1 = Φ1e + Φ1r in halfspa e 1 Φ2 in halfspa e 2. The three
and
potentials satisfy the wave equations
∇2 Φ1e,r =
1 ∂ 2 Φ1e,r , α21 ∂t2
∇2 Φ2 =
1 ∂ 2 Φ2 . α22 ∂t2
(3.91)
s2 ϕ2 , α22
(3.92)
The Lapla e transform of these equations gives
∇2 ϕ1e,r =
s2 ϕ1e,r α21
and
∇2 ϕ2 =
where ϕ1e , ϕ1r and ϕ2 are the transforms of Φ1e , Φ1r and Φ2 , respe tively, and s is the transform variable (see appendix A). We assume that the P wave starts at time t =0 at the line sour e. Therefore, the initial values of Φ1e , Φ1r and Φ2 , and their time derivatives for
t= +0, are zero outside the line sour e.
The time
derivatives have to be onsidered in the se ond derivative with respe t to
t
in
(3.91).
3.8.1 In ident ylindri al wave First, we have to study the in ident wave. Sin e the line sour e is explosive and has, therefore, ylindri al symmetry around its axis, it holds that
∇2 ϕ1e = with
R2 = x2 + (z + h)2 .
∂ 2 ϕ1e 1 ∂ϕ1e s2 + ϕ1e = ∂R2 R ∂R α21
(3.93)
The solution of (3.93), whi h an be interpreted as a
ylindri al wave in +R dire tion, is
100
CHAPTER 3.
BODY WAVES
R 1 ϕ1e = f (s) K0 ( s), s α1
(3.94)
where f (s) is the Lapla e transform of an arbitrary time fun tion K0 ( αR1 s) is one of the of zeroth order.
modied Bessel fun tions
Proof: Using the substitution x = Rα
s 1
F(t)
and
, (3.93) an be expressed as the dierential
equations of the modied Bessel fun tion
x2
dy d2 y +x − (x2 + n2 )y = 0. 2 dx dx
In the ase onsidered here, linear solutions,
K0 (x)
and
n =0.
The dierential equation has two independent
I0 (x),
respe tively. For real
x, Fig.
3.41 shows their
qualitative behaviour.
3 I0 (x)
2 1 K 0 (x) 0 1
0
Fig. 3.41: Behaviour of linear solutions
This shows that only for
x → ∞.
K0 (x)
Referen e :
(
K0 (x)
x
3
2
and
I0 (x).
is a possible solution, sin e
I0 (x)
grows innitely
M. Abramovitz and I.A. Stegun: Handbook of Mathe
mati al Fun tions, H. Deuts h, Frankfurt, 1985). Taking the inverse Lapla e transform of (3.94) in the time domain, and using the orresponden e
1 K0 s
f (s) •−◦ F (t) (F (t) ≡ 0 for t < 0) R 0 for t < R/α1 s •−◦ cos h−1 ( αR1 t ) for t > R/α1 α1
with a typi al behaviour of the solution given in Fig. 3.42;
3.8.
101
EXACT OR GENERALISED RAY THEORY  GRT
t
R/α 1 Fig. 3.42: Behaviour of the solution. the potential an be written as
Φ1e =
Z
t
R/α1
By varying
F(t),
F (t − τ ) cos h−1 (
α1 τ )dτ R
(t ≥ R/α1 ).
(3.95)
the ylindri al wave an be given dierent time dependen
ies. Equation (3.95) is the analogue to the potential
Φ1e =
1 RF
spheri al wave from an explosive point sour e. In the following, we treat the ee ts an be studied.
spe ial ase F (t) = δ(t),
t−
R α1
of a
for whi h all important
If realisti ex itations have to be treated, the results
for the potentials and displa ements of the ree ted and dira ted waves, respe tively, derived with time dependent realisti
F(t).
For
F (t) = δ(t),
have to be onvolved with
F (t) = δ(t) Φ1e = cosh−1 (
α1 t ), R
and the orresponding radial displa ement in Rdire tion is
UR =
∂Φ1e =− ∂R
t R
t2
−
Fig. 3.43: Displa ement in R dire tion.
R2 α21
1/2
(t >
R ). α1
(3.96)
102
CHAPTER 3.
BODY WAVES
The orresponding point sour e results are
Φ1e
=
UR
=
1 R H t− R α1 1 R R 1 ∂Φ1e − . δ t− = − 2H t − ∂R R α1 Rα1 α1
Fig. 3.44: Displa ement in R dire tion.
3.8.2 Wavefront approximation for UR If we write (3.96) as
UR =
R t+
R α1
−t 1/2
t−
R α1
1/2 ,
R and onsider values of t near α , we nd the approximation 1
UR ≈
−1
1/2
(2Rα1 )
·
1 t−
R α1
1/2 .
(3.97)
t is to R/α1 , and, therefore, wavefront approximation. It is more a
urate for large R, and it is, therefore, also the fareld approximation of the ylindri al wave. Within the framework of the wavefront approximation (3.97), the impulse form
This approximation is more a
urate the loser this is alled
of the ylindri al wave is independent from R, and its amplitude is proportional −1/2 to R . Both statements be ome espe ially obvious if (3.97) is onvolved with realisti ex itation fun tions integrable.
F (t).
The singularity in (3.96) and (3.97) is
3.8.
103
EXACT OR GENERALISED RAY THEORY  GRT
3.8.3 Ree tion and refra tion of the ylindri al wave The oordinates most appropriate for the study of ree tion and refra tion are the Cartesian oordinates
x
and
z.
Equation (3.92), then, takes the form
s2 ∂2ϕ ∂2ϕ + = 2 ϕ. 2 2 ∂x ∂z α 2 Appropriate elementary solutions have the form cos(kx) exp(±imz) with k + 2 2 2 m = −s /α . From these elementary solutions, more ompli ated solutions in integral form (similar to se tion 3.7.1) an be onstru ted
ϕ=
Z
∞
f (k) cos(kx)e±imz dk,
0
rstly, the potential (3.94) of the in ident se ondly, the boundary onditions for z = 0. Spe i ally, we use the
with whi h we an try to satisfy wave, and
following ansatz
ϕ1e
=
ϕ1r
=
ϕ2
m1,2
For
K0
R∞ 0
R∞ 0
R∞
=
0
s2 = −i k 2 + 2 α1,2
R α1 s
A1 (k) cos(kx)e−im1 z dk B1 (k) cos(kx)eim1 z dk A2 (k) cos(kx)e−im2 z dk
!1/2
(negative
(z > −h)
(3.98)
imaginary for positive radi ands).
in (3.94) an integral representation, similar to (3.85) for the spher
i al wave, an be found. With this
ϕ1e
(with
f(s)=1, sin e F (t) = δ(t)) it follows
that
ϕ1e
1 = s
Z
0
∞
1
k2
A omparison with
+
s2 α21
ϕ1e
1/2
1/2 # s2 2 cos(kx) exp − z + h k + 2 dk. α1 "
from (3.98) for
A1 (k) =
z>h
gives
e−im1 h . ism1
(3.99)
104
CHAPTER 3.
The boundary onditions for
BODY WAVES
z= 0 are ( ompare se tion 3.6.5)
∂ ∂Φ2 (Φ1e + Φ1r ) = , ∂z ∂z
ρ1
∂ 2 Φ2 ∂2 (Φ + Φ ) = ρ . 1e 1r 2 ∂t2 ∂t2
The Lapla e transform gives
∂ϕ2 ∂ (ϕ1e + ϕ1r ) = , ∂z ∂z
ρ1 (ϕ1e + ϕ1r ) = ρ2 ϕ2 .
From (3.98), it follows that
m1 (A1 (k) − B1 (k)) ρ1 (A1 (k) + B1 (k))
= =
m2 A2 (k) ρ2 A2 (k),
and from this
B1 (k) = Rpp (k)A1 (k)
and
A2 (k) = Bpp (k)A1 (k)
. Thus,
ρ2 m 1 − ρ1 m 2 ρ2 m 1 + ρ1 m 2
Rpp (k) = The potentials
ϕ1r
ϕ2
and
ϕ1r ϕ2
= =
and
Bpp (k) =
2ρ1 m1 . ρ2 m 1 + ρ1 m 2
(3.100)
are, therefore,
Z
∞
Z0 ∞
A1 (k)Rpp (k) cos(kx)eim1 z dk A1 (k)Bpp (k) cos(kx)e−im2 z dk.
0
w and u of the verti al and horizontal displa ement W ∂ϕ U, respe tively, an, in general, be written as w = ∂ϕ ∂z and u = ∂x
The Lapla e transforms and
and spe i ally
w1r = u1r w2 = u2
R∞ 0
R∞ 0
Rpp (k) sim1 Bpp (k) sim1
im1 cos(kx) e−im1 (h−z) dk −k sin(kx) −im2 cos(kx) e−i(m2 z+m1 h) dk. −k sin(kx)
(3.101)
3.8.
EXACT OR GENERALISED RAY THEORY  GRT
105
transformed ba k.
This is impossible
These Lapla e transforms must now be
with (A.9) in appendix A. One, rather, uses an approa h that is based in transforming (3.101) several times with fun tion theory methods until the integrals are of the form
w u The inverse
W (t) U (t)
= Z(t)
=
R∞ 0
Z(t)e−st dt.
an then be identied
(3.102)
dire tly.
positive real s, i.e., we do not onsider the whole onvergen e halfplane of the Lapla e
An important limitation has to be mentioned rst: we only onsider transform, but only the positive real axis.
This simplies the omputations
signi antly, without limitation of its generality, sin e the Lapla e transform is an analyti al fun tion.
It is, therefore, determined in the whole onvergen e
halfplane by its values on the real axis, where the integral (3.101) is real. With
cos(kx) = Re(e−ikx ) w1r = u1r w2 = u2
Re Re
R∞ 0
R∞ 0
and
sin(kx) = Re(ie−ikx ),
(3.101) an be written as
m1 e−i(kx+m1 (h−z)) dk −k −m2 Bpp (k) e−i(kx+m1 h+m2 z) dk sm1 −k
Rpp (k) sm1
The next step is a hange of the integration variables
u=
m1,2
= =
with
m1,2
u axis.
!1/2 2 1/2 s −i −s2 u2 + 2 = −is −u2 + α−2 1,2 α1,2 1/2 −s u2 − α−2 = −sa1,2 1,2 1/2
.
(3.104)
The transformed integration path is, therefore, in the sheet of the of the square root
the sheet with
The
gives
a1,2 = u2 − α−2 1,2
plane
(3.103)
ik s
so that the integration path is now along the positive imaginary transformation of the square root
.
a1,2 ,
a1,2 ≃ −u).
in whi h
a1,2 ≃ u
for
u → ∞
Riemann
holds (and not in
Introdu ing (3.104) in (3.100), gives
106
CHAPTER 3.
ρ2 a1 − ρ1 a2 ρ2 a1 + ρ1 a2
Rpp (u) =
With
k = −isu
w1r u1r w2 u2
and
dk = −isdu,
= Re
Z
and
= Re
Z
2ρ1 a1 . ρ2 a1 + ρ1 a2
(3.105)
it follows from (3.103)
+i∞
Rpp (u)
0
Bpp (u) =
BODY WAVES
+i∞
Bpp (u)
0
−i − au1
i aa21 − au1
e−s(ux−ia1 (h−z)) du
(3.106)
e−s(ux−ia1 h−ia2 z) du.
(3.107)
These expressions already have a ertain similarity with (3.102) sin e s only o
urs in the exponential term. The next step is, therefore, a new hange in the integration variable in (3.106)
t = ux − ia1 (h − z)
(3.108)
t = ux − ia1 h − ia2 z.
(3.109)
in (3.107)
From both equations,
u
has to be determined as a fun tion of
be inserted in (3.106) and (3.107), respe tively.
t
and has to
This will be dis ussed later
in more detail. At the same time, the integration path has to be transformed a
ordingly. For
u = 0,
it follows from (3.108) and (3.109), respe tively, that
t(0) = t(0) =
t0 t0
The transformed integration paths
C1
= =
h−z α1 h z α1 + α2 .
(3.110)
C2 (for u → +i∞, they
(for (3.106)) and
spe tively, start on the positive real taxis. For
(3.107)), reapproa h an
asymptote in the rst quadrant whi h passes through the entre of the oordinate system and has the slope
tan γ = x/(h + z)
tan γ = x/(h − z)
(for the ase (3.109)).
(for the ase (3.108)) and
In (3.108),
z
is always negative in
(3.109) always positive. The transformed integration paths are shown in Fig. 3.45.
3.8.
107
EXACT OR GENERALISED RAY THEORY  GRT
Fig. 3.45: Transformed integration paths.
In the ase of the ree ted wave (left in Fig. 3.45), in the ase of the refra ted wave,
C1
C1
is part of a hyperbola;
is part of a urve of higher order whi h is
similar to a hyperbola. We then get
w1r u1r
w2 u2
= Re
Z
Rpp (u(t))
(
Bpp (u(t))
(
C1
= Re
Z
C2
The last step is now to deform the paths
Cau hy's integral
C1
−i
)
du −st e dt dt
(3.111)
ia2 (u(t)) a1 (u(t)) − a1u(t) (u(t))
)
du −st e dt. dt
(3.112)
− a1u(t) (u(t))
and
C2
towards the real axis a 
ording to . The path C1,2 an now be repla ed by the path ′ C 1,2 + C1,2 in Fig. 3.46 if no poles of the integrands in (3.111) and (3.112) are lo ated inside the two paths.
Im t
C1,2 C’1,2 C1,2 t0 Fig. 3.46: Integration paths in the omplex plane.
Re t
108
CHAPTER 3.
BODY WAVES
bran h points
This is satised be ause the only singularities of a1 and a2 are the −1 and u = ±α2 , respe tively, and these bran h points are integrable ′ singularities and not poles. Finally, the ontribution of the urve C1,2 goes to
u = ±α−1 1
zero if its radius be omes innite. Thus,
w1r u1r
Z
=
"
∞ h−z α1
+
Re Rpp (u(t)) h−z α1
Z
0
w2 u2
Z
=
∞
− a1u(t) (u(t))
[0] e−st dt "
Z
+ αz 2 h z α1 + α2
)
−i
# du −st e dt dt (3.113)
Re Bpp (u(t))
h α1
+
(
(
ia2 (u(t)) a1 (u(t)) − a1u(t) (u(t))
)
#
du −st e dt dt
[0] e−st dt.
(3.114)
0
In these expressions, only the real part of the square bra kets has to be onsid−st ered sin e e is real and the integration is only over real . The addition of the
t
se ond integral with vanishing ontribution was only done for formal reasons, to allow integration over
t
from 0 to
∞ a
ording to (3.102).
Equation (3.113) and
(3.114) have, therefore, the standard form of a Lapla e transform, from whi h the original fun tion an be read
dire tly.
W1r and U1r of (h−z)/α1 . This is
The displa ements
the ree ted wave are, therefore, zero between the time 0 and not surprising sin e
(h − z)/α1
is the travel time from the sour e perpendi ular
down to the ree ting interfa e and ba k to level
z
of the sour e. This time
is, therefore, smaller, or at most equal, to the travel time of the rst ree ted onsets at this point. For
W1r U1r with
u(t)
t > (h − z)/α1
"
= Re Rpp (u(t))
h/α1 + z/α2
are zero, and for
W2 U2 u(t)
(
−i
− a1u(t) (u(t))
)
du dt
#
(3.115)
from (3.108).
Similarly, the displa ements
with
it holds that
W2 and U2 of the dira ted wave t > h/α1 + z/α2 , it holds that "
= Re Bpp (u(t))
(
ia2 (u(t)) a1 (u(t)) − a1u(t) (u(t))
)
du dt
#
for
0 ≤ t
t0
du dt . In the ase of
(from (3.110)) and the knowledge of the derivative
(3.108), this is very easy sin e here
u(t) =
x 2t R x 2t R
R = x2 + (h − z)2
− +
1/2
h−z 2 R h−z i 2 R
u
t21 − t2
an be given expli itly
21
for
1 2 2
t2 − t1
for
≤ t ≤ t1 =
R α1
(3.117)
t > t1 .
is the distan e of the sour e point from its mirror image,
i.e., from the point with the oordinates
a
ording to
h−z α1
x = 0
and
z = +h; t1 = R/α1
is,
Fermat's prin iple, the travel time of the a tual ree tion from the
interfa e. The urve of
u(t)
is given in Fig. 3.47. The derivative
omputed dire tly from (3.117). It has a singularity at
Fig. 3.47: The path of
u(t)
du/dt
in the omplex plane.
In the ase of the refra ted wave,
u(t)
has to be omputed numeri ally with a
similar urve as for the ree ted wave (Fig. 3.47).
Fig. 3.47: The path of
u(t)
an be
t = t1 .
for the refra ted wave in the omplex plane.
110
CHAPTER 3.
The numeri al omputations of a large eort.
u(t),
BODY WAVES
and its derivative, are possible without
In the following se tion, we fo us on the ree ted wave using
(3.115) and (3.117).
3.8.4 Dis ussion of ree ted wave types We assume that the
P velo ity
in the lower halfspa e is larger than that of
α2 > α1 . First, we onsider u(t1 ) = x/(Rα1 ) < α−1 2 . Sin e x/R = sin ϕ (ϕ= ∗ ∗ means that sin ϕ < α1 /α2 = sin ϕ (ϕ = riti al angle
the upper halfspa e ontaining the line sour e; re eivers
P (x, z)
for whi h
angle of in iden e), this
of in iden e), and this implies
pre riti al in iden e of the ylindri al wave.
Fig. 3.49: Sket h for line sour e and its mirror point.
In this ase,
−1 u(t) for t < t1 is smaller than α−1 2 and, thus, even smaller than α1 .
a1 (u(t)) and a2 (u(t)), a
ording to (3.104), are positive imaginary, Rpp (u(t)), a
ording to (3.105), is real. Sin e du/dt is real for the times
Therefore, and
onsidered, it follows with (3.115) that the real part of the square bra kets for all
t < t1
is zero. For
t > t1 , u(t)
be omes omplex, and the real part is nonzero.
Not surprisingly, the displa ement, therefore, starts at If
α1 < α2 ,
t = t1 .
this argument holds for arbitrary re eiver lo ations in the upper
halfspa e. If
u(t1 ) = x/(Rα1 ) > α−1 2 ,
the
head wave
angle of in iden e ϕ is larger than the angle ϕ∗ ,
as the rst onset. In this ase, a2 (u(t)) be omes u(t2 ) = α−1 2 . The same does not hold for a1 (u(t)). Thus, Rpp (u(t)) has nonzero real and imaginary parts for t > t2 , and, −1 in (3.108), it therefore, (3.115) is already nonzero for t > t2 . Putting u = α2 follows that and we expe t a
real at the time t2 whi h is dened via
t2 =
x −2 21 + (h − z)(α−2 1 − α2 ) < t1 . α2
3.8.
111
EXACT OR GENERALISED RAY THEORY  GRT
This is the arrival time of the head wave as expe ted a
ording to Fermat's prin iple for the ray path from Q to P(x,z) in Fig. 3.50.
Fig. 3.50: Path of head wave from sour e Q to re eiver P.
Considering this ase at time
t = t1 , we expe t,
due to the sudden hange in the
u(t) propagates, signi ant hanges in the displa ement (3.115), ree tion proper gives a signi ant signal, and this is something
urve on whi h i.e., that the
whi h indeed an be observed. Our derivation has shown that the head wave, and also the ree tion, an be derived from the potential
Φ1r , i.e., no separate des ription was ne essary for the
head wave. If we had studied solid media, we, possibly, ould have identied
P
a se ond arrival whi h is an additional interfa e or boundary wave (
to
S
onversion) ( ompare, e.g., the work by Garvin quoted earlier). In se tion 3.7 (see head wave in Fig.
3.37), we en ountered a similar situation in that the
head wave was in luded in the solution. Both methods (ree tivity method in se tion 3.7 and GRT this se tion) give a omplete solution unless simpli ations for numeri al reasons are introdu ed.
Fig. 3.51: Sket h of the displa ement (in horizontal and verti al dire tion) for pre riti al (left) and post riti al (right) in iden e, respe tively.
112
For
CHAPTER 3.
t → ∞,
BODY WAVES
the limit of the displa ement is nonzero, as for the ase of the in
ident wave ( ompare (3.96)). The singularities at
t = t1
are always integrable.
Therefore, onvolution with a realisti ex itation fun tion
F(t)
is always possi
ble. On the right side of Fig. 3.51, the displa ement starts before the ree tion proper arriving at
t1 .
Fermat's prin iple, therefore, does not give exa tly the
arrival time of the rst onset in the ase where the ree tion is not
arrival.
the rst
3.9.
113
RAY SEISMICS IN CONTINUOUS INHOMOGENEOUS MEDIA
Fig. 3.52: Theoreti al seismograms (verti al displa ement) for a rustal model. α1 = 6, 4km/s, α2 = 8, 2 km/s, ρ1 = 3, 0 g/ m3 , ρ2 = 3.3 g/ m3 ,
Parameters: h
= 30
km, z
= −30
km.
On the left, the exa t seismograms (GRT) and
on the right a wavefront approximation is shown.
From G. Müller:
Exa t
ray theory and its appli ation to the ree tion of elasti waves from verti ally inhomogeneous media, Geophys. Journ. R.A.S. 21, S. 261283, 1970.
For realisti
F(t),
this dis repan y is usually small. A full example using the
theory of this se tion, is given in Fig. 3.50b. For a des ription of the dierent wave types
et ., see also se tion 3.7.2.
In the ase of a layered medium with
more than one interfa e, the wave eld an
be broken into separate ray ontributions, as was done for the lamella in se tion 3.6.5. For ea h ray ontribution, a formula of the type of (3.115) or (3.116) an be given, whi h an ontain head or boundary wave ontributions. This is the reason for the name "exa t or generalised ray theory". Another ommon name is the "Cagniardde Hoopmethod".
Exer ise 3.12: Give wavefront approximations for the ree ted wave and head wave, i.e., expand (3.115) around the arrival times wave), respe tively.
t1
(of the ree tion) and
Distinguish between
slowly
(of the head
t = t1 and t = t2 , respe tively, and t − t1 , t1 − t and t − t2 , respe tively.
an be repla ed by their values for varying terms, whi h depend on
t2
varying ontributions, whi h
rapidly
3.9 Ray seismi s in ontinuous inhomogeneous media With the ree tivity method and the GRT, we have dis ussed waveseismi methods, whi h if applied in
ontinuous inhomogeneous
media (in our ase ver
ti ally inhomogeneous media), require a segmentation in homogeneous regions (in our ase homogeneous layers).
Waveseismi methods for ontinuous in
homogeneous media, without this simplied representation of real media, are often more ompli ated ( ompare example in se tion 3.10). In this se tion, we will now sket h the
rayseismi (or rayopti al) approximation
of the wave the
ory in inhomogeneous media. We will also show that it is the
approximation
high frequen y
of the equation of motion (2.20) for the inhomogeneous elasti
ontinuum. We restri t our dis ussion again to a simplied ase, namely the propagation of
SH waves
in a twodimensional inhomogeneous medium.
velo ity
β
and shear modulus
µ
displa ement omponent is in the
depend
The
y dire tion where density ρ, S only on x and z. The only nonzero
sour e is assumed to be a linesour e in the
y dire tion v = v(x, z, t).
114
CHAPTER 3.
BODY WAVES
3.9.1 Fermat's prin iple and the ray equation Fermat's prin iple states that the travel time of the (SH )wave from the sour e Q to an arbitrary re eiver P along the seismi ray is an extremum and, therefore, stationary, i.e., along ea h innitesimally adja ent path between P and Q (dashed in Fig. 3.53) the travel time is either larger or smaller.
Fig. 3.53: Ray with extremum path and innitesimally adja ent ray.
In most ases, the travel time along a seismi ray is a minimum, but there are also ases, where it is a maximum (e.g., the body waves PP, SS, PKKP). If we des ribe an arbitrary path from
{x = x(p), z = z(p)},
ds = We onsider now
many
p = p1
p = p2
at
Q
and
P
to
Q
via a
the element of the ar length
"
dx dp
2
+
Q P, respe tively.
su h paths from at
dz dp
2 # 21
to
s
parameter representation
an be written as
dp.
(3.118)
P. They all have the same value p annot be identi al to
Therefore,
s ; p ould, for example, be the angle between the line onne ting the oordinate x or zaxis. The seismi ray is the
entre to the point along the way and the path for whi h
T =
Z
p2
β
−1
(x(p), z(p))
p1
"
dx dp
2
+
dz dp
2 # 12
dp =
Z
p2
F
p1
dx dz x, z, , dp dp dp
dz dx = x′ , = z′ dp dp is an extremum.
The determination of the seismi ray has, therefore, been
redu ed to a problem of
equations
al ulus of variations.
∂F d ∂F =0 − ∂x dp ∂x′
and
This leads to the
∂F d ∂F = 0. − ∂z dp ∂z ′
EulerLagrange
3.9.
RAY SEISMICS IN CONTINUOUS INHOMOGENEOUS MEDIA
115
This gives then, for example,
1 ∂ x +z 2 ∂x ′2
Division by
′2
(x′2 + z ′2 )1/2 ,
square bra ket with
1 β(x, z)
d 1 x′ = 0. − dp β(x, z) (x′2 + z ′2 ) 12
multipli ation of nominator and denominator of the
dp, and use of (3.118) gives d ds
1 dx β ds
∂ = ∂x
1 . β
(3.119)
d ds
1 dz β ds
∂ = ∂z
1 . β
(3.120)
Similarly,
Equations (3.119) and (3.120) are the the parameter representation of the ray. With
dierential equations of the seismi ray in
{x = x(s), z = z(s)} where s
is now the ar length
dz dx = sin ϕ, = cos ϕ ds ds (ϕ= angle of the ray versus the
z dire tion), it follows that
sin ϕ β cos ϕ d ds β d ds
=
∂ ∂x
=
∂ ∂z
1 β . 1
(3.121)
β
x
1 0
Q
ϕ z
dx ds dz
11 00 00 11 00 11
P
11 00
Fig. 3.54: Ray in xz oordinate system.
These two equations an now be onverted into another form of the ray equation (show)
116
CHAPTER 3.
1 ∂ϕ = ds β This dierential equation for the ray path. The inverse of
sin ϕ
ϕ(s) is dϕ/ds
∂β ∂β − cos ϕ ∂z ∂x
BODY WAVES
.
(3.122)
well suited for numeri al omputations of
−1 ∂β ∂β ds = β sin ϕ − cos ϕ r= dϕ ∂z ∂x is the
radius of the urvature of the ray.
(3.123)
The ray is urved strongly (r is smaller)
where the velo ities hange strongly (∇β large).
Spe ial ases a)
β=
onst.
From (3.122), it follows that
b)
β = β(z)
dϕ/ds = 0.
(no dependen e on
The ray is straight.
x)
The rst equation in (3.121) gives, after integration,
sin ϕ = q = const β
Snell's law)
along the whole ray (
sour es and re eivers at the level apex
S. The ray parameter q
velo ities
β(0)
where
z=0
q
is the
(3.124)
ray parameter
is onne ted to the takeo angle
at the sour e and the turning point
q=
ϕ
zs
z Fig. 3.55: Ray in 1D medium.
β(zs ),
00000 11111 P 00 11 00 11 11 00 S 1 0
ϕ0 ,
the seismi
respe tively, via
sin ϕ0 1 = . β(0) β(zs )
11111∆ Q 00000 z=0 0 1
of the ray. For
the ray is symmetri with respe t to its
x
3.9.
RAY SEISMICS IN CONTINUOUS INHOMOGENEOUS MEDIA
zs
A turning point depth distan e
∆
117
β(z) < β(zs ) for all z < zs . For z = 0, it holds (using (3.124))
is only possible if
in whi h the ray reappears at level
that
∆(q) = 2
Z
S
dx = 2
Z
zs
Z
tan ϕdz = 2q
0
Q
zs 0
β −2 (z) − q 2
21
dz.
(3.125)
The ray's travel time is
T (q) = 2
Z
S
Q
ds =2 β
Z
0
zs
dz =2 β cos ϕ
Z
0
zs
− 1 β −2 (z) β −2 (z) − q 2 2 dz.
Equations (3.126) and (3.127) are parameter representations of the
urve of the model.
(3.126)
travel time
An example for ray paths and travel time urves in a model
with a transition zone is given in Fig. 3.56.
Fig. 3.56: Ray paths and travel time urves in a model with a transition zone.
The
slope
of the travel time urve is (show)
dT = q. d∆
)
β = a + bx + cz
(linear dependen e from
In this ase, it follows from
(3.127)
x
and
z)
118
CHAPTER 3.
BODY WAVES
c sin ϕ − b cos ϕ dϕ = ds β and by dierentiation with respe t to
s (ϕ and β are fun tions of s
x
via
and
z)
c cos ϕ + b sin ϕ dϕ c sin ϕ − b cos ϕ d2 ϕ = (b sin ϕ + c cos ϕ) = 0. − ds2 β ds β2 This implies
dϕ/ds
is a onstant and, therefore, the urvature radius along the
The ray is, therefore, a ir le, or a se tion of it.
whole ray.
from (3.123) if
β
Its radius
r
ϕ0 .
to the takeo angle
follows
ϕ
equal
M, the entre of the ir le, an be found from
as
is hosen identi al to the value at the sour e Q and
Q
shown in Fig. 3.57.
Fig. 3.57: Ray paths in a model with a linear velo ity law in x and z. The travel time from
T =
Z
P
Q
with An
ds = β
Z
Q ϕ1
ϕ0
to
P
is
dϕ = c2 + b c sin ϕ − b cos ϕ
δ = arctan(b/c).
arbitrary
velo ity law
β(x, z),
1 2 −2
i h tan (ϕ12−δ) h i ln tan (ϕ0 −δ) 2
given at dis rete points
(xn , zn )
in a model,
an linearly be interpolated pie ewise between three neighbouring points. Then the laws derived here an be applied. The ray then onsists of several se tions of
ir les, and at the transition between two regions with dierent linear velo ity laws, the tangent to the ray is ontinuous.
A orresponding travel time and
plotting program is, for many pra ti al appli ations, already su ient.
Exer ise 3.13 a) Show that point M is on the line
β = 0.
b) Derive the formula for T given above.
3.9.
RAY SEISMICS IN CONTINUOUS INHOMOGENEOUS MEDIA
119
3.9.2 High frequen y approximation of the equation of motion The equation of motion for inhomogeneous, isotropi media (2.20) an for
SH 
wave propagation in twodimensional media without volume for es, be simplied to
∂ 2v ∂ = ∂t2 ∂x
ρ
∂v ∂ ∂v µ + µ . ∂x ∂z ∂z
(3.128)
For the time harmoni ase we use the ansatz
v(x, z, t) = A(x, z) exp [iω (t − T (x, z))] . This is an ansatz for expe t that amplitude
(3.129)
high frequen ies sin e only for su h frequen ies an we A(x,z) and travel time fun tion T(x,z) to be frequen y
independent in inhomogeneous media. In homogeneous media far from interfa es, this is true for all frequen ies as long as one is a few wavelengths away from the sour e. Using (3.129) in (3.128), and sorting with respe t to powers of
ω,
it follows that
2
( "
ω A µ
+iωA
∂T ∂x
2
+
∂T ∂z
2 #
−ρ
)
2 ∂µ ∂T ∂ T ∂2T ∂µ ∂T ∂ ln A ∂T ∂ ln A ∂T + + +µ + + 2µ ∂x ∂x ∂z ∂z ∂x2 ∂z 2 ∂x ∂x ∂z ∂z 2 ∂µ ∂A ∂µ ∂A ∂ A ∂2A (3.130) − = 0. + + +µ ∂x ∂x ∂z ∂z ∂x2 ∂z 2
For su iently high frequen ies, the three terms of this equation are of
dierent
magnitudes. To satisfy (3.130), ea h term has then to be zero independently,
espe ially the rst two terms
2 2 ∂T 1 ∂T + = ∂x ∂z β2 ∂ ln A ∂T ∂µ ∂T ∂ ln A ∂T ∂µ ∂T 2µ = − + − − µ∇2 T. ∂x ∂x ∂z ∂z ∂x ∂x ∂z ∂z
(3.131)
(3.132)
Eikonal equation, and it ontains on the right side the S velo ity β = (µ/ρ)1/2 . After solving the Eikonal, T is transport equation (3.132), and ln A and A are determined.
Equation (3.131) is the lo ationdependent inserted in in the
This, in prin iple, solves the problem. The frequen yindependent third term
120
CHAPTER 3.
in (3.130) will usually not be zero with
A
from (3.132).
BODY WAVES
Solution (3.129) is,
therefore, not exa t, but be omes more a
urate, the higher the frequen y. We still have to derive the onditions under whi h the rst two terms of (3.130) are indeed of dierent order and, therefore, the separation into (3.131) and (3.132) is valid.
ea h single ea h single sum
The ondition follows from the requirement that
summand in the se ond term has to be small with respe t to mand in the rst term, for example,
2 ∂µ ∂T ω ≪ ω 2 µ ∂T . ∂x ∂x ∂x
From (3.131), it follows roughly
With
µ = ρβ 2 ,
it follows
∂T /∂x = 1/β .
Thus,
β ∂µ µ ∂x ≪ ω. β ∂ρ ∂β ≪ ω. + 2 ρ ∂x ∂x
(3.133)
Similar relations follow from the other summands in (3.130). Usually, the required onditions are formulated as follows: the high frequen y approximations (3.131) and (3.132) are valid for frequen ies whi h are large with respe t to the
velo ity gradients
ω ≫ ∇β =
"
∂β ∂x
2
+
∂β ∂z
2 # 21
.
(3.134)
Equation (3.133) shows also, that density gradients have also an inuen e. Equation (3.134) an be expressed even more physi ally: the relative hange of the velo ity over the distan e of a wavelength has to be smaller then
2π
(show).
Example We solve (3.131) and (3.132) in the simplest ase of a plane in zdire tion with the assumption that
ρ, β
and
µ
SH wave propagating z. Ansatz
depend only on
(3.129) then simplies to
v(z, t) = A(z) exp [iω (t − T (z))] . The solution of (3.131) is
(3.135)
3.9.
RAY SEISMICS IN CONTINUOUS INHOMOGENEOUS MEDIA
1 dT = , T (z) = dz β(z) where
T (z)
is the
Z
z
0
121
dζ β(ζ)
Swave travel time, with respe t to the referen e level z = 0.
Equation (3.132) an be written as
1 dµ µ dβ 2µ d ln A =− + 2 . β dz β dz β dz Thus,
−1 d ln (ρβ) 2 d ln β d ln µ = = − dz dz dz 12 ρ(0)β(0) . = A(0) ρ(z)β(z)
d ln A dz
1 2
A(z) The amplitudes of the impedan e
ρβ .
SH wave
vary, therefore, inversely proportional to the
The nal solution of (3.135) is
ρ(0)β(0) v(z, t) = A(0) ρ(z)β(z)
21
Z exp iω t −
0
z
dζ β(ζ)
.
From these results we on lude that, in the ase onsidered, an frequent
SH wave propagates without hanging its form.
(3.136)
impulsive, high
Exer ise 3.14 Vary the velo ity
z,
but via a
exa tly via
β
and the density
ρ
not
step somewhere in between. the SH refra tion oe ient
ontinuously
from depth 0 to depth
Then the amplitudes an be derived (3.40).
Show that (3.136) gives the
same results if the relative hange in impedan e is small with respe t to 1. Hint: Expansion in
both
ases.
3.9.3 Eikonal equation and seismi rays From (3.129), it follows that surfa es of onstant phase are given by
t − T (x, z) = const. In the impulse ase, these surfa es are the wavefronts separating perturbed and unperturbed regions. This is why the term wavefront is used also here. Complete dierentiation with respe t to
t
gives
122
CHAPTER 3.
BODY WAVES
− → dx ∂T dx ∂T dz + = ∇T · = 1, ∂x dt ∂z dt dt where
− → dx/dt = (dx/dt, 0, dz/dt)
(3.137)
is the propagation velo ity of the wavefront.
− → dx/dt and the ve tor ∇ T , whi h is perpendi ular to the wavefront, are parallel, sin e a
ording → − to the Eikonal equation (3.131) ∇ T  = 1/β , it holds that dx/dt = β . This The obvious interpretation of (3.137) for isotropi media is that
means that the wavefronts propagate perpendi ular to themselves with the lo al velo ity
The
β.
orthogonal traje tories of the wave are dened as seismi rays.
We still have
to show that they are the rays dened via the Fermat's prin iple. We demonstrate this by showing that the dierential equations of the seismi ray, (3.119) and (3.120), also follow from the Eikonal equation. As before, we des ribe the ray via its parameter representation
s.
Ve tor
− → dx/ds = (dx/ds, 0, dz/ds)
{x = x(s), z = z(s)}
with the ar length
is a unit ve tor in ray dire tion for whi h,
using the statements above, we an write
− → dx = β∇ T = β(∂T /∂x, 0, ∂T /∂z). ds
(3.138)
Instead of (3.119), we, therefore, have
d ds
1 dx β ds
=
d ds
∂T ∂x
.
Using (3.138) and the Eikonal equation on the right side, we derive
d ds
1 dx β ds
= = =
2 ∂ 2 T dx ∂ 2 T dz ∂ 2 T ∂T ∂ T ∂T + = β + ∂x2 ds ∂x∂z ds ∂x2 ∂x ∂x∂z ∂z " 2 2 # β ∂ ∂T 1 ∂T β ∂ = + 2 ∂x ∂x ∂z 2 ∂x β 2 β (−2) ∂β 1 1 ∂β ∂ . =− 2 = 2 β 3 ∂x β ∂x ∂x β
That is identi al to (3.119) and a similar derivation holds for (3.120). following is the
ray equation in ve tor form,
The
whi h is also valid in the three
dimensional ase
d ds
− →! 1 1 dx =∇ . β ds β
(3.139)
3.9.
RAY SEISMICS IN CONTINUOUS INHOMOGENEOUS MEDIA
This shows that
mi s.
123
ray seismi s is a high frequen y approximation of wave seiskinemati aspe ts of wave
We have, until now, limited the dis ussion on
propagation, i.e., on the dis ussion of wave paths, travel times and phases.
Dynami
parameters, espe ially amplitudes, were not dis ussed ex ept in the
simple example in se tion 3.9.2. The following se tion gives more details on this aspe t. Before doing this, we give the form of (3.139) whi h is often used in numeri al
al ulations, espe ially in three dimensions. equation of 2nd order for 1st order for
− → x
and the
the absolute value
1 β)
− → x
The
single
ordinary dierential
(3.139) is repla ed by a system of
→ p = slowness ve tor −
− →
two
equations of
1 dx β ds (ve tor in ray dire tion with
− → − → dx dp 1 → = β− p, =∇ . ds ds β Ee tive numeri al methods for the solution of systems of ordinary dierential equations of 1st order exist, e.g., the RungeKuttamethod.
3.9.4 Amplitudes in ray seismi approximation Within the framework of ray seismi s developed from Fermat's prin iple, amplitudes are usually omputed using the assumption that the energy radiated into a small ray bundle, remains in that bundle. This assumption implies that no energy exits the bundle sideways via dira tion or s attering and no energy is ree ted or s attered ba kwards. This is only valid for high frequen ies. In the following, we derive a formula whi h des ribes the hange of the displa ement amplitude along a ray radiated from a line sour e in a two dimensional inhomogeneous medium. The medium shall have no dis ontinuities. We onsider (see Fig. 3.58) a ray bundle emanating from a line sour e a width of the point
dl(M ) at the
P.
referen e point
M
lose to
1 0 Q 0 1 dl(M)
11 00 00 11 M
Q
and a width of
x z
dl(P)
P 0 1
Fig. 3.58: Ray bundle emanating from the sour e Q.
Q
with
dl(P ) near
124
CHAPTER 3.
Equation (3.129) holds for the displa ement in
M
and
BODY WAVES
P, or in real form
v = A sin [ω(t − T ] . Our aim is the determination of the amplitude ration
A(P )/A(M ).
We rst de
termine the energy density of the wave, i.e., the sum of kineti and potential en1 1 2 2 2 2 ergy per unit volume. The kineti energy density is ρv 2 ˙ = 2 ρω A cos [ω(t − T )]. Averaged over the period 2π/ω , the potential and kinemati energy density have 1 1 2 2 2 the same value 4 ρω A sin e the average of cos x is identi al to 2 . Then the energy density averaged over a period an be written as
∆E 1 = ρω 2 A2 . ∆V 2 Consider a ube with the volume
∆V = dl dy ds; its ross se tion dl dy is perds is exa tly 1 wavelength β2π/ω .
pendi ular to the ray bundle and its length The energy
∆E =
1 2 2 ρω A ∆V = πωρβA2 dl dy 2
ontained within this ube ows per period through the ross se tion
∆E
the ray bundle. Sin e no energy leaves the bundle,
M. From this the amplitude ratio
at
P
dl dy
of
is the same as at
follows as
1 A(P ) ρ(M )β(M )dl(M ) 2 . = A(M ) ρ(P )β(P )dl(P ) As in (3.136), impedan e hanges o
ur along the ray.
(3.140)
The square root of
the hange in the ross se tion is the important parameter for the amplitude variation. In the most general threedimensional ase, the ross se tion
surfa e
of the
threedimensional
dl
has to be repla ed by
ray bundle.
Equation (3.140) an be approximated by tra ing su iently many rays through the medium using the methods dis ussed previously, and then determining their perpendi ular distan es (or rossse tion surfa es in three dimensions). These methods are based mainly on the solution of the ray equation (3.139), whi h is also alled the equation of the
kinemati ray tra ing.
A more stringent approa h
to al ulate (3.140) is based on dierential equations whi h are dire tly valid for the ross se tion of a ray bundle; they are alled equations of the
ray tra ing.
dynami
Their derivation annot be treated here; for details, see the book
of ervený, Molotov and P²en ík (1977). For point
P
on a horizontal prole, e.g., at
z = 0,
a loser look at (3.140) and
Fig. 3.59 helps in understanding the physi al meaning.
3.9.
RAY SEISMICS IN CONTINUOUS INHOMOGENEOUS MEDIA
d∆
∆
P 0 1 0 ϕ1
11 00 Q 00 11
a
1M 0
ϕ0
125
x,∆ d∆ dl(P)=cos ϕ d ϕ d ϕ 0 0
dl(M)=ad ϕ 0
z Fig. 3.59: Ray paths in an inhomogeneous medium.
The horizontal distan e of
Q
to
dl(P )
P.
P
is
∆(ϕ0 ) with the takeo angle ϕ0
The distan e of the referen e point
M
from
Q
is
a.
of the ray from
With
dl(M )
and
from Fig. 3.59, it follows from (3.140) that
21 1 A(P ) a ρ(M )β(M ) 2 . · = d∆ A(M ) ρ(P )β(P ) cos ϕ dϕ 0
(3.141)
P is on, or lose to, the neighturning points of the travel time urve on the horizontal prole
This expression shows that problems o
ur if bourhood of
( ompare, e.g., the travel time urve in Fig. 3.56). At these points,
d∆/dϕ0
hanges its sign, and that an happen either with a ontinuous or non ontinuous pass through a zero. In the rst ase, innite amplitudes o
ur in ond ase, the amplitudes be ome non ontinuous. and nonphysi al.
P ; in the se 
Both ases are unrealisti
Equations (3.141) and (3.140), respe tively, an, therefore,
only be used at some distan e from the turning points ( austi s) of the travel time urves. Unfortunately, this means that the points with some of the largest amplitudes annot be treated properly under these assumptions; more sophisti ated methods (like the WKBJ method, to be dis ussed later, or the Gaussian Beam method) must be employed. Despite this disadvantage, the formulae given above (and their orresponding equations in three dimensions) are very useful in seismologi al appli ations. They an be easily extended to in lude refra tions and ree tions at dis ontinuities.
This requires the determination of hanges in the ross se tion of
the ray bundle at dis ontinuities, and the in lusion of ree tion and refra tion
oe ients. A further problem of the energy ansatz used in this se tion is that it only gives the amplitudes of a seismi ray but no information on its phase hanges that o
ur in
addition
always su ient to add
to the phase hanges in the travel time term. It is not
exp [iω(t − T )] to (3.140) and (3.141), respe tively, e.g.,
on retrograde travel time bran hes in the ase that the velo ities are only a
126
CHAPTER 3.
fun tion of
z.
BODY WAVES
The WKBJ method and the Gaussian Beam method solve this
problem by tra king an additional parameter, the KMAH index named after Keller, Maslov, Arnold and Hormander, whi h ounts the austi s en ountered along the ray. In the following se tion, the WKBJ theory for verti al inhomogeneous media, whi h avoids some of the ray theory problems dis ussed, is presented; it ontains more wave seismi elements.
Exer ise 3.15 Use the ray parameter
q instead of the takeo angle ϕ0 in the amplitude formula
(3.141) in the ase of a verti ally inhomogeneous medium and then use (3.127). What is the relation between the amplitudes and the travel time urve
T (∆)?
WKBJ method
3.10
Now we will onsider total ree tion at a verti ally inhomogeneous medium using the WKBJ method.
3.10.1 Harmoni ex itation and ree tion oe ient We onsider a medium whose velo ity for
z > 0
β(z)
in the medium. A
for
z≤0
is
β(0),
i.e., onstant and
ontinuous fun tion of z, i.e., no dis ontinuities exist plane SH wave may propagate obliquely in the lower half
an be any
z < 0 with the horizontal wavenumber k = ω sin ϕ/β(0) (ray parameter q = sin ϕ/β(0) (ϕ = angle of in iden e). Note the dieren e to the example in
spa e
hapter 3.9.2 with verti al propagation. Then the ray seismi s of the verti ally inhomogeneous medium, se tion 3.9.1, suggests inserting
∂T /∂x = q = constant
in the Eikonal equation (3.131) ( om
pare with (3.127)). This then, gives
T (x, z) = qx +
Z
0
z
−2 1 β (ζ) − q 2 2 dζ,
travel time of the Swavefront from the interse tion of the origin x,z ). The rest of the dis ussion is as in se tion 3.9.2., and leads to
and, thus, the to the point (
v0 (x, z, t) = s(ζ)
=
1 Z z µ(0)s(0) 2 A(0) exp iω t − qx − s(ζ)dζ (3.142) µ(z)s(z) 0 1 −2 β (ζ) − q 2 2 . (3.143)
3.10.
q
127
WKBJ METHOD
is the horizontal and
s
is identi al to (3.136).
S wave.
is the verti al
slowness
Equation (3.142) is the
q= 0, (3.142) WKBJ approximation of the
of the wave. For
It is a useful high frequen y approximation, as long as
β −1 (z) > q ,
i.e., as long as the seismi ray whi h an be asso iated with the wave
propagating horizontally.
is not
If the velo ity, e.g., with in reasing depth de reases or −1 , (3.142) is appli able . if it in reases, but does not rea h the value q
for all z
For ases of interest and a medium with in reasing velo ities for in reasing −1 depth, a depth zs is rea hed where β(zs ) = q . At this depth, where the ray propagates horizontally, (3.142) diverges. Equations (3.129), (3.131) and (3.132) are
insu ient
for the des ription of the waveeld near the turning point of
rays. If (3.142) is onsidered for
z > zs
with
β(z) > β(zs ),
i.e., the velo ity
ontinues to in rease, a stable result an, again, be obtained. The integral in the exponential term from
zs
to
z
is imaginary, thus, giving an exponential
z SH wave de reases as expe ted.
de ay of the amplitudes with in reasing , i.e., below the ray's turning point the amplitude of the
For
z < zs ,
the waveeld is
insu iently des ribed by (3.142) sin e (3.142) represents only the downward
propagating in ident SH wave. A similar equation an be given for the ree ted SH wave upward propagating from the turning point
µ(0)s(0) v1 (x, z, t) = RA(0) µ(z)s(z)
21
Z exp iω t − qx +
z
s(ζ)dζ 0
.
(3.144)
That this wave propagates upwards an be seen from the positive sign before the integral in the exponent.
R is, as an be seen by the sele tion z=0
(3.144), the amplitude ratio wave; in other words, the
v1 (x, 0, t)/v0 (x, 0, t)
in (3.142) and
of the ree ted to the in ident
ree tion oe ient of the inhomogeneous halfspa e
is z>0.
Its determination requires a quantitative onne tion of the whole eld
v0 + v1 z > zs .
To tie these two solutions together is, as mentioned before, not possible
for
z < zs
with the already mentioned exponentially de aying eld for
with the high frequen y approximation of the equation of motion used until now. The required onne tion be omes possible with another high frequen y approximation of (3.128), namely a wave equation with depth dependent velo ity
∇2 V v
= =
∂2V ∂x2 1 1
2
+ ∂∂zV2 = V
1 ∂2V β 2 (z) ∂t2
µ 2 (z)
)
.
(3.145)
This high frequen y approximation is valid under ondition (3.134), as an be shown by inserting in (3.128). For plane waves, it follows from the ansatz
V (x, z, t) = B(z) exp [iω(t − qx)]
128
CHAPTER 3.
via (3.145) an ordinary dierential equation for
BODY WAVES
B(z)
B ′′ (z) + ω 2 β −2 (z) − q 2 B(z) = 0.
This equation has now to be solved for large ω . ω 2 f (z)B(z) = 0 in the neighbourhood of a zero of
alled For
(3.146)
The solutions of
B ′′ (z) +
f(z) and for large ω is generally WKBJ solution after the authors  Wentzel, Kramers, Brillouin, Jereys.
z < zs ,
the previously dis ussed superposition of (3.142) and (3.144) of the
in ident and ree ted wave of
v
results. For
z
z > zs
the exponentially de aying
is in the immediate neighbourhood of zs 2 2 has to be examined in more detail. We approximate the oe ient ω s (z) 2 of B(z) (with from (3.143)) and get, with s (zs ) = 0, β(zs ) = q −1 and β ′ (zs ) > 0,
solution follows.
The ase that
s(z)
linearly
B ′′ (z) − 2ω 2 q 3 β ′ (zs )(z − zs )B(z) = 0.
(3.147)
This equation an, with the substitution,
1 y(z) = 2ω 2 q 2 β ′ (zs ) 3 (z − zs )
be transformed into the dierential equation of the
(3.148)
Airy fun tions
C ′′ (y) − yC(y) = 0. The solution of interest to us,
C(y) = Ai(y),
is dis ussed in appendix E (more
on Airy fun tions an be found in M. Abramovitz and I.A. Stegun: Handbook of Mathemati al Fun tions, H. Deuts h, Frankfurt, 1985). The depths
z < zs
Ai(y). From Fig. E.2, it follows that the transition from the os illatory solution B(z) of (3.147) with z < zs to the exponentially damped solution for z > zs is without singularity. This then, also holds for the displa ement v, in ontrast to what one would (z
> zs )
orrespond to arguments
y < 0 (y > 0)
of
expe t from (3.142) and (3.144). The os illatory behaviour of
v0
and the ree tion
v1 ,
B(z)
for
z < zs
indi ates that the in ident wave
build a standing wave with nodes of the displa ement
at depths whi h orrespond to the zeros of the Airy fun tion. The ree tion
oe ient R in (3.144) is now determined in su h a way, that the superposition of (3.142) and (3.144), in the term that depends on z, is identi al to the Airy fun tion. Due to the high frequen y assumption, the asymptoti form of Ai(y) for large negative y an be used
Ai(y) ≃ π Furthermore, for
z < zs
− 12
− 41
y
sin
π 2 23 y + 3 4
.
(3.149)
3.10.
129
WKBJ METHOD
v0 + v1
h i1 Rz µ(0)s(0) 2 = A(0) µ(z)s(z) exp iω t − qx − 0 s s dζ Rz Rz Rz · exp iω z s s dζ + R exp 2iω 0 s s dζ · exp −iω z s s dζ .
(3.150)
The
z dependen e
v0 + v1
of
is given by the urved bra ket.
It will now 2 2
be determined in approximation. With the approximation (3.147) 2ω 2 q 3 β ′ (zs )(zs − ζ), it follows
ω
with
y = y(z)
Z
zs
s dζ z
1 = ± 2ω 2 q 3 β ′ (zs ) 2
Z
zs
z
ω s (ζ) =
1
(zs − ζ) 2 dζ
1 2 3 = ± 2ω 2 q 3 β ′ (zs ) 2 (z − zs ) 2 3 2 3 = ± y 2 3 = ±Y
from (3.148).
The positive (negative) sign holds for positive
(negative) frequen ies. Thus, for the urved bra kets in (3.150)
{· · ·} = e±iY + Ze∓iY with the abbreviation
Z Z = R exp 2iω
0
With
Z = ±i
for
ω
{· · ·} = =
zs
s dζ .
(3.151)
0), it follows that
1 π (1 ± i)(cos Y + sin Y ) = 2 2 (1 ± i) sin Y + 4 1 π 2 23 y + 2 2 (1 ± i) sin , 3 4
and, therefore, the required agreement with the main term in (3.149). in (3.151) gives then the
Rz ω exp −2iω 0 s s(ζ)dζ = i ω 1 s(ζ) = β −2 (ζ) − q 2 2 β(0) β(zs ) = q −1 = sin ϕ. R
Z = ±i
ree tion oe ient in the WKBJapproximation
(3.152)
130
CHAPTER 3.
total ree tion ).
Its absolute value is 1 (
i.e., a onstant phase shift of
±π/2
for
BODY WAVES
It des ribes only the
ω < 0 (ω > 0)
phase shifts,
is added to the phase
shifts due to the travel time in the exponential term of
R.
Compared to the
ree tion oe ients of layered media, derived earlier without approximation; the form of (3.152) is
simple.
It is only valid for su iently high frequen ies
( ondition (3.134)) and for angle of in iden e
ϕ with
total ree tion. Ree tion
oe ients of the type of (3.152) are useful in seismology but even more so for the propagation of sound waves in o eans or the propagation of radio waves in the ionosphere ( ompare, e.g., Budden (1961) and of Tolstoy and Clay (1966)). The
ree ted SHwave
observed at the oordinate entre follows, then, by in
serting (3.152) in (3.144)
v1 (0, 0, t) = A(0)i
Z zs ω s(ζ)dζ . exp iω t − 2 ω 0
(3.153)
Then
τ (q) = 2
Z
zs
s(ζ)dζ
0
is the
delay time,
, i.e., the time between the interse tion of the in ident and
the ree ted wave with the oordinate entre. This time delay orresponds to the ray segments
AC, BD, or OP.
Fig. 3.60: Constru tion of a austi from the envelopes of rays.
Note that the wavefronts are
z< 0 are the wavefronts plane.
Fig.
urved
3.60 shows also that the line
for
z>0;
z = zs
only in the homogeneous region
is the envelope of all rays.
Su h
envelopes are alled austi s, and they are hara terised by large energy on entrations. Within the amplitude approximation formula (3.140), and due to
3.10.
131
WKBJ METHOD
dl(P)= 0
z = zs , innite amplitudes would result ± π2 , as dis ussed before, an be interpreted
for a point P on the austi
there. The additional phase shift of
physi ally as the ee t of the strong intera tion of ea h ray with its neighbouring rays in the vi inity of the austi . More ompli ated austi s o
ur in verti ally
inhomogeneous media, if the in ident wave is from a point or line sour e, reπ spe tively. In su h ases, the phase shift per austi en ountered, is ± . More 2
ompli ated austi s are en ountered in two and three dimensional media.
3.10.2 Impulsive ex itation and WKBJseismograms If an
impulsive wave,
produ ing a displa ement
v0 (0, 0, t) = F (t)
at the oor
dinate entre, instead of a harmoni wave is in ident, it follows from (3.153)
F(t)
that the orresponding ree tion is the time delayed Hilbert transform of ( ompare se tion 3.6.3)
v1 (0, 0, t) = FH (t − τ (q)). This means that the ree tion for or horizontal slowness
q)
all
(3.154)
angles of in iden e
ϕ
(or ray parameters
have the same form ex ept for the time delay
τ (q).
This is, therefore, dierent from the results for a dis ontinuity of rst order in
hapter 3.6.3 (there the impulse form hanged in the ase of total ree tion also with the angle of in iden e When
ylindri al waves
ϕ).
are onsidered, the prin iple of superposition is used.
line sour e in the oordinate entre, is represented by many plane waves with radiation angles ϕ from 0 to π/2. This orresponds to positive values of q. The
The ylindri al wave, assumed to originate from an isotropi ally radiating
ree tions are superimposed similarly. First, (3.154) is generalised for arbitrary
x>0 v1 (x, 0, t) = FH (t − τ (q) − qx) = FH (t) ∗ δ (t − τ (q) − qx) . Then these plane waves are integrated over
WKBJseismogram
at distan e
v(x, 0, t) = FH (t) ∗ The
Z
0
π 2
x
ϕ
from 0 to
π/2
and, thus, the
from the line sour e is derived
δ (t − τ (q) − qx) dϕ = FH (t) ∗ I(x, t).
(3.155)
impulse seismogram I(x, t) an now be derived numeri ally via I(x, t) ti
Usually, the
P = i δ(t − ti )∆ϕi = τ (qi ) + qi x , qi =
sin ϕi β(0) .
ϕi are hosen equidistant (∆ϕi = ∆ϕ = const).
(3.156)
The delta fun tions
are shifted from the times ti to their immediate neighbouring time points
I(x, t)
132
CHAPTER 3.
BODY WAVES
and possibly amass there, i.e., in the dis retised version of I(x, t), multiples of ∆ϕ o
ur there. I(x, t) is then onvolved with FH (t). The most time onsuming part is the omputation of the delay time τ (qi ); on the other hand, e ient rayseismi methods exist for that. In omparison to the ree tivity method and the GRT, the WKBJmethod is signi antly faster. There are also other numeri al realisations of this method than (3.155) and (3.156). WKBJseismograms have
other phase relations and impulse forms
pe ted from (3.152) and (3.154), respe tively.
than ex
This is due to the summation
of many plane waves. Prograde travel time bran hes (see Fig. 3.56) show no phase shift, i.e., the impulse form of the in ident ylindri al wave is observed there. Phase shifts and impulse form hanges only o
ur on retrograde travel time bran hes. Furthermore, the seismogram amplitudes are nite in the vi inity of the turning points of travel time urves, i.e., the WKBJmethod is valid at austi s. WKBJseismograms for a simple rustmantle model and a line sour e at the Earth's surfa e are shown in Fig. 3.61. The omputations were performed with a program for
SH waves; even so, the velo ity model (Fig.
3.62) is valid for P
waves. An a ousti Pwave omputation would give, in prin iple, the same result for
pressure.
The travel time urve of the ree tion
PM P
from the rustmantle
boundary (Moho) is retrograde and the travel time bran h of the refra ted wave
Pn
from the upper mantle is prograde. Times are redu ed with 8 km/s.
Fig. 3.61: WKBJseismograms for a simple rustmantle model and a line sour e at the Earth's surfa e.
3.10.
133
WKBJ METHOD
Fig. 3.62: Velo ity model used in Fig. 3.61.
The impulse forms of these waves are as expe ted: radiated wave, whereas
PM P
Pn
has the form of the
is roughly the Hilbert transform of it. At distan es
smaller than the riti al distan e ( a. 100 km), the amplitudes in rease strongly. At the riti al point, whi h is lo ated on a austi within the rust, the wave eld remains nite. In seismologi al appli ations of WKBJseismograms, their approximate nature, due to the high frequen y approximation (3.152) for the ree tion oe ient, should be kept in mind.
This approximation is insu ient in regions of the
Earth where wave velo ity and density hange rapidly with depth, e.g., at the
oremantle boundary or at the boundary of the inner ore of the Earth.
134
CHAPTER 3.
BODY WAVES
Chapter 4
Surfa e waves 4.1 Free surfa e waves in layered media 4.1.1 Basi equations In addition to the body waves that penetrate to all depths in the Earth, another type of wave exists whi h is mostly limited to the neighbourhood of the surfa e of the Earth alled
surfa e waves.
These waves propagate along the surfa e of the
Earth, and their amplitudes are only signi ant down to the depth of a few wave lengths. Below that depth, the
displa ement is negligible.
Be ause surfa e waves
are onstrained to propagate lose to the Earth's surfa e, their amplitude de ay as a fun tion of sour e distan e is smaller than for body waves, whi h propagate in three dimensions.
This is why surfa e waves are usually the dominating
signals in the earthquake re ord. Another signi ant property is their
dispersion,
i.e., their propagation velo ity is frequen y dependent. Therefore, the frequen y within a wave group varies as a fun tion of time ( ompare example in 4.1.4).
These are some of the main observations and explanations for surfa e waves. The rst s ientists to study surfa e waves (Rayleigh, Lamb, Love, Stoneley et al.) found the theoreti al des riptions explaining the main observations. The fundamental tenet of this approa h is the des ription of surfa e waves as an
eigenvalue problem but omitting the sour e of the elasti waves.
We onsider a
layered halfspa e with parameters as given in Fig. 4.1 with a free surfa e. 135
136
CHAPTER 4.
SURFACE WAVES
Fig. 4.1: Layered halfspa e with free surfa e.
We work with Cartesian oordinates
x,y,z
and assume
independen e from y.
Then, the equations from se tion 3.6.2 an be used whi h have been derived to des ribe ree tion and refra tion. The separation of the displa ement into
P, SV and SH  ontributions holds as before, as does the fa t, that the PSVSHwaves. Surfa e waves of the PSVtype are alled Rayleigh waves. They are polarised in the xzplane
ontributions propagate independently from the
The potentials
Φ
horizontal displa ement
u
=
∂Φ ∂x
−
∂Ψ ∂z
verti al displa ement
w
=
∂Φ ∂z
+
∂Ψ ∂x .
and
Ψ,
in ea h layer, satisfy the wave equation
∇2 Φ = Surfa e waves of the
(4.1)
1 ∂2Φ , α2 ∂t2
∇2 Ψ =
1 ∂ 2Ψ . β 2 ∂t2
(4.2)
SH type are alled Love waves. They are polarised in y v in ea h layer the following wave equation
dire tion, and for the displa ement holds
∇2 v =
1 ∂2v . β 2 ∂t2
The boundary onditions are given in se tion 3.6.2. The
(4.3)
ansatz for Φj , Ψj
in the jth layer of the model for harmoni ex itation (ω
Φj Ψj
=
Aj (z) Bj (z)
> 0)
and
x i Aj (z) = exp [i (ωt − kx)] exp iω t − Bj (z) c h
vj
is
(4.4)
4.1.
137
FREE SURFACE WAVES IN LAYERED MEDIA
and
vj
h x i = Cj (z) exp [i (ωt − kx)] . = Cj (z) exp iω t − c
(4.5)
Consider the onditions for whi h a plane wave exists, whi h propagates in xdire tion with the is
?
phase velo ity , where
Then the fun tions
Aj (z), Bj (z)
and
is identi al in all layers. How large
Cj (z)
exist, so that
An (z) Bn (z) lim = 0. z→∞ Cn (z)
This problem is an
(4.6)
eigenvalue problem, and and the wavenumber k = ω/c, eigenvalues. In many ases (for xed ω ), a
respe tively, are the orresponding nite number
(≤ 1)
of eigenvalues exist. The problem is omparable to that of
natural os illations (or free os illations) of nite et .) ( ompare exer ise 4.1).
determining the frequen ies of bodies (beams, plates, bodies
Here we are only interested in the ase where the eigenvalues are real (>0). This has the largest pra ti al appli ation. The orresponding surfa e waves are alled
normal modes.
There exist also waves whi h an be des ribed with omplex
k = k1 − ik2 (k1,2 > 0).
These surfa e waves are alled
their amplitude de reases exponentially with is
exp(−k2 x).
leaking modes
k:
sin e
Their phase velo ity
ω/k1 .
The ansatz with plane waves negle ts the inuen e of ex itation. This, then, leads to a major simpli ation of the problem. Su h surfa e waves are alled in ontrast to
for ed
free
surfa e waves whi h are ex ited by spe i sour es. The
analogy to the free and for ed resonan es of limited bodies is also helpful here in the ontext of ex itation. The treatment of free surfa e waves is an important requirement for the study of for ed surfa e waves ( ompare se tion 4.2). We will soon show that the dispersive properties of both wave types are identi al. Sin e this property depends on the medium, they an be used to determine medium parameters. This is why the study of the dispersion of free surfa e waves is of great pra ti al importan e.
Exer ise 4.1 P
The radial os illations of a liquid sphere with velo ity α are des ribed by Φn (r, t) ∼ (eiωn t/r ) sin(ωn r/α) (n=1,2,3...). Determine the eigen
the potential frequen ies
ωn
from the ondition that the surfa e of the sphere at
r=R is stress
free (prr from exer ise 3.4); give the radial displa ement. Where are the nodal planes?
138
CHAPTER 4.
SURFACE WAVES
4.1.2 Rayleigh waves at the surfa e of an homogeneous halfspa e The halfspa e (z>0) has the velo ities
α and β
for
P  and S waves, respe tively.
Inserting the ansatz
Φ = A(z) exp [i (ωt − kx)] into the wave equation (4.2) with equations for
A(z) and B(z), e.g.,
and
Ψ = B(z) exp [i (ωt − kx)]
∇2 = ∂ 2 /∂x2 + ∂ 2 /∂z 2 ,
A′′ (z) + k 2
(4.7)
gives the dierential
c2 − 1 A(z) = 0. α2
The general solution of this equation is
A(z) = A1 e Due to (4.6),
δ
−ikδz
+ A2 e
ikδz
with
δ=
has to be purely imaginary. Then,
δ , a rst statement on the phase possible: c < α. It also holds that
the properties of be omes
c2 −1 α2
A2 = 0
21
.
has to hold. From
velo ity of the Rayleigh wave
A(z) = A1 e−ikδz , and, similarly, it follows that
B(z) = B1 e−ikγz with
γ = c2 /β 2 − 1
1/2
(negative imaginary). This further limits
: c < β.
The potential ansatz (4.7) an now be written as
Φ = A1 exp [i (ωt − kx − kδz)] , Inserting the
Ψ = B1 exp [i (ωt − kx − kγz)] .
boundary onditions pzz = pzx = 0 for z= 0 with pzz
(3.29), it follows that
−
ω2 λ = A1 + 2 −k 2 δ 2 A1 − k 2 γB1 α2 µ −2k 2 δA1 + −k 2 + k 2 γ 2 B1 =
0 0.
and
pzx
(4.8)
from
4.1.
139
FREE SURFACE WAVES IN LAYERED MEDIA
Division by
−k 2
λ/µ = (α2 − 2β 2 )/β 2
and use of
c2 α2 − 2β 2 +2 α2 β2
c2 −1 α2
gives
A1 + 2γB1 c2 2δA1 + 2 − 2 B1 β
= 0 = 0.
This leads to
− 2 A1 + 2γB1 2 −2δA1 + βc 2 − 2 B1 c2 β2
0
= =
0.
This system of equations only has nontrivial solutions minant is zero.
This leads to an equation for :
In the range of interest
c2 −2 β2
0 < c < β,
c2 −2 β2
2
2
A1
(4.9)
and
B1 ,
if its deter
+ 4δγ = 0.
we have
1 1 c2 2 c2 2 1− 2 . =4 1− 2 β α
Squaring this gives
c4 β 2 c2 β2 c2 c6 = 0. − 8 4 + 24 − 16 2 − 16 1 − 2 β2 β6 β α β2 α Solution
c=0
(4.10)
is not of interest, therefore, only the terms in the bra ket have
to be examined. For
c = 0,
it is negative, and for
c = β , positive. Therefore, at β . The eigenvalue problem
least one real solution of (4.10) exists between 0 and
has, thus, a solution, i.e., along the surfa e of a homogeneous halfspa e a wave
an propagate, the amplitudes of whi h de ay with depth. In this simple ase
no dispersion
o
urs and
In the spe ial ase
λ=µ
ω. √ α = β 3), c = 0, 92β .
is independent of
(i.e.,
wave is only slightly slower than the We now examine the
displa ement
S wave.
In general, the Rayleigh
of the Rayleigh wave
u = −ikA1 e−ikδz + ikγB1 e−ikγz exp [i (ωt − kx)] .
140
CHAPTER 4.
With
γB1 = − c2 /2β 2 − 1 A1
u = −ikA1 e
−ikδz
e
+
−ikδz
+
SURFACE WAVES
(from (4.9)), it follows
c2 −ikγz exp [i (ωt − kx)] −1 e 2β 2
(4.11)
c2 − 1 e−ikγz =: a(z). 2β 2
Similarly,
w = −ikδA1 e−ikδz − ikB1 e−ikγz exp [i (ωt − kx)] .
With
B1 = δA1 c2 /2β 2 − 1 "
+
e−ikδz +
A1
−ikδz
(from (4.9)), it follows
c2 −1 2β 2
w = −ikδA1 e
We assume that
−1
−1
e
−1
e−ikγz =: b(z).
c2 −1 2β 2
−ikγz
#
exp [i (ωt − kx)]
(4.12)
is positive real and onsider the real parts of (4.11) and
(4.12)
u = w = For
z= 0,
kA1 a(z) sin(ωt − kx) −δkA1 b(z) cos(ωt − kx).
it holds that a(0)>0 and b(0) β2 , γ1
η2 + 1 µ1 γ1 + µ2 γ2 = η2 − 1 µ1 γ1 − µ2 γ2
is negative imaginary and
(4.18)
c < β2 .
is also negative imaginary. Then, the right side of (4.18) is real
and larger than 1. The left side is also real, but smaller than 1. Therefore, no real solution The
S velo ity
layer, i.e.,
of (4.18) exists in this ase. in the halfspa e, therefore, has to be larger than that of the
β2 > β1 .
In this ase, we an ex lude values of
with the same arguments as for and
β2 .
In this ase
γ1
β1 > β2 .
between 0 and
This leaves values for
between
β1 β1
is positive real, and both sides in (4.18) are omplex.
The absolute value of both sides is 1. Thus, the mat hing of the phases gives the
dispersion equation of the Love waves
−2kγ1 d1 = −2 arctan From this, it follows with
µ2 γ2  . µ1 γ1
ω = kc
µ2 c−2 − β2−2
µ1 β1−2 − c−2
12 12
h 1 i = tan ωd1 β1−2 − c−2 2 .
This trans endent equation is solved by sele ting
(4.19)
with β1 < c < β2 and inver
sion of the radi and, thus, giving the orresponding ω (or the the orresponding ω ′ s). For a dis ussion, we introdu e a new variable
general
x
146
CHAPTER 4.
x(c)
= ωd1 β1−2 − c−2
c(x)
=
ωd1
ω 2 d2 β2 1
1
−x2
21 .
12
µ1
(4.20)
The left side of (4.19) an then be written as
µ2 c−2 − β2−2
SURFACE WAVES
21
1 µ2 2 2 −2 ω d1 β1 − β2−2 − x2 2 = f (x, ω). 1 = µ1 x β −2 − c−2 2 1
Equation (4.19) an then be expressed (see also Fig. 4.6) as
f (x, ω) = tan x.
Fig. 4.6:
f (x, ω)
f (x, ω)
and tan x.
is real between
x= 0 and its zero x0 = ωd1 β1−2 − β2−2
This zero moves to the right as a fun tion of se tions of in
f (x, ω)
with
tan x.
ω
The interse tion
21
.
(4.21)
and reates, thus, more inter
xi = xi (ω)
gives, substituting
(x) from (4.20), the dispersion relation of the ith mode (i = 1, 2, . . .) ci (ω) =
ωd1 ω 2 d21 β12
−
12 .
x2i (ω)
(4.22)
4.1.
147
FREE SURFACE WAVES IN LAYERED MEDIA
The ith mode o
urs only if be omes
νi =
x0 > (i − 1)π .
uto frequen y
i−1 ωi = 1 . 2π 2d1 β1−2 − β2−2 2
The ith mode exists only for frequen ies mode,
With (4.21), its
ν > νi .
(4.23)
The rst mode (fundamental
i= 1) exists, due to ν1 = 0, at all frequen ies.
The phase velo ity of ea h mode at its uto frequen y is most simply derived from the fa t that at this point the tangent is zero in (4.19)
ci (ωi ) = β2 . For
ω → ∞, xi → (i−1/2)π and the tangent in (4.19) approa hes ∞. lim ci (ω) = β1 .
ω→∞
(4.24) Therefore,
(4.25)
An upperlimiting frequen y does not exist, and the velo ities of the layer and the halfspa e are the limiting values of the phase velo ity. A al ulated example for the dispersion urves of the rst three Love modes of a rustmantle model, onsisting of a halfspa e with a layer above, is given in Fig. 4.7.
Fig.
4.7:
model.
Dispersion urves of the rst three Love modes of a rustmantle
148
CHAPTER 4.
In addition to the urves of the phase velo ities
,
the
SURFACE WAVES
group velo ities U
are
shown
U=
dω dc dc c c =c+k =c−Λ = = ω dc dk dk dΛ 1 − c dω 1 + Tc
(Λ= wavelength,
T = period).
dc dT
.
(4.26)
The group velo ity ontrols, as we will see, the
propagation of an impulse from the sour e to a re eiver, i.e., ea h frequen y travels with its group velo ity from the sour e to the re eiver,
not
with its
phase velo ity. Phase and group velo ities an be determined from observations (see se tion 4.1.4 and 4.1.5). Thus, both an be used for interpretation.
Nodal planes and eigen fun tions Finally, we examine the
nodal planes
of the Love modes for the ase just dis
ussed, i.e., the surfa es on whi h the displa ement is zero. From (4.13) with (4.14) and (4.15), it follows that
v1 v2 where
v1
and
v2
= =
2E1 cos(kγ1 z) exp [i (ωt − kx)] E2 exp (−ikγ2 (z − d1 )) exp [i (ωt − kx)]
are for the layer and the halfspa e, respe tively. These expres
sions also hold for ea h individual mode, and the dispersion relation (4.22) has to be used. The relation between
E1
E2
and
is
E2 = 2E1 cos(kγ1 d1 ); E1
an be
hosen arbitrarily. This means that at the surfa e
z = 0,
the
maximum displa ement
is always
observed, and that nodal planes an only o
ur in the layer, but not in the halfspa e. Their position is determined by the zeros of the osine. For the ith mode they an be derived via the equation
kγ1 z = ω β1−2 − c−2 i where
Ni
21
z = (2n − 1)
is their number determined by
For the lower frequen y limit
ω = ωi ,
(i − 1)π
π 2
(n = 1, 2, . . . , Ni ),
z ≤ d1 .
from (4.23), it follows due to (4.24)
π z = (2n − 1) . d1 2
This is satised for
z=
2n − 1 d1 2i − 2
with
(4.27)
n = 1, 2, . . . , Ni = i − 1.
4.1.
For
149
FREE SURFACE WAVES IN LAYERED MEDIA
ω = ωi ,
therefore,
i−1
nodal planes exist. Their spa ing is
d1 (i = 3, 4, . . .). i−1
∆z =
i= 1) has no nodal plane, neither for ω = ω1 = 0 nor for nite
The rst mode (
ω > 0. The other extreme on the frequen y s ale of ea h mode is dis ussion of the behaviour of
lim ω
ω→∞
β1−2
−
ci (ω)
c−2 i
12
for
ω→∞
xi (ω) = = lim ω→∞ d1
Equation (4.27) leads to the fa t that all satisfy the relation
ω = ∞.
From the
(see (4.25)), it follows that
z ≤ d1
1 π i− . 2 d1
have to be determined whi h
1 z π i− π = (2n − 1) . 2 d1 2 These are the values
z=
For
ω = ∞,
therefore,
i
2n − 1 d1 2i − 1
with
n = 1, 2, . . . , Ni = i.
nodal planes exist with the spa ing
∆z =
d1 i − 12
(i = 2, 3, . . .).
The hange relative to the situation where
ω = ωi
holds, is rst, the de rease in
the spa ing of the nodal planes, se ond a general move to shallower depth and nally, the addition of the ith nodal plane
z = d1 .
This means that for
the halfspa e remains at rest.
ω=∞
Fig. 4.8 shows quantitatively the amplitude behaviour of the rst three modes as a fun tion of depth for the rustmantle model used for Fig. 4.7. The period is also indi ated. Su h amplitude distributions are alled
eigen fun tions
modes. They follow from the zdependent part of the displa ement dis ussed above.
v1
of the
and
v2
150
CHAPTER 4.
SURFACE WAVES
Fig. 4.8: Quantitative amplitude behaviour of the rst three modes as a fun tion of depth for the rustmantle model used for Fig. 4.7.
Exer ise 4.3 Derive the dispersion equation for free Rayleigh waves in a liquid medium onsisting of a layer over a halfspa e and ompare this to (4.19). Sket h a gure similar to Fig. 4.6. What is the dieren e, espe ially for the rst mode?
4.1.4 Determination of the phase velo ity of surfa e waves from observations In the last se tion, we saw how the phase velo ity of a Love mode in a layered halfspa e an be determined if the halfspa e is known. In the same way (but with more ompli ations), the same an be done for Rayleigh waves. We now dis uss the derivation of the phase velo ity from only
one
mode is present.
observations. We assume that lters have to be used to
If that is not the ase,
separate the dierent modes. Sin e these are sometimes ompli ated methods, they are not dis ussed here. An overview, and appli ations for surfa e waves,
an be found, e.g., in Aki and Ri hards, Dahlen and Tromp, and Kennett. We work here with plane surfa e waves, i.e., we negle t the sour e term. The method for the determination of the phase velo ity thus derived, is su iently a
urate for pra ti al purposes.
The
modal seismogram
of any displa ement
omponent at the Earth's surfa e an be written for propagation of the mode in
x dire tion as
4.1.
151
FREE SURFACE WAVES IN LAYERED MEDIA
1 u(x, t) = 2π
Z
+∞
−∞
x dω. A(ω) exp iω t − c(ω)
(4.28)
This is a superposition of the previously dis ussed harmoni surfa e waves with
c(ω) u(x, t).
the aid of the Fourier integral. derived from the re ordings of mograms for dierent at (arbitrary)
x= 0.
Φ(ω) = arg A(ω),
x
are dierent.
is the phase velo ity of the mode to be Sin e
A(ω)
is frequen y dependent, the seis
is the spe trum of the displa ement
The amplitude spe trum is
i.e.,
A(ω)
and the phase spe trum
A(ω) = A(ω)eiΦ(ω) . We assume that
u(x, t)
x = x1
is known for
and
x = x2 > x1
and apply a
Fourier analysis to the seismograms
u(x1,2 , t) = where
G1,2 (ω)
1 2π
Z
+∞
G1,2 (ω)eiωt dω,
(4.29)
−∞
is the spe trum of the seismogram for
x = x1,2 .
The omparison
of (4.29) with (4.28) gives
G1,2 (ω) = G1,2 (ω)eiϕ1,2 (ω) = x1,2 x1,2 = A(ω) exp i Φ(ω) − ω . A(ω) exp −iω c(ω) c(ω) If the time
t= 0 is identi al for both seismograms, and, if possible, jumps of ±2π ϕ1,2 (ω) G1,2 (ω)...)
have been removed from the numeri ally determined phases
−π and +π ), the phases of the top ( A(ω)...) an be mat hed and give
between bottom
ϕ1,2 (ω) = Φ(ω) − ω Subtra ting
ϕ1 (ω)
from
ϕ2 (ω),
(observation
(usually and the
x1,2 . c(ω)
the unknown phase spe trum
Φ(ω)
of
u(0, t)
an els and the following result for the phase velo ity is left
c(ω) =
(x2 − x1 )ω . ϕ1 (ω) − ϕ2 (ω)
(4.30)
For pra ti al appli ations of this method, the surfa e waves have to be re orded at two stations whi h are on a great ir le path with the sour e. In ase
three
stations are available, this requirement an be ir umvented by onstru ting a triangle between the stations. In both approa hes, the phase velo ities derived are representative for the region
between
the stations.
152
CHAPTER 4.
SURFACE WAVES
Fig. 4.9: Example for seismograms of Rayleigh waves from L. Knopo, et al., 1966. The tra es have been shifted.
The interpretation based on the phase velo ity from Fig. 4.9 is given in Fig. 4.10. It shows short period group velo ity observations from near earthquakes as well as phasevelo ity measurements for the region of transition (Central Alps to northern Foreland, Fig. 4.9) (from L. Knopo, St. Müller and W.L. Pilant: Stru ture of the rust and upper mantle in the Alps from the phase velo ity of Rayleigh waves, Bull. Seism. So . Am.
Fig. 1966.
56, 10091044, 1966).
4.10: Short period group velo ity observations from L. Knopo, et al.,
4.1.
153
FREE SURFACE WAVES IN LAYERED MEDIA
4.1.5 The group velo ity Seismograms with dispersion, as shown in Fig. variation of frequen y with time, so that
one
4.9, often have a very slow
frequen y an be asso iated with
a ertain time. If this is done for two dierent distan es
x1
and
x2 ω
ir le through the sour e), and if the times, for whi h frequen y are t1 (ω) and
t2 (ω),
(on a great is observed
it follows that
U (ω) =
x2 − x1 . t2 (ω) − t1 (ω)
(4.31)
U (ω) is the velo ity with whi h this frequen y, or a wave group with a small frequen y band ∆ω around the frequen y ω , propagates. U is, therefore, alled the group velo ity. The theory of surfa e waves from point sour es in se tion 4.2 gives the even simpler formula seismogram;
r
U (ω) = r/t(ω), whi h only requires one t(ω) is relative to the time
is the distan e from the sour e, and
when the wave started (sour e time).
In pra ti e, this seismogram is ltered
in a narrow band with the entral frequen y
ω.
The arrival time
t(ω)
is at
the maximum of the envelope of the ltered seismogram. The group velo ity
an, therefore, in prin ipal be determined without di ulty from observations. Another question is, how the group velo ity is onne ted with the phase velo ity and, thus, with the parameters of the Earth, i.e., the velo ity of
P  and S waves
and density, as a fun tion of depth. To study this relation, we start from the des ription of the modal seismogram (4.28) and express it using the wavenumber
u(x, t) =
For su iently large
x
1 2π
Z
k = ω/c
as
+∞
−∞
A(ω) exp [i (ωt − kx)] dω.
(4.32)
t
and, therefore, also large , the phase
ϕ(ω) = ωt − kx is
rapidly varying ompared to A(ω).
ω ϕ(ω)
has hanged by
2π ,
whereas
This means, for example, that for hanging
A(ω)
is pra ti ally un hanged, and the
ω
interval, therefore, does usually not ontribute to the integral (4.32). This is espe ially true if
ω = ω0 ,
ϕ(ω) an be approximated linearly. This no longer holds for ϕ(ω) has an extremum (see Fig. 4.11), i.e., when it be omes
for whi h
stationary.
154
CHAPTER 4.
Fig. 4.11: Extremum of This frequen y
ω0
SURFACE WAVES
ϕ(ω).
depends on
t
and follows from
ϕ′ (ω0 ) = t − x
dk ω=ω0 = 0. dω
(4.33)
The frequen y ω0 , whi h satises (4.33), dominates at time t in the modal seismogram. Sin e for plane surfa e waves lo ation and time origin are arbitrary, (4.33) an be written for two distan es and
t2 (ω0 ),
lo ity
x1
and
x2 and orresponding times t1 (ω0 )
respe tively. Subtra tion gives the
basi formula for the group ve
x2 − x1 dω = U (ω0 ) = ω=ω0 . t2 (ω0 ) − t1 (ω0 ) dk ω
and
k
are onne ted via the phase velo ity
c = ω/k .
(4.34)
Using this, the expli it
group velo ity (4.26) an be derived dire tly from (4.34) (see exer ise 4.4). The arguments sket hed here are the entral ideas of the
phase.
method of stationary
We will use it later to al ulate integrals of the form (4.32) approxi
mately; here, it was only used to demonstrate that it is the group velo ity whi h determines the sequen e and possible interferen e in the modal seismogram. To make this statement more obvious, we onsider an arbitrary mode of the
Rayleigh waves of a liquid halfspa e with a layer at the top.
Its dispersion urves
for phase and group velo ity look like those in Fig. 4.12 ( ompare exer ise 4.3).
Fig. 4.12: Dispersion urves for a liquid halfspa e with a top layer.
4.1.
155
FREE SURFACE WAVES IN LAYERED MEDIA
Instead of (4.31), or the left of (4.34), we use
U (ω) = r/t(ω),
whi h refers to
the sour e, to interpret the urve of the group velo ity. From the trend in
Un , we on lude that for the mode onsidered, the frequen ies ωn arrive at an arbitrary
in the neighbourhood of the lower limiting frequen y distan e
r
from the sour e.
This assumes that su h frequen ies are a tually
ex ited at the sour e. Their group velo ity is is
r/α2 .
α2 ,
and their group travel time
For later times, whi h are still smaller than
r/α1 ,
the frequen y of the
arriving os illations slowly in reases, orresponding to the steep trend in the
urve of
Un .
greater then
This wave train is alled the
r/α1
there are
two
frequen ies,
fundamental wave. ω′
and
ω ′′ ,
At later times
whi h ontribute to the
seismogram. This has the ee t that the higher frequen y waves (water waves) ride on the fundamental waves. The frequen ies of the two waves approa h ea h
r/Unmin .
other for in reasing time and be ome identi al at time
Here,
is the minimal group velo ity. The orresponding wave group is the and it onstitutes the end of the modal seismogram.
Unmin
Airy phase,
Exa t omputations of
modal seismograms, dis ussed later, onrm these qualitative statements. The example in Fig. 4.13 shows the behaviour of pressure of the fundamental mode (from C.L. Pekeris: Theory of propagation of explosive sound in shallow water, Geol. So . Am. Memoir No. 27, 1948).
Fig. 4.13: Pressure of the fundamental mode from C.L. Pekeris, 1948. 3 3 1500 m/s , α2 = 1650 m/s, ρ1 = 1 g/ cm , ρ2 = 2 g/ cm , d1 = 20 m.
The dispersion of the fundamental wave of the example in Fig.
α1
=
4.13, shown
for a sour e distan e of 9200 m , i.e., de reasing group velo ity with in reasing frequen y (in rease of frequen y with time), is alled
regular dispersion.
For the
water wave, the group velo ity grows with the frequen y (frequen y in reases with in reasing time); this is alled
inverse
dispersion.
The notation regular
and inverse dispersion should not be onfused with the expressions
normal
and
156
CHAPTER 4.
anomal
SURFACE WAVES
dispersion, whi h express, that the group velo ity is larger or smaller
than the phase velo ity, respe tively. Fig. 4.14 is a sket h to demonstrate the basi propagation properties of a dispersive wave train. It assumes that the sour e radiates an impulse with onstant spe trum in the frequen y band only regular dispersion.
ω1 ≤ ω ≤ ω2
and that the medium produ es
The larger the propagation distan e, the longer the
wave train be omes. At the same time, the amplitudes de rease (not shown in Fig. 4.14).
Fig. 4.14: Basi propagation properties of dispersive wave trains.
Constant frequen ies with a slope of
dr/dt
o
ur on
straight lines
through the origin (r
e.g., a ertain maximum or an interse tion with zero, are
lines,
= 0, t = 0)
Constant phases, situated on urved
that is identi al to the group velo ity.
and the frequen y varies along these urves.
The lo al slope
dr/dt
of
these urves is the phase velo ity for the dominating frequen y at this time.
Exer ise 4.4 Derive (4.26) for the group velo ity. How does abnormal dispersion, respe tively?
dc/dω
behave for normal and
4.1.
157
FREE SURFACE WAVES IN LAYERED MEDIA
Exer ise 4.5 a) What is the form of the most general phase velo ity velo ity
U (ω)
c(ω) for whi h the group
is onstant? Interpret the orresponding seismogram (4.28).
c1 (ω) and 1/U = dk/dω =
b) What is the most general onne tion between phase velo ities
c2 (ω) with identi al d(ω/c)/dω .
group velo ities
U1 (ω)
U2 (ω)?
and
Use
4.1.6 Des ription of surfa e waves by onstru tive interferen e of body waves Up to this point, surfa e waves have been treated for the most part theoreti ally, namely based on an ansatz for the solution of dierential equations. Input in these equations have been the on entration of the wave amplitude near the surfa e, propagation along the surfa e and dispersion. We have not rea hed a physi al understanding how these waves an be onstru ted. In this se tion, we will show for the simple example of Love waves in a halfspa e with one layer at the top, that surfa e waves an basi ally be understood as arising from onstru tive interferen e of body waves whi h are ree ted between the interfa es. We onsider SH body waves whi h propagate up and down in the layer with an angle of in iden e and ree tion
ϕ.
The ree tion at the surfa e is loss
free; the ree tion oe ient a
ording to (3.39) equals +1. During ree tion at the lower boundary of the layer (z = d1 ), energy loss through ree tion ∗ o
ur, as long as ϕ < ϕ = arcsin(β1 /β2 ). In this ase, the amplitude of the ree ted waves de rease with the number of ree tions at the lower boundary. ∗ If on the other hand, ϕ > ϕ , the ree tion oe ient at the lower interfa e has the absolute value of +1 (see (3.42)).
Thus, no wave in the lower half
spa e exists transporting energy away from the interfa e and the amplitude of ea h
single
multiple ree tion is preserved.
The wave eld in this layer is
then basi ally ontrolled by the interferen e of
ertain values of
ϕ
all
multiple ree tions.
For
there will be onstru tive interferen e and for other values
there will be destru tive interferen e, respe tively. We try to determine those
ϕ
whi h show onstru tive interferen e. To a hieve this, we approximate the
momentary waveeld pi ture of Fig. 4.15, lo ally, by plane parallel wavefronts with a orresponding angle of in iden e
ϕ
(see Fig. 4.16).
Fig. 4.15: Pi ture of momentary wave eld.
158
CHAPTER 4.
SURFACE WAVES
Fig. 4.16: Approximation of Fig. 4.15 by plane, parallel wavefronts. This approximation is better at larger distan es from the sour e. The limitation on
plane
free
waves means that we onsider
normal modes.
The phases of neighbouring wavefronts 1, 2, 3 are and
Φ 3 = Φ1 + ǫ 1 + ǫ 2 ,
A and B, respe tively.
ree tions in
the phase dieren e travel time
respe tively, where
ωs/β1
Φ3 − Φ1
s=2
SH waves.
follows for
ǫ1
and
To ensure that wave 1 and 3 are
2π .
in phase,
With
onstru tive interferen e
2ωd1 cos ϕ + 2nπ, β1
ǫ2
Φ1 (arbitrary), Φ2 = Φ1 + ǫ1 ǫ2 are the phase shifts of the
d1 sin ϕ = 2d1 cos ϕ, tan ϕ
we derive the following ondition for
The phase shifts;
and
has to be equal to the phase dieren e due to the
plus a multiple of
ǫ1 + ǫ2 =
ǫ1
n = 0, 1, 2, . . . .
(4.35)
are the phases of the ree tion oe ients for plane ϕ > ϕ∗ = arcsin(β1 /β2 ), it
Sin e we only onsider post riti al
ǫ1
from (3.42) and with
ω>0
12 2 β2 2 sin ϕ − 1 −ρ β b 2 2 β12 ǫ1 = −2 arctan = −2 arctan . a ρ1 β1 cos ϕ
For the ree tion at
B,
it holds that
ǫ2 = 0 ,
sin e a
ording to (3.39), the
ree tion oe ient at the free surfa e is always +1. Substituting all of the above into (4.35), we get an equation for in iden e
ϕ,
whi h produ e a for given
ρ 2 β2 arctan
β22 β12
sin2 ϕ − 1
ρ1 β1 cos ϕ
12
ω,
those
angles of
onstru tive interferen e
ωd1 cos ϕ + nπ, = β1
n = 0, 1, 2, . . .
(4.36)
4.1.
159
FREE SURFACE WAVES IN LAYERED MEDIA
In this equation, we introdu e the apparent velo ity
c=
β1 sin ϕ
(4.37)
horizontal
with whi h the wavefronts propagate in a dire tion. With cos ϕ = 2 and ρ1,2 β1,2 = µ1,2 , we derive by reversing (4.36), an equation
β1 (β1−2 − c−2 )1/2 for
µ2 c−2 − β2−2 µ1 β1−2 − c−2
21
h 1 i −2 −2 2 . 12 = tan ωd1 β1 − c
This equation is identi al with the dispersion equation (4.19) of Love waves. We would have found the same equation, if we had onsidered waves whi h propagate upwards (and not downwards) in Fig. 4.16. The superposition of both
verti al wavefronts.
is not only an apparent velo ity, but also the phase velo ity of this
groups of waves gives, for reason of symmetry, a wave with Therefore,
resulting wave. We also see that the Love waves in the halfspa e with a top layer are produ ed by
onstru tive interferen e of body waves
whi h have a
post riti al
angle of
in iden e. For these angles of in iden e, no energy is lost from the layer into the halfspa e. The energy remains in the layer whi h a ts as a perfe t
guide.
wave
This is generally true for normal modes of Love and Rayleigh waves in
horizontally layered media, in whi h ase that normal modes exist. From this, we an also on lude that the phase velo ity of normal modes an, equal to the
S velo ity of the halfspa e under the layers
at most, be
c ≤ βn . If it were larger, energy would be radiated into the halfspa e in the form of an
S wave. Leaking modes
also o
ur by onstru tive interferen e of body waves.
In this ase, the angles of in iden e are larger than
βn .
pre riti al,
and the phase velo ity is
Thus, radiation into the lower halfspa e o
urs and the wave
guide is not perfe t. Finally, it should be noted that the explanation of surfa e waves via onstru tive interferen e of body waves annot be applied to the
Rayleigh waves.
fundamental mode of
The Rayleigh wave of the homogeneous halfspa e, for example,
exists without additional dis ontinuities at the surfa e. No simple explanation exists for the fundamental mode of Rayleigh waves.
Exer ise 4.6 Determine the dispersion urves for a liquid layer whose boundaries are (1) both free, (2) both rigid, (3) one rigid and one free, respe tively. Use the arguments of se tion 4.1.6 and ompare with the solution of the orresponding eigenvalue problem. Give the group velo ity and sket h the pressuredepth distributions.
160
CHAPTER 4.
SURFACE WAVES
4.2 Surfa e waves from point sour es 4.2.1 Ideal wave guide for harmoni ex itation Expansion representation of the displa ement potentials We study the propagation of mono hromati sound waves from an explosive point sour e in a liquid layer with a free surfa e situated above a
rigid
half
spa e.
Fig. 4.17: Explosive point sour e in a liquid layer with a free surfa e atop a rigid halfspa e.
This is an ideal wave guide sin e no waves an penetrate the halfspa e. For su h a s enario, the key features of surfa e waves from point sour es an be studied without too mu h mathemati al eort. For the displa ement potential
Φ
in the layer, we assume the following integral
ansatz, using the analogy to (3.84) and applying (3.85) for the potential of the iωt spheri al wave from the sour e. In the following, the timedependent term e is omitted
Φ
=
Z
∞
ω2 − k2 α2
J0 (kr)
0
l
J0 (kr)
=
12
k −ilz−h e + A(k)e−ilz + B(k)eilz dk il
(4.38)
.
is the Bessel fun tion of rst kind and zeroth order,
k
and
zontal and verti al wavenumber, respe tively. The square root
l
are the hori
l is, as in se tions
3.6.5 and 3.7, either positive real or negative imaginary. It an be shown that
Φ
is a solution of the wave equations in ylindri al oordinates. The rst term in (4.38) is the wave from the sour e, the se ond and third orrespond to the waves
4.2.
161
SURFACE WAVES FROM POINT SOURCES
z dire tion,
propagating in positive and negative
respe tively.
A(k)
and
B(k)
follow from the boundary onditions for the interfa es
z=0:
stress
pzz = ρ
z=d:
normal displa ement
∂2Φ = −ρω 2 Φ = 0 ∂t2
or
Φ=0
∂Φ = 0. ∂z
This gives
A(k) + B(k)
=
A(k) − e2ild B(k)
=
k − e−ilh il k − eilh . il
The solution of this system of equation is (please he k)
k cos [l(d − h)] il cos(ld) k sin(lh) −ild e . l cos(ld)
A(k)
= −
B(k)
=
Inserting them into (4.38) gives
0≤z≤h:
Φ
=
h≤z≤d:
Φ
=
∞
sin(lz) cos [l(d − h)] dk l cos(ld) Z0 ∞ sin(lh) cos [l(d − z)] 2 kJ0 (kr) dk. l cos(ld) 0 2
Z
kJ0 (kr)
(4.39)
(4.40)
Before these expressions are solved with methods from omplex analysis, it should be noted that an ex hange of sour e and re eiver does not hange the value of
Φ.
Displa ement and pressure are also the same for this ase. This is
an example for
re ipro ity relations, whi h is important in the theory of elasti
waves. The poles
kn
of the integrand in (4.39) and (4.40) are determined via
dln = d This gives
ω2 − kn2 α2
12
π = (2n − 1) , 2
n = 1, 2, 3 . . .
162
CHAPTER 4.
kn =
ω2 (2n − 1)2 π 2 − 2 α 4d2
12
SURFACE WAVES
.
(4.41)
The innite number of poles are situated on the real axis between
+ω/α and
−ω/α
and
on the imaginary axis, respe tively. The number of poles on the real
axis depends on
ω.
Due to these poles, the integration path in (4.39) and (4.40)
have to be spe ied in more detail.
We hoose path
C1
in Fig.
4.18 whi h
ir umvents the poles in the rst quadrant.
Fig. 4.18: Integration path
C1
whi h ir umvents the poles in the rst quadrant.
In the following, we dis uss only (4.39) in detail. Equation (4.40) an be solved similarly. We use the identity
J0 (kr) =
where
(1)
H0 (kr)
fun tions )
and
(2)
H0 (kr)
i 1 h (1) (2) H0 (kr) + H0 (kr) , 2
Hankel
are Bessel fun tions of the third kind (=
and zeroth order (Appendix C, equations (C.2) and (C.3), respe 
tively). Then,
Φ=
h i sin(lz) cos [l(d − h)] (1) (2) dk. k H0 (kr) + H0 (kr) l cos(ld) C1
Z
Using relation (C.6) from appendix C,
(1)
(2)
H0 (x) = −H0 (−x)
(4.42)
4.2.
163
SURFACE WAVES FROM POINT SOURCES
the rst part of the integral in (4.42) an be written as
−
Z
(2)
C1
kH0 (−kr)
sin(lz) cos [l(d − h)] dk = l cos(ld)
(2)
C2
C2
where u=k is used. The integration path
C1
Z
uH0 (ur)
sin(lz) cos [l(d − h)] du l cos(ld)
is pointsymmetri al to the path
with respe t to the oordinate entre, but it goes from
k
this in (4.42) and with onsistent use of
Φ=
Z
(2)
C
kH0 (kr)
Fig. 4.19: Integration path
C
sin(lz) cos [l(d − h)] dk = l cos(ld)
from
real axis.
Integration path
C,
−∞
−∞
to 0. Inserting
as integration variable, gives
to
+∞
Z
I(k)dk.
ir umventing the poles on the
therefore, is, as indi ated in Fig. 4.19, from
and ir umventing the poles on the real axis. Despite the nonuniqueness of the square root fun tion of
k.
The reason for this is that
sign of the square root of
l
does
not
e.g., if the halfspa e is not rigid, more ompli ated (introdu tion of Now we apply the
I(k)
l
in (4.43),
I(k)
−∞
to
+∞
is a unique
l
is an even fun tion of , thus, the
matter. For more ompli ated wave guides,
I(k)
is not unique and the theory be omes
bran h uts ).
remainder theorem on the losed integration path shown in C and a half ir le with innite radius. The
Fig. 4.20 whi h onsists of path
only
(4.43)
C
singularities in luded are the poles of
I(k).
164
CHAPTER 4.
Fig. 4.20: Integration path
SURFACE WAVES
C
and half ir le with innite radius.
Z
= −2πi
Then
Z L
+
L
C
∞ X
ResI(k)k=kn .
n=1
indi ates the lo kwise integration in the lower plane of Fig. 4.20, and ea h
term in the sum is the residue of If the
asymptoti representation
I(k)
at the pole
k = kn .
of the Hankel fun tion is used, it follows for
large arguments (Appendix C, equation (C.4))
(2) H0 (kr)
≃
2 πkr
12
e−i(kr− 4 ) . π
(4.44)
We see that their values on the semi ir le in the lower halfplane, where a negative imaginary part, be omes zero (for R going to
∞).
k
has
The orresponding
integral also goes to zero, and we have found a representation of the potential
Φ
as an innite sum of residuals. The determination of the residue of the quotient
f1 (k)/f2 (k)
if
kn
at the lo ation
is a pole of
rst
kn
with
f2 (kn ) = 0
is done with the formula
f1 (k) f1 (kn ) Res , = ′ f2 (k) k=kn f2 (kn)
order. In our ase,
kn
follows from (4.41) and is either
positive real or negative imaginary. The orresponding
ln = Furthermore, it holds here, that
l
is
(2n − 1)π . 2d
f2 (k) = cos(ld),
f2′ (kn ) = d
kn sin(ln d). ln
f2 (kn ) = cos(ln d) = 0
and
4.2.
165
SURFACE WAVES FROM POINT SOURCES
Finally,
cos [ln (d − h)] = cos(ln d) cos(ln h) + sin(ln d) sin(ln h) = sin(ln d) sin(ln h). Thus,
Res I(k)k=kn = We, therefore, get the following
sion, for whi h eiωt Φ=−
1 (2) H (kn r) sin(ln z) sin(ln h). d 0
representation of the potential Φ as an expan
has now to be added again for ompleteness
∞ h 2πi X πz i (2) πh sin (2n − 1) H0 (kn r)eiωt . sin (2n − 1) d n=1 2d 2d
This expression not only holds for
0 ≤ z ≤ h
but also for arbitrary depth,
sin e a
ording to (4.40) the same expansion an be found. displa ement omponents
an be derived.
(4.45)
From (4.45) the
∂Φ/∂r and ∂Φ/∂z and the pressure p = −pzz = ρω 2 Φ
For the ideal wave guide the eld an be onstru ted
solely
from the ontribu
tions from the poles, ea h of whi h represents a mode, as will be shown later. For ompli ated wave guides, ontributions in the form of urve integrals in the
omplex
k plane have to be
added to the pole ontributions. These additional
ontributions orrespond mostly to body waves.
Modes and their properties Ea h term in (4.45) represents a
mode.
This is only a denition, but it ts well
into the mode on ept introdu ed in the previous se tions for
free surfa e waves. r, we an use
If we onsider, for example, the terms in (4.45) for large distan es (4.44)
(kn r > 10)
√ π ∞ h −2 2πiei 4 X πz i 1 πh i(ωt−kn r) Φ= . sin (2n − 1) sin (2n − 1) 1 e d 2d 2d 2 (k r) n n=1
(4.46)
kn . Their numω . They orrespond to waves with ylindri al wavefronts whi h propagate in +rdire tion with the frequen y dependent phase
The most important terms in (4.46) are those with positive ber is nite and in reases with velo ity
− 1 ω (2n − 1)2 π 2 α2 2 cn (ω) = =α 1− . kn 4d2 ω 2
(4.47)
166
CHAPTER 4.
If the eigenvalue problem for
SURFACE WAVES
free surfa e waves in the same wave guide is solved nth free normal mode
( ompare exer ise 4.6), it follows for the
h πz i iω Φn = A sin (2n − 1) e 2d with
cn (ω)
t− cnx(ω)
,
(4.48)
from (4.47). Furthermore, the terms in (4.46) and (4.48) agree, that
des ribe the
zdependen e agree. It, therefore, makes sense to name the single nth for ed normal mode, if kn > 0. The dieren e
terms in (4.45) and (4.46) the
with respe t to (4.48) is in the amplitude redu tion proportional to the addition of a term that depends on the sour e depth the
ex itation fun tion of the mode.
h.
r−1/2
and in
This term is named
If the the sour e is lo ated at a nodal plane
of the free mode (4.48), the ex itation fun tion is zero, and the mode is not ex ited. Maximum ex itation o
urs, if the sour e is at a depth where the free mode has its maximum. From the omparison of the free and the for ed normal modes, the importan e of the study of free modes be omes obvious. It des ribes the dispersive properties and the amplitudedepth distributions (eigen fun tions) of the for ed normal modes and, therefore, their most important property.
ompli ated wave guides.
The terms in (4.46) with negative imaginary
kn
This also holds for more
are not waves but represent
os illations with amplitudes that de rease exponentially in
r dire tion.
They
only ontribute to the wave eld near the sour e, where (4.45) has to be used for ompleteness.
The fareld is dominated by normal modes.
The number of
nodal planes of the nth mode is n, and their spa ing is 2d/(2n−1)
(n = 2, 3 . . .).
The potential
have a node for
z=0
Φ,
horizontal displa ement
and a maximum for
z = d,
The opposite is true for the verti al displa ement
∂Φ/∂r
and pressure
p
respe tively (see Fig. 4.21).
∂Φ/∂z .
Fig. 4.21: Modes and nodal planes, n = 1, 2, 3, 4.
The phase velo ity (4.47) of the nth mode an be written as
4.2.
SURFACE WAVES FROM POINT SOURCES
167
ω 2 − 21 n cn (ω) = α 1 − ω
(4.49)
with the lower frequen y limit
ωn =
(2n − 1)πα . 2d
Innitely high phase velo ities an o
ur. A
ording to (4.26), the group velo ity is
Un (ω) = α 1 −
ω 2 12 n
ω
.
(4.50)
Fig. 4.22: Group and phase velo ities. An important property of the ideal wave guide with a rigid and a free interfa e is that the angular frequen ies and waves length
Λ > 4d,
ω < ω1 = πα/2d
respe tively)
(or the frequen ies
ν < α/4d
annnot propagate undamped.
This no
longer holds for the ideal wave guide with two rigid walls ( ompare exer ise 4.6). In this ase, an additional fundamental mode exists, in addition to the modes dis ussed before, but with dierent limiting frequen ies. That mode an o
ur at all frequen ies, and its phase velo ity is frequen y independent and equal to
α.
4.2.2 The modal seismogram of the ideal wave guide In this se tion, we will ompute the orresponding arbitrary summand in (4.45), exa tly. method of stationary phase to do this.
modal seismogram
for an
In the next se tion, we will use the
168
CHAPTER 4.
SURFACE WAVES
The potential (4.45) orresponds to time harmoni ex itation, i.e., for the potential of the explosion point sour e
Φ0 = it holds that
F (t) = eiωt .
1 F R
R , t− α
(4.51)
From this mode of ex itation, we now will move to
the ex itation via a delta fun tion, F (t) = δ(t). Multiplying (4.45), without the iωt fa tor e , with the spe trum of the delta fun tion F (ω) = 1, gives the Fourier transform of the displa ement potential. Finally, the result is transformed ba k into the time domain. These modal seismograms an then be onvolved with realisti ex itation fun tions stood for
F (t),
but the basi features an already be under
F (t) = δ(t).
For this, we onsider the nth mode in the expansion (4.45). Its Fourier transform for ex itation via a delta fun tion is, ex ept for geometry fa tors, equal (2) 2 2 1/2 /α. We now use the Lapla e to H n (ω) = iH0 (kn r) with kn = (ω − ωn ) transform ( ompare se tion A.1.7)
(2)
hn (s) = H n (−is) = iH0 and the relation
(2)
H0 (−ix) =
h r 1 i −i s2 + ωn2 2 α
2i K0 (x) π
between the Hankel fun tion and the modied Bessel fun tion
K0 (x) (see se tion
3.8). This gives then
hr 1 i 2 hn (s) = − K0 s2 + ωn2 2 . π α
The original fun tion an then be found in tables of the Lapla e transform. It r r is zero for t < α , and for t > α it holds that
Hn (t) = −
2 · π
cos ωn t2 − t2 −
r2 α2
21
r2 α2
12
.
Thus the nth mode of the potential an be written as
Φn =
0
4 d
for
h
cos πz sin (2n − 1) πh sin (2n − 1) 2d 2d
2
ωn t
r2 −α 2 2
r t2 − α 2
21
12 i
for
t< t>
r α r α. (4.52)
4.2.
169
SURFACE WAVES FROM POINT SOURCES
That is a normal mode for
all n sin e the delta fun tion ontains arbitrarily high
frequen ies, ensuring that the lower limiting frequen y of ea h normal mode an be ex eeded. For simpli ity, we limit the dis ussion in the following to the potential All on lusions also hold for displa ement and pressure.
Φn .
The seismogram in
t = r/α with a singularity that is integrable. Then the 1/t, for times large ompared to r/α, while os illating. The most important feature of Φn is its frequen y modulation or dispersion. The frequen ies de rease from large values to the limiting frequen y ωn of the mode
Fig. 4.23 starts at time
amplitudes de rease with
onsidered. The dispersion in the example shown is, therefore, inverse.
Fig. 4.23: Seismogram showing frequen y modulation (dispersion).
What we have learned about the group velo ity in se tion 4.1.5 an be onrmed with (4.52). We rst ask whi h frequen ies
ω
dominate at a ertain time
t0
in
the modal seismogram. Outside the singularity, the dis ussion an be limited
f (t) near t = t0 to be able cos [f (t)] in the neighbourhood of t = t0 by the mono hromati os illation cos [ϕ0 + ω(t0 )t]. Here ϕ0 is a phase that is independent of t, and ω(t0 ) is the instantaneous angular frequen y required. This leads
to the osine fun tion in (4.52). We plan to linearise to approximate the fun tion
to
cos [f (t)] ≈ cos [f (t0 ) + f ′ (t0 )(t − t0 )] . From this, it follows that
ω(t0 ) = f ′ (t0 ).
If applied to (4.52), it follows
− 21 r2 2 ω(t0 ) = ωn t0 t0 − 2 . α
170
CHAPTER 4.
From this, we derive the quotient of frequen y
ω(t0 )
SURFACE WAVES
r/t0 , i.e., the velo ity with whi h a wave group
propagates from the sour e to the re eiver, and we get (with
ω0 = ω(t0 )) " 2 # 21 r ωn = Un (ω0 ) =α 1− t0 ω0 with
Un (ω0 )
from (4.50), i.e., exa tly the
group velo ity
of the nth mode. We,
therefore, onrm the statement from se tion 4.1.5: that ea h frequen y that is radiated from the sour e propagates to the re eiver with the group velo ity. The
omplete seismogram
in the wave guide is produ ed by onvolving the
modal seismogram (4.52) with a realisti ex itation fun tion
F (t), the spe trum
of whi h has an upper limiting frequen y, and sum. Only those normal modes (4.52) ontribute signi antly to the seismogram whi h have lower limiting frequen ies that are smaller than the upper limiting frequen y of
F (t).
Often the
response of hydrophones and seismometers, together with the dissipative me hanisms in the wave guide, redu e the number of modes. In pra tise usually only a few modes ontribute to the observed surfa e waves.
4.2.3 Computation of modal seismograms with the method of stationary phase The omputation of modal seismograms is only possible
exa tly
guides (with rigid and/or free boundaries, respe tively).
for ideal wave
In the following, an
approximation is dis ussed and demonstrated, whi h gives the modal seismogram for the fareld form of a normal mode of the type of (4.46). This is the
method of stationary phase
mentioned before.
Multiplying a normal mode in (4.46) with the spe trum fun tion
F (t)
F (ω)
of the ex itation
in (4.51), then transforming ba k into the time domain, gives
the modal seismogram as a Fourier integral. To avoid integration over negative frequen ies, we use the fa t that
real
fun tions
f (t)
annot only be represented
as
Z
+∞
1 f (t) = Re π
∞
f (t) =
1 2π
f (ω)eiωt dω
−∞
but also as
This gives the modal seismogram as
Z
0
f (ω)eiωt dω.
4.2.
171
SURFACE WAVES FROM POINT SOURCES
Φn = Re
(
) √ Z ∞ π h πh πz i 1 F (ω) i(ωt−kn r) −2 2iei 4 √ √ √ e sin (2n − 1) sin (2n − 1) dω . πd 2d 2d r 0 kn (4.53)
The approximate omputation of the integral over on the fa t that at times
t
ω
is based, as in se tion 4.1.5,
to be onsidered, the phase
ϕ(ω) = ωt − kn r is usually
rapidly varying
ompared to fun tion
tribute little to the integral in (4.53).
(4.54)
F (ω).
Su h frequen ies on
This is dierent for frequen ies with
stationary phase values. Su h a frequen y
ω0
follows from the equation
dkn =0 ϕ (ω0 ) = t − r dω ω=ω0 ′
and depends on
t.
This means that the frequen y
ω0 ,
for whi h the group
velo ity is
dω r Un (ω0 ) = = , dkn ω=ω0 t
t. From this follows the prin iple of determining the group velo ity from an observed modal seismogram. For a given time t, relative to the sour e time, the moment frequen ies and the orresponding group velo ities, using t and sour e distan e r, are determined. The sour e time and epi entre of the earthquake,
dominates the modal seismogram at time
therefore, have to be known. This gives a pie e of the group velo ity dispersion
urves. One has now to verify this pie e of the urve via forward modelling. The asso iation of a ertain frequen y to a ertain time, ne essary here, is in prin iple not unique, but the error asso iated with it an be estimated. With this method, applied to surfa e waves of earthquakes, several important results on the stru ture of the Earth were found, for example, the average rustal thi kness in dierent parts of the Earth is shown in the dierent bran hes in Fig. 4.10. A disadvantage of this method is that the result is only an average over the
whole
region between sour e and re eiver. Therefore, today several stations
are used in the interpretation of the phase velo ity ( ompare se tion 4.1.4). The
omputation of the modal seismogram requires then the following additional
steps: for given time t, we expand the phase (4.54) at the frequen y
ω0 ,
whi h
is determined by
Un (ω0 ) =
r t
(4.55)
172
CHAPTER 4.
= ϕ(ω0 ) + 21 ϕ′′ (ω0 )(ω − ω0 )2
ϕ(ω) ϕ′′ (ω0 ) An
0
=
important requirement
∞
Z
F (ω) iϕ(ω) √ e dω kn
≈ ≈
SURFACE WAVES
r
2 (ω ) Un 0
is that
dUn dω (ω0 )
ϕ′′ (ω0 ) 6= 0.
.
(4.56)
Then,
ω0 +∆ω
1 ′′ F (ω) 2 √ exp i ϕ(ω0 ) + ϕ (ω0 )(ω − ω0 ) dω 2 kn ω0 −∆ω Z F (ω0 )eiϕ(ω0 ) ω0 +∆ω i ′′ p ϕ (ω0 )(ω − ω0 )2 dω. exp 2 kn (ω0 ) ω0 −∆ω
Z
Here, we limited our dis ussion to the neighbourhood of the frequen y
ϕ(ω) is stationary. The other frequen ies 1/2 x = (ω − ω0 ) 12 ϕ′′ (ω0 ) , we get ω0 +∆ω
Z
exp
ω0−∆ω
i ′′ 2 ϕ (ω0 )(ω − ω0 ) dω 2
=
≈
with
2 ′′ ϕ (ω0 )
12 Z
+∆ω
2 ϕ′′ (ω0 )
12 Z
+∞
2π ϕ′′ (ω0 ) 1 R +∞ R +∞ π 2 2 2 cos x dx = sin x dx = . 2 −∞ −∞ =
ω0 , where
do not ontribute signi antly. With
−∆ω
ϕ′′ (ω0  2
ϕ′′ (ω0 ) 2
2
1/2
e±ix dx
−∞
12
π
e±i 4
The positive and the negative sign in the exponential term hold, if
ϕ′′ (ω0 ) >
0 and < 0, respe tively. The extension of the limits √ of the integration to x = ±∞ is possible, sin e they are proportional to r and r is very large.
Furthermore, signi ant ontributions to the integral ome only from relatively small values of x ( a.
x ≤ 5).
Putting all this together, the modal seismogram
for the ideal wave guide in the approximation given by the method of stationary ′′ phase (with ϕ (ω0 ) > 0) an be written as
Φn
= ·
√ π h πz i πh −2 2iei 4 √ sin (2n − 1) sin (2n − 1) Re 2d 2d πd ) 1 2 2π π F (ω0 )eiϕ(ω0 ) ei 4 . ′′ rkn (ω0 ) ϕ (ω0 )
2
1/2 e±ix dx
(
(4.57)
4.2.
173
SURFACE WAVES FROM POINT SOURCES
Next, one uses
ϕ(ω0 ) =
ω0 t − kn (ω0 )r, " 2 # 12 ωn ω0 1− , α ω0 " 2 # 21 ωn Un (ω0 ) = α 1 − ω0
kn (ω0 ) = r t
=
−1/2 ω0 = ω0 (t) = ωn t t2 − r2 /α2 from (4.57). After some al ulations, and for the assumption F (ω0 ) = 1, whi h orresponds to the ex itation fun tion F (t) = δ(t), the following modal seismogram and
ϕ′′ (ω0 )
from (4.56) and deletes
an be derived (please onrm)
Φn =
0
4 d
for
h
cos πz sin (2n − 1) πh 2d sin (2n − 1) 2d
ωn t2 − 2 t2 − r 2 α
12 i r2
α2
21
t
r α r α. (4.58)
We, therefore, get the stringent results of (4.52). This is surprising, onsidering the approximations used. From this we an draw the general on lusion that the method of stationary phase is a good approximation for normal modes even for more ompli ated wave guides.
dUn dω (ω0 ) = 0 and with stationary values of the group velo ity (whi h do not o
ur for ideal wave guides), the For frequen ies
ω0
with
ϕ′′ (ω0 ) = 0,
i.e., with
expansion in (4.56) has to be extended by one additional term. The treatment of the al ulations following is, therefore, slightly dierent (see, for example, Appendix E). It leads to the
behaviour of Airy phases
and shows that they
are usually the dominating parts of the modal seismograms ( ompare also Fig. 4.13).
4.2.4 Ray representation of the eld in an ideal wave guide In the last two se tions we have learned that the waveeld in an ideal wave guide is omposed of for ed normal modes. Furthermore, we found in se tion 4.1.6, that free normal modes are omposed of multiple ree ted plane body waves in the wave guide. This raises the question, an the eld of a point sour e in an ideal wave guide also be represented by the superposition of multiple ree tions? In other words, in this ase is there also a ray representation of the wave eld? In addition, it is interesting to see if mode and ray representations of the wave eld are then also equivalent.
174
CHAPTER 4.
SURFACE WAVES
We rst examine the ree tion of the spheri al wave
1 Φ0 = F R0
R0 t− α
(4.59)
at the interfa e of the wave guide in the neighbourhood of the point sour e, e.g., the free surfa e.
Q1 000000 111111 0000000000 1111111111 h 0000000000 1111111111 z=0 0000000000 r 1111111111 0000000000 1111111111 R1 h 111111 000000 000000000 111111111 0000000000 1111111111 1 0 000000000 111111111 0 1 0000000000 1111111111 Q0 000000000 111111111 0 1 R0
z=d
P(r,z)
1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 01 1 0 1 z
Fig. 4.24: Ree tion of a spheri al wave at the interfa e of the wave guide in the neighbourhood of the point sour e.
The potential of the ree tion is
−1 Φ1 = F R1 R1
is the distan e between
P (r, z)
R1 . t− α
and the
image sour e Q1 .
(4.60)
As long as the
ree tion from the lower (rigid) interfa e of the wave guide has not rea hed
Φ0 + Φ1 and, Φ0 + Φ1 satisfy, therefore, for su h times, the 2 2 surfa e z = 0 (pzz = ρ∂ (Φ0 + Φ1 )/∂t ).
the surfa e, the potential in the neighbourhood of the surfa e is therefore, zero on the surfa e.
ondition of no stress at the
Similarly, if we onsider the ree tion of the spheri al wave from interfa e
z = d,
Q0
at the
the potential of the wave ree ted there an be written as
1 F Φ2 = R2
R2 , t− α
(4.61)
R2 now has to be determined for a new image sour e with the z  oordinate d + (d − h) = 2d − h. That Φ0 + Φ2 satises the boundary onditions ∂(Φ0 + Φ2 )/∂z = 0 for z = d (zero normal displa ement), an be seen easily, sin e for where
points in that interfa e
4.2.
175
SURFACE WAVES FROM POINT SOURCES
R0
= R2
∂R0 ∂z
=
∂R2 ∂z
=
z − h R0
d−h R0 z=d d−h ∂R0 2d − h − z =− =− − . R2 R2 ∂z z=d =
With the two previously onsidered ree tions of the spheri al wave originating from
Q0
at the interfa es of the wave guide, boundary onditions an only be
satised for ertain times, e.g., only as long as the ree tions
Φ1
and
Φ2
have
rea hed the opposite interfa e, respe tively. Sin e they are of the same form as
Φ0 ,
higher order ree tions an be onstru ted in the same way. Ea h ree tion
and multiple ree tion seems to ome from an image sour e, whi h was reated by the appli ation of multiple mirror images of
Q0
at the interfa es (Fig. 4.25).
The sign of the orresponding potential is negative if the number of ree tions at the surfa e is odd, otherwise it is positive. Ea h image sour e orresponds to a
ray
from the sour e to the re eiver whi h has undergone a ertain number of
ree tions.
Fig. 4.25: Image sour es for multiple ree tions in the wave guide.
176
CHAPTER 4.
The ray on ept an, without di ulties, be generalised for
SURFACE WAVES
solid
media (in lud
ing Swaves), but this is not true for the on ept of the image sour e. This is why, in general, and also in the ase presented here, we speak of a
tation of the wave eld.
Φ
ray represen
It an be expressed as
∞ X
1 Rj1 Rj2 1 + F t− F t− (−1) − Rj1 α Rj2 α j=0 1 1 Rj3 Rj4 + − F t− F t− Rj3 α Rj4 α
=
j
(4.62)
with
2
2 Rj1
=
(2jd + h + z) + r2
2 Rj2
=
2 Rj3
(2jd − h + z)2 + r2
=
2 Rj4
=
2
(2(j + 1)d − h − z) + r2 2
(2(j + 1)d + h − z) + r2 .
Only those terms in (4.62) are nonzero, for whi h the argument is positive, and for whi h
F (t)
at
t = 0
is not zero.
The number of su h terms is nite and
in reases with time. For the ideal wave guide, the determination of the ontribution of a ray is simple, sin e it follows the same time law as the exiting spheri al wave. For other wave guides, methods like those presented in se tion 3.8 have to be used. The resulting numeri al eort is then signi antly greater and seems only justied if not too many rays have to be summed up, but that is ne essary for large horizontal distan es from the sour e where the paths of many rays be ome very similar. In this ase, the representation of the wave eld as a sum of only a few normal modes is
signi antly more e ient.
The ray representation is most suited for
su h distan es from the sour e where the typi al normal mode properties of the wave eld have
not yet
developed.
Finally, we will show that (4.62) for tively, are
F (t) = eiωt
and (4.39) and (4.40), respe 
two dierent representations of the same wave eld.
dis ussion rst to the ase
h ≤ z ≤ d.
If
F (t) = eiωt
We limit our
is inserted into (4.62),
and the Sommerfeld integral (3.85) for a spheri al wave is used for ea h term, iωt it follows (the fa tor e is again omitted)
Φ
=
Z
0
∞
J0 (kr)
∞ kX (−1)j [− exp (−il 2jd + h + z) il j=0
4.2.
177
SURFACE WAVES FROM POINT SOURCES
+ exp (−il 2jd − h + z)
+ exp (−il 2(j + 1)d − h − z) − exp (−il 2(j + 1)d + h − z)] dk. z ≥ h, the ontributions are exp(−i2ljd) an be separated If
Φ
=
Z
0
∞
equal to the arguments
everywhere.
Then
∞ X k j J0 (kr) −e−2ild [− exp (−il (h + z)) il j=0
+ exp (−il (−h + z))
+ exp (−il (2d − h − z)) − exp (il (2d + h − z))] dk. The expansion in the rst square bra ket has a sum of −ild
an be extra ted giving se ond square bra ket, e
Φ
=
1/(1 + e−2ild ).
From the
∞
k [− exp (il (d − h − z)) 2il cos(ld) 0 + exp (il (d + h − z)) + exp (−il (d − h − z)) Z
J0 (kr)
− exp (−il (d + h − z))] dk.
The remaining square bra ket is equal to
−2i sin [l (d − h − z)] + 2i sin [l (d + h − z)] = 4i cos [l (d − z)] sin (lh) . Thus,
Φ=2
Z
0
∞
J0 (kr)
k cos [l(d − z)] sin(lh) dk, l cos(ld)
(4.63)
whi h agrees with (4.40). If sour e and re eiver are ex hanged in (4.62), the potential of a single ray is un hanged sin e it depends only on the path travelled. Therefore, ex hanging
z
with
h
in (4.63) gives the potential for
Thus, the proof of the identity of (4.62) (for
0 ≤ z ≤ h whi h leads to (4.39). F (t) = eiωt ) with (4.39) and (4.40),
respe tively, is omplete. Finally, we would like to reiterate ( ompare se tion 4.2.1) that a representation of the wave eld
by normal modes alone
is only possible for ideal wave guides
178
CHAPTER 4.
SURFACE WAVES
whi h have upper and lower boundaries that are ompletely ree ting for all angles of in iden e. In other media, additional ontributions (body waves, leaky modes) o
ur whi h are not due to the poles in the omplex plane like the normal modes.
Exer ise 4.7 Study the polarisation of the displa ement ve tor of the se ond free normal
n
mode of the ideal wave guide ( =2 in (4.48)) as a fun tion of depth.
Exer ise 4.8 An explosive point sour e is lo ated at depth halfspa e. The displa ement potential
Φ
h below the free surfa e of a liquid
is the sum of the potentials (4.59) of
the dire t wave and (4.60) for the ree tion. Give an approximation for
Φ whi h
holds under the following onditions (dipole approximation) : a) The dominant period of the ex itation fun tion travel time
h/α
b) The distan e
F (t)
is mu h larger than the
from the sour e to the surfa e.
r
to the re eiver is mu h larger than
Introdu e spheri al oordinates
R
and
ϑ
h.
relative to the point
( ompare Fig. 4.24). What is the result, if the surfa e is not free but rigid?
r = 0,
and
z=0
Appendix A
Lapla e transform and delta fun tion A.1 Introdu tion to the Lapla e transform A.1.1 Literature Spiegel, M.R.
: Lapla e Transformation, S haum, New York, 1977
Riley, K.F., Hobson, M.P. and Ben e, J.C.
: Mathemati al methods for
physi s and engineering, A omprehensive guide, Cambridge University Press, Cambridge, 2nd edition, 2002
A.1.2 Denition of the Lapla e transform The Lapla e transform asso iates a fun tion transforms a fun tion
f (s) =
Z
0
F (t)
=
f (s) =
∞
F (t)
into the fun tion
f (s) with f (s).
the fun tion
F (t),
or it
e−st F (t)dt = L {F (t)}
original fun tion image fun tion
Symboli notation:
(Lapla e − transform,
f (s) •−◦ F (t) (•−◦ = symbol t = s =
real variable
σ + iω
(of
abbreviated L
of asso iation) time)
omplex variable
179
− transform)
180
APPENDIX A.
LAPLACE TRANSFORM AND DELTA FUNCTION
A.1.3 Assumptions on F(t) 1. F(t) is usually a real fun tion 2.
F (t) ≡ 0
3.
F (t)
for
t < 0 (satised
for many physi al parameters  ausality)
should be integrable in the interval
[0, T ],
and for
t> T
it should
hold that
F (t) < eγt
with real
γ.
f (s) of F (t) s > γ ( onvergen e halfplane ). All limited fun tions (α > 0), sin βt et . have an Ltransform but also nonlimited −1/2 n , t and eαt (n, α > 0). Note assumption 2. Many fun tions fun tions as t −1 t2 in physi s also have an Ltransform. The fun tions t and e do not have an These are su ient onditions for the existen e of the Ltransform
s
for omplex −αt as, e.g., e
with Re
Ltransform.
A.1.4 Examples a)
0 for t < 0 H(t) = Heaviside step fun tion (unit step) 1 for t ≥ 0 Z ∞ ∞ 1 1 f (s) = e−st dt = − e−st 0 = for Re s > 0 ( onverg. halfplane) s s 0 1 H(t) ◦−• s F (t)
=
b)
F (t) f (s)
= =
Z
0
eδt · H(t) For
δ=0
◦−•
0 eδt ∞
for for
t δ
H(t)
)
sin at 1 . · H(t) ◦−• 2 a s + a2 Tables of many more orresponden es an be found in the literature given.
A.1.
181
INTRODUCTION TO THE LAPLACE TRANSFORM
A.1.5 Properties of the Lapla e transform Similarity theorem a > 0 F (at) ◦−•
Z
∞
e
−st
F (at)dt =
0
Therefore,
Z
∞
e
s −a at
0
d(at) 1 F (at) = a a
Z
∞
0
1 s F (at) ◦−• f . a a
Thus, only the Ltransform of
F (t)
s
e− a τ F (τ )dτ
(A.1)
has to be known.
Example: A
ording to se tion A.1.4
et · H(t) ◦−•
1 . s−1
With the similarity theorem, it follows that
eat · H(t) ◦−•
1 a
s a
1 1 = , −1 s−a
i.e., the result of the dire t omputation in se tion A.1.4.
Displa ement theorem
Fig. A.1: Displa ement theorem.
F (t − ϑ) ◦−•
Z
F (t − ϑ) ◦−• e
∞
0 −ϑs
e
−st
F (t − ϑ)dt =
f (s)
Z
∞
e−s(τ +ϑ) F (τ )dτ = e−ϑs f (s)
0 (A.2)
Damping theorem (α arbitrary omplex) e
−αt
F (t) ◦−•
Z
∞
e−(s+α)t F (t)dt = f (s + α)
0
e−αt F (t) ◦−• f (s + α)
(A.3)
182
APPENDIX A.
LAPLACE TRANSFORM AND DELTA FUNCTION
Dierentiation theorem F ′ (t) ◦−•
∞
Z
e−st F ′ (t)dt = e−st F (t)∞ 0 +s
0
Z
∞
e−st F (t)dt
0
The rst term is zero at its upper limit due to the assumption 3 from hapter A.1.3. Thus,
F ′ (t) ◦−• sf (s) − F (+0)
F (+0)
=
lim t→0 t>0
right
F(t) is the limit from the
side.
(A.4)
Generalisation:
F ′ (t) F ′′ (t) . . .
F
(n)
◦−• sf (s) − F (+0) ◦−• s2 f (s) − sF (+0) − F ′ (+0)
(A.5)
(t) ◦−• s f (s) − s F (+0) − s F (+0) (n−2) (n−1) − . . . − sF (+0) − F (+0) n
n−1
n−2
′
Integration theorem G(t) =
Z
0
t
1 F (τ )dτ ◦−• f (s) s
(A.6)
From this, it follows that
G′ (t) = F (t) ◦−• f (s) − G(+0) = f (s).
Convolution theorem Z
0
t
F1 (τ )F2 (t − τ )dτ ◦−• f1 (s)f2 (s)
The integral is alled onvolution of Furthermore, it holds that
F1
with
F2 ,
symboli notation
(A.7)
F1 ∗ F2 .
A.1.
183
INTRODUCTION TO THE LAPLACE TRANSFORM
Z
0
t
F1 (τ )F2 (t − τ )dτ =
Z
t
0
F1 (t − τ )F2 (τ )dτ
or
F1 ∗ F2 = F2 ∗ F1 , i.e., the onvolution is ommutative. Further elementary properties of the Ltransform are that it is rstly homogeneous and linear, i.e., it holds that
a1 F1 (t) + a2 F2 (t) ◦−• a1 f1 (s) + a2 f2 (s), and se ondly, that from
F (t) ≡ 0
it follows that
f (s) ≡ 0
and vi e versa.
An important property, whi h follows from the denition of the Ltransform, is
lim
Re s→+∞ Only
then
is a fun tion
f (s)
f (s) = 0.
(A.8)
an Ltransform and an be transformed ba k (see
next hapter).
A.1.6 Ba ktransform −1
F (t) = L
1 {f (s)} = 2πi
Z
c+i∞
ets f (s)ds
(A.9)
c−i∞
Fig. A.2: Convergen e halfplane of Lapla e inversetransform. The integration path is parallel to the imaginary axis and has to be situated in the onvergen e halfplane of the integration path, the left.
f (s)
f (s),
otherwise,
is arbitrary. To the right of
annot have singularities, but it an have them to
The integration path an be deformed in a
ordan e with Cau hy's
integral law and the remainder theorem.
184
APPENDIX A.
LAPLACE TRANSFORM AND DELTA FUNCTION
A.1.7 Relation with the Fourier transform A ommon representation of the Fourier transform
F (ω) =
+∞
Z
F (ω)
of a fun tion
F (t)e−iωt dt.
F (t)
is
(A.10)
−∞
F (ω)
is also alled the
omplex spe trum
of
F (t); ω
is the angular frequen y.
The inversetransform is given by
1 F (t) = 2π
Z
+∞
F (ω)eiωt dω.
(A.11)
−∞
This equation an be interpreted as the superposition of harmoni os illations. This is the reason why the Fourier transform is often used in physi s. If the behaviour of a system, whi h an be des ribed by linear dierential equations, is known for harmoni ex itation, its behaviour for impulsive ex itation an be determined via (A.11).
To do this, the ex itation has to be broken into its
spe tral omponents a
ording to (A.10). Then, the problem is solved for ea h spe tral omponent, and, nally, all spe tral solutions are superimposed via (A.11). Su h an approa h is used in se tion 3.6.3 in the study of the ree tion of impulsive waves at an interfa e. It is often less physi al, but often more elegant and simple, to use the Ltransform. The onne tion between
F (t)
that are zero for
t r1 follows from τ = 0. Then The retardation the sphere
Z τ r1 1 1 − rα (τ −ϑ) 1 dϑ . U (r, t) = U1 (τ ) + α − U1 (ϑ)e r r r1 0
A.3 The delta fun tion δ(t) A.3.1 Introdu tion of δ(t) We examine the result (A.17) for the me hani al resonator
A.3.
δ(T )
THE DELTA FUNCTION
Y (t) =
for the for e
1 ωm
K(t) = Iδǫ (t),
Z
0
t
195
K(τ )e−αω0 (t−τ ) sin ω(t − τ )dτ,
where
I
is a onstant with the dimension of for e
time (=dimension of an impulse) and
δǫ (t) =
0
1 ǫ
0
for for for
×
t