Troyan V., Kiselev Y. Statistical Methods of Geophysical Data Processing
Подождите немного. Документ загружается.
The relationships of sounding signal and medium parameters 133
(∆l) in a direction n is connected, at the first, with an absorption, which one is
proportional to an intensity (ϕ), secondly, with a dissipation of the radiation ϕ
from a direction n on all by another directions (in Fig. 4.7 the examples of angular
distributions of the dissipation of the radiation intensity by a small (a), large (b)
and super-large (c) particles) are given:
∼
1
4π
Z
q(n
0
← n)ϕ(x, n)dn
0
= ϕ(x, n)
1
4π
Z
q(n
0
← n)dn
0
,
Thirdly, with the contribution to the intensity of a flux with a direction of the
radiation n, dissipated in all other directions:
∼
1
4π
Z
q(n ← n
0
)ϕ(x, n
0
)dn
0
,
Fourthly, with a radiation source. Therefore the transport equation takes a form
Fig. 4.7 Examples of angular distributions of a dissipation of a radiation intensity by small (a),
large (b) and super-large (c) particles.
(with taking into account of ∇
n
= (n,
∇
∇) and α
s
+ α
a
∆
= α)
(n ·∇
∇
∇)ϕ(x, n) + αϕ(x, n) −
1
4π
Z
q(x, n, n
0
)ϕ(x, n
0
)dn
0
= s(x, n). (4.72)
And in the operator form
Lϕ = s, L = (n ·
∇
∇) + α
ˆ
I −
1
4π
< q(n, n
0
)|, (4.73)
where here symbol < q(n, n
0
)| means
R
q(·)dn
0
; ϕ = ϕ(x, n) is a field of the sounding
signal; s = s(x, n) is a signal source. Substituting the operator L from the equation
(4.73) as L = L
0
− S (S is a scattering operator),
S = α
ˆ
I +
1
4π
< q(n, n
0
)|,
we write down the solution of the transport equation as an expansion on the multiple
scattering:
ϕ = ϕ
0
+ L
−1
0
SL
−1
0
s + L
−1
0
SL
−1
0
SL
−1
0
s + . . . ,
where
L
−1
0
s = ϕ
0
, ϕ = ϕ
0
exp
−
Z
α(x)dl
.
Let’s note, that the equation of the carrying (4.72) is used for the solution of
the wide range of problems, for example, at studying of the transport of monochro-
matic neutrons (in this case it is accepted to call it as the Boltzmann equation), in
problems of X-ray and spectral tomography, at laser sounding etc.
This page intentionally left blankThis page intentionally left blank
Chapter 5
Ray theory of wave field propagation
The exact solutions of problems of sounding signals propagation are constructed
only for limited number of media models. As a rule, this set of media includes
uniform medium, layered homogeneous medium, uniform medium with the inclu-
sions of high symmetry. For the interpretation of real geophysical fields (seismic,
acoustic, electromagnetic) it is necessary to construct the approximate solutions for
the wave propagation in non-uniform media. So, it exists in the Earth not only the
interfaces, on which one the elastic properties vary by jump, but also areas, inside
which there is a smoothly varying variation of elastic properties. From the physical
point of view the ray theory is interpreted as follows: the waves propagate with
local velocities along ray pathways and arrive in the observation point with ampli-
tudes described by a geometrical spreading of rays from a source to the receiver
point. At an enunciating of this chapter we shall follow the description introduced
in (Ryzhikov and Troyan, 1994).
5.1 Basis of the Ray Theory
One of the most common method of the solution of the equations of wave propa-
gation is the method of geometrical optics (Babic and Buldyrev, 1991; Babic et al.,
1999; Kravtsov, 2005; Bleistein et al., 2000). This method is a shortwave asymp-
totic of a field in weak non-uniform, slow non-stationary and weak-conservative
media: the sizes of the inhomogeneity are much greater than the wavelength and
time intervals of the non-stationarity much more than the period of oscillation. The
shortwave asymptotic allows to consider the medium locally as homogeneous and
stationary and is based on the assumption of a wave field in a form (ϕ) a “quick”
phase and a “slow” amplitude multiplier factors.
Let us consider a formal scheme of the space-time ray method. Let a field ϕ sat-
isfies homogeneous linear equations, which we shall write down, following (Kravtsov
and Apresyan, 1996), in a form of the integral equation:
Lϕ =
Z
L(x, x
0
)ϕ(x
0
)d
4
x
0
= 0 (5.1)
135
136 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING
(in this section we put x = (x
1
, x
2
, x
3
, t)). Let us assume that the operator L
contains a part L
0
, carrying the responsibility for the description of the field ϕ in
the conservative medium, and the part L
0
with a small parameter α:
L = L
0
+ αL
0
. (5.2)
We represent the field ϕ as an asymptotic expansion (the Debye expansion):
ϕ = A(y, α) exp
i
α
τ(y)
, (5.3)
where y = αx is “ a slow” argument; τ is a phase (or eikonal); A is a wave amplitude.
The expansion of A over power of α has a form
A = A
0
(y) + αA
1
(y) + O(α
2
).
Let L is a near differential operator, i.e. L(x, x
0
) 6= 0 under small ∆x = x
0
− x.
Then the expansion of the phase function exp{(i/α)τ(αx)} in the point x
0
accurate
to O(α
2
) can be represented as:
exp
i
α
τ(α(x + ∆x))
= exp
i
α
τ(αx
0
) + α(∆x · ∂
y
)τ +
1
2
α
2
(∆x · ∂
y
)
2
τ + O(α
2
)
= exp
i
α
τ(y)
× exp{i(p · ∆x)}
1 +
i
2
α(∆x · ∂
y
)
2
τ
+ O(α
2
),
where p
∆
= ∇∇
y
τ. The expansion of the amplitude accurate to O(α
2
) has a form:
A(αx
0
) = A(α(x + ∆x)) = (1 + α(∆x · ∂
z
))A(z),
where z = αx. Finally, the wave field in a point x
0
is given by the formula
ϕ(x
0
) = exp
i
α
τ(y)
exp{i(p · ∆x)}
×
1 + α
(∆x · ∂
z
) +
i
2
(∆x · ∂
y
)
2
τ(y)
× A(z)
z=y
+ O(α
2
). (5.4)
Regarding this expansion and keeping linear terms on α, we can rewrite (5.1) on
the next form:
Z
L(x, x
0
)ϕ(x
0
)d
4
x
0
= exp
i
α
τ(αx)
×
Z
d
4
x
0
L(x, x
0
) exp{ip(x
0
− x)}
×
1 + α)
(∆x · ∂
z
) +
i
2
(∆x · ∂
y
)
2
τ(y)
A(z)
= exp
i
α
τ(αx)
1 − α
i(∂
p
∂
z
) +
i
2
(∂
p
∂
y
)
2
τ(y)
×
Z
d
4
x
0
L(x, x
0
) exp[ip(x
0
− x)]A(z).
Ray theory of wave field propagation 137
Introducing the denotation
R
d
4
x
0
L(x, x
0
) exp[ip(x
0
− x)] ≡
ˆ
L(p), let us represent
(5.1) as
1 − iα
(∂
p
∂
z
) +
1
2
(∂
p
∂
y
)
2
τ(y)
y=αx
ˆ
L(p)A(z)
z=αx
= 0.
Finally, expanding the amplitude in the expression for A(z) in a point x by power
α and representing L(p) in the form L = L
0
+ αL
0
(and equate the terms of the
equal order on α), we obtain the system of recurrence equations to determine the
amplitude (now we put the parameter α being equal an identity element):
ˆ
L
0
(p)A
0
(x) = 0, (5.5)
ˆ
L
0
(p)A
1
(x) =
i
(∂
p
∂
x
) +
1
2
(∂
p
∂
y
)
2
τ(y)
y=x
ˆ
L
0
(p) −
ˆ
L
0
(p)
A
0
(x) (5.6)
etc.
If the field is multicomponent, for example, the displacement field in the Lame
equation, electromagnetic field, then the amplitude A is a vector function, and the
operator
ˆ
L(p) is a matrix operator. The resolvability condition (5.5) is reduced to
the requirement of the equality to zero of the determinant:
det L
0
(p) = 0. (5.7)
If L
0
corresponds to the conservative medium, then the operator
ˆ
L
0
is the Hermitian
operator and the next canonical representation of the matrix L
0
takes place:
ˆ
L
0
=
X
i
λ
i
e
i
e
T
i
, (5.8)
where {λ
i
} are eigenvalues; {e
i
} are orthonormal eigenvectors of the matrix L
0
, i.e.
L
0
e
i
= λ
i
e
i
, e
T
i
e
j
= δ
ij
.
The resolvability condition (5.7) subject to a canonical form (5.8) gives the
dispersion equations:
λ
i
= λ
i
(x, p) = 0, i = 1, 2, . . . , (5.9)
which determine the relationship between space wave vectors ∂τ/∂x
1
, ∂τ /∂x
2
,
∂τ/∂x
3
) and the frequency ∂τ/∂x
0
∆
= ∂τ/∂t, as well as the phase (eikonal) τ space
and time behavior. Because p = ∂τ/∂x, then the dispersion equation is the first-
order differential equation, which is usually solved by the method of characteristics.
It is usual to write the equation for the characteristic in the Hamilton form, intro-
ducing the parameter l:
d
l
x = ∂
p
λ
i
,
d
l
p = −∂
x
λ
i
. (5.10)
The equation (5.10) describes the rays corresponding to the different types of waves,
while i takes the different values.
138 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING
5.2 Ray Method for the Scalar Wave Equation
The scalar wave equation represents the mathematical model underlying the de-
scription of the set of the physical processes of the signal propagation. Therefore
enunciating of the concrete examples to use of the space-time ray method (Babic
et al., 1999) we shall begin from the wave equation
ϕ = ∆ϕ −
1
c
2
∂
2
∂t
2
ϕ = 0. (5.11)
In accordance with the statement form of the operator equation of the sounding
signal propagation Lϕ = s, in this case L is
L = ∆ −
1
c
2
∂
2
∂t
2
; ϕ = ϕ(x, t); s = 0
(we consider the wave field in an area without any sources).
Let us write the general form of an asymptotic expansion (5.3), extract the time
in the explicit form, as stated below
ϕ(x, t) = A(α, αx, αt) exp[i/ατ(αx, αt)].
Let’s note, that the wave operator is a local operator (differential opera-
tor), i.e. in representation of the operator L as an integral operator (Lϕ =
R
L(x, t; x
0
, t
0
)ϕ(x
0
, t
0
)dx
0
dt
0
), operator kernel is singular, i.e. the sign of integra-
tion get lost.
Approximation of the phase part of the wave field ϕ in the vicinity of x
0
, t
0
(expansion (5.4)) accurate to O(α) looks like
ϕ(x, t) ∼ exp[i(p∆x − p
0
∆t)],
where p =
1
α
∇
∇τ
x
0
,t
0
, p
0
=
1
α
∂
∂t
τ
x
0
,t
0
,
∆t = t − t
0
, ∆x = x − x
0
,
i.e. the approximation in the vicinity of x
0
, t
0
appears to be a plane wave. The
resolvability equation (5.5) implies that
(p, p) −
1
c
2
(x, t)
p
2
0
= 0 (5.12)
(here L
0
= ∆ − [c
2
(x, t)]
−1
∂
2
/∂t
2
), therefore
L
0
(p) = exp[−i(px − p
0
t)]
∆ −
1
c
2
∂
2
∂t
2
exp[i(px − p
0
t)].
The resolvability condition (5.12), which we have written in the form of the dis-
persion equation, gives a connection of the frequency p
0
∆
= ω and the wave vector
p
∆
= k (spatial frequency):
ω = ω(k) = c|k|, (5.13)
where c = c(x, t); ω is the solution of the dispersion equation |L
0
(ω(k))| = 0.
Ray theory of wave field propagation 139
The ray approximation of the solution of (5.11) could be represented in the next
form:
ϕ(k, ω) = A(k)δ(L
0
(ω(k))) = A
k
δ(ω − ω(k)), (5.14)
where A
k
∆
= A|∂ωL(k, ω)|
−1
(a consequence of the formula δ(f (x)) = |∂
x
f|
−1
×
×δ(x − x
0
), where x
0
is a root of the equation f(x) = 0 ). The presence in the
representation of (5.14) the Dirac delta function δ(|L
0
|) means, that the amplitude
of a field ϕ can be differ from zero only for wave vectors k, fitting to the dispersion
equation, i.e. a constraint equation for the frequency and the modulus of the wave
vector (5.13).
Let us note, that the expression (5.14) corresponds to the single type of
wave (single-mode propagation in conservative medium), that means the equation
|L
0
(ω(k))| = 0 has a one real root.
In the general case of the multimode propagation, the wave field ϕ can be rep-
resented by the superposition of the several types of waves:
ϕ =
X
µ
A
k
δ(ω − ω
µ
(k)).
In the space-time representation the solution ϕ is described as a set of the travelling
waves:
ϕ(x, t) =
1
(2π)
3/2
Z
A
k
exp[i(kx − ω(k)t)]dk.
Writing the inverse Fourier transform of ϕ(x, t), it is possible to find out that
A
k
e
iω(k)t
=
1
(2π)
3/2
Z
ϕ(x, t)e
−ikx
dx,
i.e. the space Fourier operator of the free wave field is an oscillating one.
The resolvability equation (5.5) written for a phase τ is called a general eikonal
equation:
(∇
x
τ, ∇
x
τ) =
1
c
2
(x, t)
∂
∂t
τ
2
. (5.15)
This equation is a generalisation of the eikonal equation for the known case of the
stationary media (c = c(x), τ = τ − t, c
−2
(x)
∆
= n
2
(x) plays usually the role of the
refractive index).
We may name the equation (5.15) as a characteristic equation of the wave equa-
tion (5.11). The hypersurfaces τ(x, t) = const are the characteristics of the wave
equation. The equation (5.15) can be represented as the Hamilton–Jacobi equation:
∂τ
∂t
+ H(
∇
∇
x
τ, x, t) = 0, (5.16)
where H = c(x, t)|
∇
∇τ|. The primary method to solve this equation is the charac-
teristics method. Let us remind the algorithm of this method for the general form
of the equation (5.16):
H
x
µ
,
∂τ
∂x
µ
= 0, µ = 0, . . . , M, (5.17)
140 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING
i.e. we consider the general partial equation of the second order (index µ = 0 is
connected with the time parameter). Let supplement (5.17) by the initial conditions
defined on hypersurface Q, we then start to solve the Cauchy problem.
We consider Q as a parameter defined value, i.e.
Q = {x
µ
: x
µ
0
(γ
1
, . . . , γ
m
)}, γ ∈ Γ,
where the vectors {∂x/∂γ}are linear independent (surface Q is regular, im kleinen).
Let the determinant D is nonsingular on the surface Q:
D =
∂H/∂τ
0
. . . ∂H/∂τ
m
∂x
0
/∂γ
1
. . . ∂x
m
/∂γ
1
. . . . . . . . .
∂x
0
/∂γ
m
. . . ∂x
m
/∂γ
m
6= 0, (5.18)
where
τ
µ
∆
=
∂
∂γ
µ
τ.
To set the initial conditions on Q as
τ
Q
= τ
0
(γ),
∂τ
∂x
µ
Q
= τ
0
µ
(γ).
The characteristic system of equations for the equation (5.17) is the system of equa-
tions
dx
µ
/ds = ∂H/∂τ
µ
,
dτ
µ
/ds = −∂H/∂x
µ
,
dτ/ds = −
X
µ
τ
µ
∂H/∂τ
µ
, (5.19)
where
s :
∂(x
0
, . . . , x
m
)
∂(s, γ
1
. . . , γ
m
)
s=0
6= 0
owing to (5.18), i.e. locally (in the vicinity of s) the coordinate frame (s, γ
1
, . . . , γ
m
)
is non-degenerate.
Let us present the method of derivation of the system of characteristic equations.
Let us consider (5.17) as a constraint equation of the independent variables {x
µ
}
and {p
µ
∆
= ∂τ/∂x
µ
}, i.e. let’s find the solution in phase space X × P, where
X = {x = (x
0
, x
1
, . . . , x
m
)}, P = {p = (p
0
, p
1
, . . . , p
m
)}, H(x
µ
, p
µ
) = 0. Owing to
this constraint, the differential H must be identically equal to zero:
dH =
X
µ
∂H
∂x
µ
dx
µ
+
∂H
∂p
µ
dp
µ
= 0.
Hence, it is obviously, the formal vectors ∂H/∂x, ∂H/∂p and correspondingly vec-
tors dx, dp are orthogonal. For the orthogonality condition to be certainly satisfied,
it is necessary to put the first components of dx, dp to be proportional to the second
Ray theory of wave field propagation 141
components of ∂H/∂x, ∂H/∂p, whereas the relation of other components has an
opposite sign of the constant of proportionality.
Introducing the formal parameter s ( i.e. we let p = p(s), x = = x(s)), we can
write
dx
µ
/ds = ∂H/∂p
µ
,
dp
µ
/ds = −∂H/∂x
µ
.
To exclude the parameter s we should write the third equation, which determines
dτ/ds under the condition p = ∂τ/∂x:
dτ
ds
=
X
µ
∂τ
∂x
µ
dx
µ
ds
=
X
µ
p
µ
∂H
∂p
µ
.
Solving the system (5.19) under initial conditions, we obtain x
µ
= x
µ
(s, γ
1
, . . . γ
m
),
τ
µ
= τ
µ
(s, γ
1
, . . . γ
m
), τ = τ(s, γ
1
, . . . γ
m
). Using the initial coordinates
(x
0
, . . . , x
m
), we find the functions τ = τ (x) and τ
µ
= τ
µ
(x), that are the solu-
tions of the Cauchy problem, which satisfies both (5.17) and initial conditions.
Returning to the equation (5.16), let us write the characteristic equations (x
0
≡
t):
dt
ds
= 1,
dx
1
ds
=
∂H
∂τ
1
. . .
dx
m
ds
=
∂H
∂τ
m
,
dτ
1
ds
= −
∂H
∂x
1
. . .
dτ
m
ds
= −
∂H
∂x
m
,
dτ
ds
= τ
0
+
m
X
µ=1
ˆτ
µ
∂H
∂τ
µ
. (5.20)
Here, in the last equation the addend τ
0
= ∂τ/∂t is emphasis ed and we take into
account that ∂H/∂τ
0
= 1.
Let us remind an analogy between the ray in optics and the trajectory in me-
chanics. Writing the eikonal equations as
H(x, p) =
1
2
(p, p) − n
2
(x)
= 0,
we can obtain the equations for the rays in the Hamilton form:
dx/ds = p,
dp
ds
=
1
2
∇n
2
(x).
Following to ∂τ /∂s = (p, p) we have the next equation for the eikonal τ
τ = τ
0
+
s
Z
s
0
p
2
ds = τ
0
+
s
Z
s
0
n
2
(x(s))ds.
The parameter s is connected with an arc length of the ray (dl) by virtue of
(dl)
2
= (dx)
2
=
∂H
∂p
2
ds
2
142 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING
(in the last equality we take into account the first equations of the characteristic
system). In this case
ds = dl|∂H/∂p|
−1
= dl|p|
−1
= dl/n(x),
and τ depends on arc length l as
τ = τ
0
+
B
Z
A
n(x)dl.
The system of Hamilton equations relatively to x and p, excluding p, can be
written as
d
2
x
ds
2
=
1
2
∇n
2
(x).
Let us compare this equation with the dynamic one for the motion of a particle
in the potential field (F = −∇ϕ):
d
2
x/dt
2
= −
∇
∇ϕ.
This analogy allows the results obtained in optics, to use in mechanics, and on the
contrary. Essential difference of posing of mechanical and wave problems is the
assignment of the initial data: by virtue of the specificity of the wave problems the
initial conditions determine not an alone trajectory, as it is accepted in mechanics,
but a set (congruence) of “trajectories” — rays.
For wave problems are of interest (exotic from the point of view of mechanics)
the next example
H(kτ
i
, t, x) = kH(τ
i
, t, x),
where i = 1 ÷ m, k > 0, i.e. the last equation of the system (5.20) has a form
dτ/ds = τ
0
+ H = ∂τ/∂t + H ≡ 0,
at that τ = const along the curve, that is the solution of the system (5.20) (loss of
the dispersion).
Let’s note, that the solution of the characteristic system (5.19) as the system
of the ordinary differential equations with given initial conditions, by virtue of the
uniqueness theorem represents a family of curves without intersections in the phase
space.
It is accepted in geometrical optics to call a spatially-time projection of the char-
acteristic of the eikonal equation as the spatially-time ray, x
µ
= x
µ
(s, γ
1
, . . . , γ
m
),
which one can be obtained or as the solution of the system of equations (5.19), or
as an extreme value of the functional
τ =
B
Z
A
X
µ
τ
µ
∂H
∂τ
ds
on every possible pathways pairing a dots A and B. (In Fig. 5.1 the projection
of a characteristic, given in a phase space, on the physical space x
1
x
2
with the