Parabolic and white-noise approximations for elastic waves in random media

Wave Motion 46 (2009) 237–254

Contents lists available at ScienceDirect

Wave Motion journal homepage: www.elsevier.com/locate/wavemoti

Parabolic and white-noise approximations for elastic waves in random media Josselin Garnier a,*, Knut Sølna b a b

Laboratoire de Probabilités et Modèles Aléatoires & Laboratoire Jacques-Louis Lions, Université Paris VII, 2 Place Jussieu, 75251 Paris Cedex 05, France Department of Mathematics, University of California, Irvine, CA 92697, USA

a r t i c l e

i n f o

Article history: Received 31 July 2008 Received in revised form 5 February 2009 Accepted 12 February 2009 Available online 25 February 2009

Keywords: Elasticity Paraxial Wave Averaging Stochastic flow

a b s t r a c t In this paper we carry out an asymptotic analysis of the elastic wave equations in random media in the parabolic white-noise regime. In this regime, the propagation distance is much larger than the initial beam width, which is itself much larger than the typical wavelength; moreover, the correlation length of the random medium is of the same order as the initial beam width, and the amplitude of the random fluctuations is small. In this distinguished limit we show that wave propagation is governed by a system of random paraxial wave equations. The equations for the shear waves and pressure waves have the form of Schrödinger equations driven by two correlated Brownian fields. The diffraction operators can be expressed in terms of transverse Laplacians. The covariance structure of the Brownian fields is determined by the two-point statistics of the density and Lamé parameters of the random medium. Ó 2009 Elsevier B.V. All rights reserved.

1. Introduction The parabolic approximation for wave propagation in heterogeneous media is a model used in a vast number of applications, for instance in communication and imaging [27]. It usually has the form of a Schrödinger equation that describes waves propagating along a privileged propagation axis determined by the source. It is very simple compared to the full three-dimensional wave equation, both from the theoretical and numerical viewpoints, and it enables the analysis of many important phenomena, such as laser beam propagation [10,26], time reversal in random media [3,9,22,23], underwater acoustics [27] or migration problems in geophysics [4]. The parabolic approximation for scalar waves in homogeneous and random media is rather well understood. It is even valid in some scaling regimes in which the medium fluctuations are rapid and can be approximated by a white-noise term. The main motivations for studying the white-noise paraxial wave equation are (i) it appears as a very natural model in many situations where the correlation length of the medium is relatively small, in particular much smaller than the propagation distance, (ii) it allows for the use of Itô’s stochastic calculus, which in turn enables the closure of the hierarchy of moment equations and the statistical analysis of important wave propagation problems, such as scintillation [12]. When the paraxial approximation and the white-noise approximation can be justified simultaneously for scalar waves, then the limit equation takes the form of the random Schrödinger equation driven by a Brownian field studied in particular in [6]. The proof of the convergence of the solution of the wave equation in random media to the solution of the white-noise paraxial equation was obtained for stratified weakly fluctuating media in [1] and recently for three-dimensional random media in the context of acoustic waves in [14].

* Corresponding author. Tel.: +33 1 44 27 86 93; fax: +33 1 44 27 76 50. E-mail addresses: [email protected] (J. Garnier), [email protected] (K. Sølna). 0165-2125/$ - see front matter Ó 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.wavemoti.2009.02.002

238

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

A parabolic approximation for linear elastic waves is of great interest for solving migration problems in geophysics [8]. In the theory of elastic waves the field is described by the displacement vector. It is not straightforward to extend the parabolic approximation to elastic waves because two different wavenumbers coexist for plane waves in a given direction in a homogeneous medium, one for pressure waves and one for shear waves. The effective behavior of the coupling terms between shear waves and pressure waves in heterogeneous media is not easy to incorporate in the parabolic approximation. The analysis of the effective coupling is one of the main objectives of this work. In this paper we consider three-dimensional linear elastic waves. Motivated by problems in geophysics, we assume that the privileged propagation axis is the vertical direction. This choice explains the terminology ‘‘downward-going” waves and ‘‘upward-going” waves that we use in the context of the wave decomposition that we set forth. In the literature it is possible to find parabolic approximation models for the linear elastic wave equations which are obtained by suppressing the reflection terms. In [19] this idea is used to transform the second-order boundary value problem for two-dimensional elastic waves in laterally inhomogeneous media into a first-order initial value problem in space, whose numerical solutions can be computed efficiently. In [21] a parabolic approximation is constructed for three-dimensional elastic waves where the wave field disturbance is characterized by two almost plane waves travelling with wave speeds associated with pressure and shear waves. In [15] the author uses a multi-scale expansion and deals separately with the shear waves and the pressure waves. In [5] the reflection terms are suppressed based on the argument that the phases of the coupling terms between the downward-going and the upward-going waves vary rapidly making these terms negligible in the lowest order approximation. All these derivations are formal and give quantitatively different predictions [30]. Our goal is to provide a rigorous derivation of the paraxial wave equations for elastic waves in homogeneous media and, moreover, in random media. (1) In the case of a homogeneous medium, we consider the distinguished limit in which the propagation distance is much larger than the initial beam width, which is itself much larger than the typical wavelength. We apply an invariant imbedding theorem and an averaging theorem for rapidly oscillating differential equations in order to prove the convergence of the solution of the wave equation to the solution of a system of Schrödinger equations. The invariant imbedding allows us to transform the boundary value problem (with radiation conditions) into an initial value problem. The averaging theorem allows us to obtain effective paraxial wave equations for elastic waves. In particular, the coupling terms between downward-going waves and upward-going waves, and between shear waves and pressure waves, have rapid phases. By using an averaging theorem for highly oscillatory differential equations, we show that these coupling terms average out to give non-zero effective terms (in the form of Lie brackets) for downward-going shear waves and for the downward-going pressure waves. There is no effective energy transfer between the downward-going waves and the upward-going waves, neither between pressure waves and shear waves. However, the effective coupling between the wave modes affects the propagation of the downward-going shear waves and the one of the downward-going pressure waves. It is therefore important to take these effects into account. It turns out that the explicit computation of the effective Lie brackets gives very simple forms for the parabolic approximations of the pressure waves and the shear waves. (2) In the case of a random medium, we assume additionally that the correlation length of the random medium is of the same order as the initial beam width, and that the amplitude of the random fluctuations is small. By applying diffusion approximation theorems we prove the convergence of the solution of the random elastic wave equations to the solution of a system of Itô-Schrödinger equations. These Schrödinger equations are driven by two correlated Brownian fields, whose covariance function depends on the two-point statistics of the fluctuations of the density and Lamé parameters of the medium. This result shows that it is possible to justify both the parabolic approximation and the white-noise approximation in the distinguished limit considered in this paper. The limit system permits easy numerical simulations on the one hand, and exact computations of moments using Itô’s formula on the other hand. The paper is organized as follows. In Section 2 we describe the linear elastic wave equations in a random medium. In Section 3 we present the main results, i.e. the random paraxial wave equations. In Section 4 we introduce the fundamental wave decomposition that we use in the proofs of the main results. The proofs of the results are given in Section 5 when the medium is homogeneous and in Section 6 when the medium is randomly heterogeneous. 2. Elastic waves in a random medium We consider linear elastic waves propagating in a three-dimensional medium with heterogeneous and random fluctuations. The linearized elastic wave equations in an isotropic solid have the form of the hyperbolic system [7] 3 @vk X @ skl ¼ fk ; k ¼ 1; 2; 3;  @t @xl l¼1   X 3 1 @vj @vi @ skl  þ Sijkl ¼ hij ; 2 @xi @xj @t k;l¼1

q

ð1Þ i; j ¼ 1; 2; 3;

ð2Þ

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

239

where skl ðt; xÞ is the dynamic stress (symmetric), vk ðt; xÞ is particle velocity, fk ðt; xÞ is volume source density of force, hij ðt; xÞ is volume source density of deformation rate (symmetric), qðxÞ is volume density of mass, and Sijkl ðxÞ ¼ KðxÞdij dkl þ MðxÞðdik djl þ dil djk Þ is the compliance, with

KðxÞ ¼ 

kðxÞ ; ð3kðxÞ þ 2lðxÞÞ2lðxÞ

MðxÞ ¼

1 ; 4lðxÞ

and kðxÞ and lðxÞ are the Lamé coefficients of the medium. We also use the Kronecker symbol dij ¼ 1 if i ¼ j and 0 otherwise. Eqs. (1) and (2) correspond respectively to the equation of motion and the deformation equation for the small amplitude variations around the equilibrium distributions. If we introduce the pressure p and stress tensor component v that are given respectively by

" # 3 X @uj pðt; xÞ ¼ kðxÞ ðt; xÞ ; @xj j¼1



vkl ðt; xÞ ¼ lðxÞ

 @ul @uk ðt; xÞ þ ðt; xÞ ; @xk @xl

k; l ¼ 1; 2; 3;

ð3Þ

where uðt; xÞ ¼ ðuj ðt; xÞÞj¼1;2;3 is the displacement vector, then the 3  3 stress tensor s and the velocity vector v have the following form in an isotropic solid [7]:

skl ðt; xÞ ¼ pðt; xÞdkl þ vkl ðt; xÞ; k; l ¼ 1; 2; 3;

ð4Þ

@uk vk ðt; xÞ ¼ ðt; xÞ; @t

ð5Þ

k ¼ 1; 2; 3:

Regarding the stress (a symmetric tensor of rank two) we first remark that the vector sn can be interpreted as the traction vector exerted on an elementary part of the surface bounding a volume element with (outer) unit normal n. The integral of the traction vector over the surface then gives the net force acting on the volume element due to surface forces of contact type. In addition to surface forces there may be volume forces whose magnitude scales with the volume of the volume element, examples of which are gravity and electromagnetic forces. In most applications, for instance in geophysics, the influence of volume forces are negligible. The volume density of volume forces is the source term f in the equation of motion (1). The equation of motion follows from Newton’s law for an elementary volume element upon an application of Gauss integral P theorem to convert the traction surface integral to a volume integral with kernel l @ xl skl . Regarding the stress we second comment that the time-derivative of the stress tensor component v in (3) is a scaled version of the deformation rate (symmetric tensor of rank two) defined by ð@ xi vj þ @ xj vi Þ=2. This tensor specifies the rate of change of deformation that the solid undergoes during the passage of the elastic wave. For an infinitesimal oriented line segment given by the vector n the deformation rate matrix times n gives the rate of change of the vector that comes from material deformation. The external contribution to the deformation rate is h (symmetric tensor of rank two), and comes from the action of external sources that impress a deformation rate to the solid. The contribution to the deformation rate that is induced by the stress is defined via the fourth-rank tensor in (2), the compliance S, whose parametric dependence on K and M follows from the assumed constitutive relations. The form of the second term in the left-hand side of (2) follows from the linearization in typical scaling regimes, as in seismics. In this paper we focus our attention on propagation through a random section occupying the region x3 2 ðL; 0Þ with the sources f and h located outside of the random section at the surface x3 ¼ 0. The parameterization is motivated by waves probing for instance the heterogeneous earth and we may think of x3 as the main probing direction. We shall refer to waves propagating in a direction with a negative (resp. positive) x3 -coordinate as downward-going (respectively upward-going) waves. A main ingredient in our analysis will thus be to decompose the wave field in components that can be interpreted as being locally up- and down propagating modes. This analysis is facilitated by using the velocity, v, and stress tensor component, v, as the dependent variables. We can then write (1) and (2) in the form 3 @ vjk @vj @p X þ þ fj ðt; xÞ; j ¼ 1; 2; 3; ¼ @xj k¼1 @xk @t   @ vjk @vj @vk þ ¼ lðxÞ þ g jk ðt; xÞ; j; k ¼ 1; 2; 3; @t @xk @xj " # 3 X @p @vl ; ¼ kðxÞ @t @xl l¼1

qðxÞ

ð6Þ ð7Þ ð8Þ

where the 3  3 source term g is given by

g jk ðt; xÞ ¼ 2lðxÞhjk ðt; xÞ  kðxÞ

" 3 X

# hll ðt; xÞ djk ;

j; k ¼ 1; 2; 3:

l¼1

If the shear modulus l ¼ 0 and the deformation source h = 0, then v ¼ 0 and we get the acoustic wave equations for the pressure p and the velocity v. The acoustic equations in the context of random media are studied in detail in for instance [11,13,14] and we shall here generalize the analysis to the above system for elastic waves. Elastic wave propagation in layered media is studied in for instance [16,17,28], while elastic waves in the transport regime is studied in [24].

240

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

We consider in this paper the situation in which the random medium constitutes a random section occupying the region x3 2 ðL; 0Þ that is sandwiched in between two homogeneous half-spaces. The medium parameters are:

(

l1 if x3 6 L or x3 P 0; 0   x  3 if x3 2 ðL; 0Þ; l1 1 þ e m l 2 0 e ( q0 if x3 6 L or x3 P 0;    qðxÞ ¼ q0 1 þ e3 mq ex2 if x3 2 ðL; 0Þ; 1 ¼ lðxÞ

(

kðxÞ ¼

k0 if x3 6 L or x3 P 0;    if x3 2 ðL; 0Þ; k0 1 þ e3 mk ex2

ð9Þ ð10Þ ð11Þ

with e a small parameter. The parameter e2 characterizes the ratio of the correlation length of the random medium fluctuations to the thickness of the random section. The random processes ml ðxÞ; mq ðxÞ, and mk ðxÞ describe the medium fluctuations. We assume that they are stationary and zero-mean and that they satisfy strong mixing conditions in x3 . We also assume that the amplitude of the random fluctuations is of order e3 . This is the interesting regime in which the medium fluctuations give rise to effective terms of order one when e goes to zero, as we will show in this paper. Recall also that we let the source be located at the surface in the plane x3 ¼ 0. We consider a scaling regime in which the spatial support of the source (in the transverse direction), or equivalently the initial beam width, is of order e2 . This means that the transverse scale of the source and the one of the spatial fluctuations of the medium are of the same order. Remember that the Rayleigh length for a beam with initial beam width r0 and carrier wavelength k0 is of the order of r 20 =k0 in absence of random fluctuations (the Rayleigh length is the distance from beam waist where the beam area is doubled by diffraction). In order to get a Rayleigh length of order one, we assume that the carrier wavelength of the source is of order e4 . Therefore, in this regime the source terms have the form

f ðt; xÞ ¼ f

ð0Þ



 t x? ; 2 dðx3 Þ; 4

e e

hðt; xÞ ¼ h

ð0Þ



 t x? ; 2 dðx3 Þ; 4

e e

ð12Þ

where x ¼ ðx? ; x3 Þ; x? ¼ ðx1 ; x2 Þ. These source terms model volume forces and deformation rates acting on a thin region in the x3 ¼ 0 plane. They will give jump conditions on stress and velocities across the interface. The source then generates an elastic wave beam, with shear and compressional components, that propagate mainly along the x3 -axis as we explain in detail in Section 3. We remark than an explosive source occasionally is treated as a point source. Here we model in terms of a ‘‘disk” source which may reflect for instance a seismic exploration source configuration where a set of explosive sources are treated collectively to form a probing beam or the situation with propagation from an explosion in a thin cavity [29]. The sources generate a power flow into the medium and we finish this section by an expression for this quantity which we will use in Section 3. The elastodynamic power flow across an interface @D with unit normal n is

P @D ðtÞ ¼

Z

nðxÞ  pðt; xÞdrðxÞ;

@D

where p ¼ sv is the elastodynamic Poynting vector [7]. Therefore the power flow across the horizontal plane at depth x3 in the negative x3 -direction (which will correspond to the beam direction in our set up) is

Pðt; x3 Þ ¼

Z X 3

vl sl3 ðt; x? ; x3 Þdx? ;

l¼1

and the total energy flux in this direction is

F ðx3 Þ ¼

Z Z X 3

vl sl3 ðt; x? ; x3 Þdx? dt:

ð13Þ

l¼1

3. The paraxial wave equations In this section we state the main results of the paper. Proposition 1 describes the paraxial equations for elastic waves when the medium is homogeneous and it is proved in Section 5. Proposition 2 describes the white-noise paraxial equations for elastic waves when the medium is randomly heterogeneous and it is proved in Section 6. The proofs are based on the wave mode decomposition set forth in Section 4. Proposition 1. Let us assume that the medium is homogeneous ml ¼ mq ¼ mk ¼ 0. The transmitted wave consists of the succession of two fields that emerge at x3 ¼ L around time t ¼ L=cP and t ¼ L=cS , respectively:

  L e!0 ðv; v; pÞ t ¼ þ e4 s; e2 x? ; x3 ¼ L ! ðvP ; vP ; pP Þðs; x? Þ; cP   L e!0 ðv; v; pÞ t ¼ þ e4 s; e2 x? ; x3 ¼ L ! ðvS ; vS ; pS Þðs; x? Þ: cS

ð14Þ ð15Þ

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

241

The first field ðvP ; vP ; pP Þ is a pressure wave of the form

vPjk ðs; x? Þ ¼ 0; ðj; kÞ–ð3; 3Þ;

ð16Þ

Z 1 vP33 ðs; x? Þ ¼ v P33 ðx; x? Þ expðixsÞdx; 2p vPj ðs; x? Þ ¼ 0; j ¼ 1; 2;

ð17Þ ð18Þ

f1 P

vP33 ðs; x? Þ; k0 pP ðs; x? Þ ¼ v ðs; x? Þ: 2l0 P33

vP3 ðs; x? Þ ¼

ð19Þ ð20Þ

The second field ðvS ; vS ; pS Þ is a shear wave of the form

Z 1 v Sj3 ðx; x? Þ expðixsÞdx; j ¼ 1; 2; 2p vSjk ðs; x? Þ ¼ 0; ðj; kÞ R fð1; 3Þ; ð2; 3Þ; ð3; 1Þð3; 2Þg;

vSj3 ðs; x? Þ ¼ vSj ðs; x? Þ ¼

f1 S

vSj3 ðs; x? Þ; j ¼ 1; 2;

ð21Þ ð22Þ ð23Þ

vS3 ðs; x? Þ ¼ 0;

ð24Þ

pS ðs; x? Þ ¼ 0:

ð25Þ

Here cS and cP are the shear-wave and pressure-wave velocities, and fS and fP are the shear-wave and pressure-wave impedances of the homogeneous medium:

rffiffiffiffiffiffi l0 ; cS ¼

cP ¼

q0

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2l0 þ k0

q0

;

pffiffiffiffiffiffiffiffiffiffiffi f S ¼ l0 q 0 ;

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 4l20 q0 fP ¼ : 2l0 þ k0

ð26Þ

 P33 and v  Sj3 ; j ¼ 1; 2, are given by The transmitted time-harmonic fields v

v Sj3 ðx; x? Þ ¼ v P33 ðx; x? Þ ¼

Z

Z

 inc;j ðx; x0? Þdx0? ; T S ðx; x?  x0? ; 0Þv

j ¼ 1; 2;

 inc;3 ðx; x0? Þdx0? ; T P ðx; x?  x0? ; 0Þv

where the incoming wave fields generated by the source at x3 ¼ 0 are

fS ð0Þ 1 ð0Þ g ðx; x? Þ þ f j ðx; x? Þ; j ¼ 1; 2; 2 2l0 3j fP l0 ð0Þ ð0Þ g ðx; x? Þ þ v inc;3 ðx; x? Þ ¼ f ðx; x? Þ; 2ð2l0 þ k0 Þ 33 2l0 þ k0 3 Z f ð0Þ ðx; x? Þ ¼ f ð0Þ ðs; x? Þ expðixsÞds; Z ð0Þ ðx; x? Þ ¼ gð0Þ ðs; x? Þ expðixsÞds; g " # 3 X ð0Þ ð0Þ ð0Þ g jk ðs; x? Þ ¼ 2l0 hjk ðs; x? Þ  k0 hll ðs; x? Þ djk ; j; k ¼ 1; 2; 3:

v inc;j ðx; x? Þ ¼

ð27Þ ð28Þ ð29Þ ð30Þ ð31Þ

l¼1

The transmission operators T P ðx; x? ; x3 Þ and T S ðx; x? ; x3 Þ satisfy the Schrödinger equations for x3 2 ðL; 0Þ:

@ T P ðx; x? ; x3 Þ i ¼ Dx T P ðx; x? ; x3 Þ; T P ðx; x? ; x3 ¼ LÞ ¼ dðx? Þ; @x3 2kP ? @ T S ðx; x? ; x3 Þ i ¼ Dx T S ðx; x? ; x3 Þ; T S ðx; x? ; x3 ¼ LÞ ¼ dðx? Þ: @x3 2kS ? 2

ð32Þ ð33Þ

2

@ @ Here Dx? ¼ @x 2 þ @x2 is the transverse Laplacian and kS ¼ x=c S and kP ¼ x=c P are shear-wave and pressure-wave wavenumbers. 1

2

Eqs. (27) and (28) describe the downward-going waves generated by the sources. In particular, both the volume source density of force f and the volume source density of deformation rate h generate pressure and shear waves. Eqs. (32) and (33) describe how the downward-going pressure and shear waves propagate into the homogeneous medium. The form of this Schrödinger system is simple but its derivation is not obvious. One has to show that the forward-scattering approximation is valid. One has to take into account the coupling between shear and pressure wave modes and the coupling between upward- and downward-going modes. As we will see, the brutal suppression of the upward-going (i.e. reflected) modes from the equations leads to wrong results. The reflected modes are vanishing in the limit e ! 0, which means that there is no transfer of energy from the downward-going modes to the upward-going modes, but the averaging of the

242

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

coupling terms between upward- and downward-going modes gives effective terms that contribute to the expressions of the diffraction operators. The overall result is that the diffraction operators for the pressure and shear waves take the form of transverse Laplacians. We now give the result when the medium has random fluctuations. Proposition 2. In the limit e ! 0, the transmitted wave consists of the succession of two fields that emerge at x3 ¼ L around time t ¼ L=cP and t ¼ L=cS as described by (14) and (15). The first field is a pressure wave of the form (16)–(20). The second field  Sj3 ; j ¼ 1; 2 are given by  P33 and v comprises two shear waves of the form (21)–(25). The time-harmonic fields v

v P33 ðx; x? Þ ¼ v Sj3 ðx; x? Þ ¼

Z

Z

 inc;3 ðx; x0? Þdx0? ; T P ðx; x? ; x0? ; 0Þv  inc;j ðx; x0? Þdx0? ; T S ðx; x? ; x0? ; 0Þv

j ¼ 1; 2:

The operators T P and T S are the solutions of the following Itô-Schrödinger diffusion models for x3 2 ðL; 0Þ:

i ikP  Dx0 T P ðx; x? ; x0? ; x3 Þdx3 þ T P ðx; x? ; x0? ; x3 Þ  dBP ðx0? ; x3 Þ; 2kP ? 2 i ikS  T S ðx; x? ; x0? ; x3 Þ  dBS ðx0? ; x3 Þ; Dx0 T S ðx; x? ; x0? ; x3 Þdx3 þ dT S ðx; x? ; x0? ; x3 Þ ¼ 2kS ? 2

dT P ðx; x? ; x0? ; x3 Þ ¼

ð34Þ ð35Þ

with the initial conditions

T P ðx; x? ; x0? ; x3 ¼ LÞ ¼ T S ðx; x? ; x0? ; x3 ¼ LÞ ¼ dðx?  x0? Þ:

ð36Þ

Here the symbol  stands for the Stratonovich stochastic integral, BP ðx? ; x3 Þ and BS ðx? ; x3 Þ are two correlated Brownian fields with the covariance matrix



E½BP ðx? ; x3 ÞBP ð; x0? ; x03 Þ E½BP ðx? ; x3 ÞBS ð; x0? ; x03 Þ

E½BS ðx? ; x3 ÞBP ð; x0? ; x03 Þ E½BS ðx? ; x3 ÞBS ð; x0? ; x03 Þ   C 0;PP ðx? Þ C 0;PS ðx? Þ C0 ðx? Þ ¼ ; C 0;PS ðx? Þ C 0;SS ðx? Þ



¼ minfx3 ; x03 gC0 ðx?  x0? Þ;

ð37Þ ð38Þ

where we have used the notations

C PP ðxÞ ¼ E½mP ðx0 þ xÞmP ðx0 Þ; Z 1 C 0;PP ðx? Þ ¼ C PP ðx? ; x3 Þdx3 ;

ð39Þ ð40Þ

1

C SS ðxÞ ¼ E½mS ðx0 þ xÞmS ðx0 Þ; Z 1 C 0;SS ðx? Þ ¼ C SS ðx? ; x3 Þdx3 ;

ð41Þ ð42Þ

1

C PS ðxÞ ¼

1 1 E½mP ðx0 þ xÞmS ðx0 Þ þ E½mS ðx0 þ xÞmP ðx0 Þ; 2 Z 2

ð43Þ

1

C PS ðx? ; x3 Þdx3 ;

C 0;PS ðx? Þ ¼

ð44Þ

1

with

mS ðxÞ ¼ ðml þ mq ÞðxÞ;

ð45Þ

k0 mP ðxÞ ¼ ðml þ mq ÞðxÞ  ðml þ mk ÞðxÞ: 2l0 þ k0

ð46Þ

The Itô-Schrödinger equations (34) and (35) describe how the pressure and shear waves propagate into the random medium. The important qualitative results are that the forward scattering approximation is valid in the scaled regime considered in the paper, and that the random fluctuations of the medium can be modeled by a pair of correlated Brownian fields that only depends on the two-point statistics of the random fluctuations of the parameters of the medium. In order to give a complete characterization of the transmitted field, we add that there is no other transmitted wave in the sense that for any time t0 R fL=cP ; L=cS g, e!0

ðv; v; pÞðt ¼ t 0 þ e4 s; e2 x? ; x3 ¼ LÞ ! ð0; 0; 0Þ: There is no reflected wave in the sense that for any time t0 > 0,

243

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

e!0

ðv; v; pÞðt ¼ t 0 þ e4 s; e2 x? ; x3 ¼ 0Þ ! ð0; 0; 0Þ: In [6] the existence and uniqueness has been established for each of the two random processes

V P ðx; x? ; x3 Þ ¼

Z

T P ðx; x? ; x0? ; x3 Þ/P ðx0? Þdx0? ;

V S ðx; x? ; x3 Þ ¼

Z

T S ðx; x? ; x0? ; x3 Þ/S ðx0? Þdx0? ;

for any test functions /S and /P with unit L2 ðR2 ; CÞ-norm. Furthermore, it is shown that each of the two processes V S ðx; x? ; x3 Þ and V P ðx; x? ; x3 Þ is a continuous Markov diffusion process on the unit ball of L2 ðR2 ; CÞ. The moment equations moreover satisfy a closed system at each order [12]. It is straightforward to generalize these results for the joint distribution of the pair of processes ðV S ; V P Þ. The interested reader can find explicit calculations for the mean intensity, the autocorrelation function, and other relevant quantities of the transmitted wave in [14]. The fact that the solutions of Itô-Schrödinger equations have constant L2 -norms shows that the sum of the energies of the transmitted pressure wave (14) and of the transmitted shear waves (15) (in homogeneous or random media) is equal to the energy of the incoming wave. Using (13) we indeed find that the total energy flux generated by the source and entering the random section at x3 ¼ 0 is (in the limit of e small)

F inc ¼

e4 2p

Z Z

 inc;1 ðx; x? Þj2 jv  inc;2 ðx; x? Þj2 a0 jv  inc;3 ðx; x? Þj2 jv þ þ dxdx? ; fS fS fP

for a0 ¼ ðk0 þ 2l0 Þ=ð2l0 Þ. Moreover, the corresponding total energy flux of the transmitted pressure wave (14) exiting at x3 ¼ L is

F tr;P ¼

e4 2p

Z Z

Z

a0  fP 

2   inc;3 ðx; x0? Þdx0?  dxdx? ¼ T P ðx; x? ; x0? ; 0Þv

e4 2p

Z Z

a0 jv inc;3 ðx; x? Þj2 fP

dxdx? ;

and the total energy flux of the transmitted shear wave (15) exiting at x3 ¼ L is

F tr;S ¼

2 e4 X 2p j¼1

Z Z

2 Z 2 Z Z   inc;j ðx; x? Þj2 1   e4 X jv 0 0 0   ð x ; x ; x ; 0Þ v ð x ; x Þdx d x dx ¼ dxdx? : T S ? ? inc;j ? ? ?  fS 2p j¼1 fS

This shows that we have F inc ¼ F tr;P þ F tr;S , which illustrates the fact that no energy is reflected and all the energy is transmitted in the form of the two wave fields described in the propositions.

4. The mode decomposition The purpose of the following sections is to prove Propositions 1 and 2. We rescale the transverse spatial variables so as to observe the wave at the scale of the source and we take a scaled Fourier transform in time:

ej ðx; xÞ ¼ v

Z

vj ðe4 t; e2 x? ; x3 Þ expðixtÞdt;

v ejk ðx; xÞ ¼

Z

vjk ðe4 t; e2 x? ; x3 Þ expðixtÞdt:

 e13 ; v  e23 ; v  e33 Þ satisfies e1 ; v e2 ; v e3 ; v Except at x3 ¼ 0 (where the source is) the six-dimensional vector ðv

ej @v

e3 ix e 1 @v   ¼ 4 v ; j ¼ 1; 2; @x3 e l j3 e2 @xj e3 @v ix e  ; ¼ 4 v @x3 2e l 33            e13 e1 e2 @v iqx e @ k þ 2l @ v @ l @ ve1 @ k @v @ l @ ve2 1 @ k e 1 þ ¼ 4 v v 33 ; þ þ þ  2 @x1 @x2 ix @x2 @x1 ix @x2 @x2 ix @x1 @x3 e ix @x1 e @x1 2l            e23 e2 e2 e1 e1 @v iqx e @ k þ 2l @ v @ l @v @ k @v @ l @v 1 @ k e 2 þ þ þ þ  2 ¼ 4 v v 33 ; @x2 @x1 ix @x1 @x2 ix @x1 @x1 ix @x2 @x3 e ix @x2 e @x2 2l " !  e    e  2 e 2 e e e    1 @ v e2   @ v33 2l iqx e k @ v3 @ v3 1 k þ l @ v13 @ v23 k @v 2 @ 3   ¼  4 v þ þ e þ þ ix @x21 @x3 ix @x1 @x2 @x3 k þ 2l e e2 l @x1 @x2 @x22        @ k k @ 1 k @ 1  v e  v e  v e ; @x3 2l 33 e2 @x1 l 13 e2 @x2 l 23 where

ð47Þ ð48Þ ð49Þ ð50Þ

ð51Þ

l; q and k are evaluated at ðe2 x? ; x3 Þ. The four other dependent variables v e11 ; v e12 , v e22 , and pe are given by:

v ejl ¼ 

e2 l @ vej e2 l @ vel  ; j; l ¼ 1; 2; ix @xl ix @xj

¼ p

e2 k @ ve1 e2 k @ ve2 k  e  þ v : ix @x1 ix @x2 2l 33

ð52Þ

244

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

When the medium is homogeneous ðml ¼ mq ¼ mk ¼ 0Þ, then there exist plane wave solutions that depend only on x3 and satisfy:

0 e 1 0 e 1 v1 v1 e C e C Bv Bv B  13 C B  13 C B e C B e C  C  C Bv v d B B 2 C ¼  i x MB 2 C ; Bv e C e C 4   dx3 B e v B 23 C B 23 C B e C B e C 3 A 3 A @ v @v v e33 v e33

0

0 Bq B 0 B B 0 B M¼B B 0 B B 0 @ 0

l1 0

0

0

0

0 0

0 0

0

l1 0

0 0

0

q0

0

0

0

0

0

0

0

2l0 q0 k0 þ2l0

0

0

0

1

C C C C C C: 0 C C 1 C ð2l0 Þ A 0 0 0

The diagonalization of the 6  6 matrix M gives the general form of the plane wave solutions

     kx  ðxÞ exp i kS x3 j ðxÞ exp i S 4 3 þ b ej ðx; x3 Þ ¼ f1=2 ; j ¼ 1; 2; a v j S 4 e e      kx  ðxÞ exp i kS x3 j ðxÞ exp i S 4 3 þ b v ej3 ðx; x3 Þ ¼ f1=2 a ; j ¼ 1; 2; j S e e4      k x 3 ðxÞ exp i kP x3 3 ðxÞ exp i P 4 3 þ b e3 ðx; x3 Þ ¼ f1=2 a ; v P e e4      k x 3 ðxÞ exp i kP x3 3 ðxÞ exp i P 4 3 þ b v e33 ðx; x3 Þ ¼ f1=2 a : P 4

e

ð53Þ ð54Þ ð55Þ ð56Þ

e

 ; j ¼ 1; 2 j ; j ¼ 1; 2 are upward-going shear waves, the mode a 3 is an upward-going pressure wave, the modes b The modes a j 3 is a downward-going pressure wave. are downward-going shear waves, and the mode b When the medium is heterogeneous and described by the model (9)–(11), then we introduce the generalized upwarde defined by ej and downward-going modes b going modes a j

     kx e ðx; xÞ exp i kS x3 ej ðx; xÞ exp i S 4 3 þ b ej ðx; xÞ ¼ f1=2 ; j ¼ 1; 2; a v j S e e4      kx e ðx; xÞ exp i kS x3 ; j ¼ 1; 2; ej ðx; xÞ exp i S 4 3 þ b v ej3 ðx; xÞ ¼ f1=2 a j S e e4      k x e ðx; xÞ exp i kP x3 ; e3 ðx; xÞ exp i P 4 3 þ b e3 ðx; xÞ ¼ fP1=2 a v 3 4 e e      k x k x P 3 1=2 v e33 ðx; xÞ ¼ fP ae3 ðx; xÞ exp i 4 þ be3 ðx; xÞ exp i P 4 3 ;

e

ð57Þ ð58Þ ð59Þ ð60Þ

e

where cS ; cP ; fS , and fP are the shear-wave and pressure-wave velocities and impedances of the background homogeneous medium defined by (26). The radiation conditions at þ1 and 1 impose that there is no downward-going wave in the homogeneous region x3 > 0 and no upward-going wave in the homogeneous region x3 < L:

ej ðx; x? ; x3 ¼ ðLÞ Þ ¼ 0; a

e ðx; x? ; x3 ¼ 0þ Þ ¼ 0; b j

j ¼ 1; 2; 3;

ð61Þ

see Fig. 1. The jump conditions across the source plane x3 ¼ 0 are obtained by integrating (6)–(8) across x3 ¼ 0: þ

ej 00 ¼  ½v þ

e3 00 ¼  ½v

1

l

ð0Þ gj3 ðx; x? Þ;

ð0Þ @ g33 k0 ; x 2l0 þ k0 @xj

þ 2 i  ej3 00 ¼ f ð0Þ ½v j þe

1 ð0Þ g ðx; x? Þ; 2l0 þ k0 33

þ

 e33 00 ¼  ½v

j ¼ 1; 2;

ð62Þ !

ð0Þ ð0Þ 2l0 ð0Þ i 2k0 @ g31 @ g32 : f 3 þ e2 þ 2l0 þ k0 x 2l0 þ k0 @x1 @x2

Fig. 1. Boundary conditions for the wave modes with a source at x3 ¼ 0 and radiation conditions at 1.

ð63Þ

245

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

 33 that is localized in In order to be complete, we mention the existence of a singular component to the stress tensor entry v ð0Þ ð x ; x Þdðx Þ, moreover, singular components to the stress tensor entries the plane x3 ¼ 0 and that is of the form e4 xð2lik0þk0 Þ g ? 3 33 0 v jk ; j; k ¼ 1; 2 that are also localized in the plane x3 ¼ 0 and that are of the form e4 igð0Þ ðx; x? Þdðx3 Þ=x. jk From the jump conditions (62) and (63) and the radiation conditions at þ1 we obtain

 ðx; x? Þ; e ðx; x? ; x3 ¼ 0 Þ ¼ b b inc;j j

j ¼ 1; 2; 3;

ð64Þ

with

! ð0Þ 1=2 2 @ g33 k0  ðx; x Þ ¼ fS gð0Þ ðx; x Þ þ 1 f ð0Þ  ie ðx; x? Þ; b inc;j ? ? j 1=2 2l0 j3 x 2l0 þ k0 @xj 2fS  b inc;3 ðx; x? Þ ¼

f1=2 1 ð0Þ P g ðx; x? Þ þ 1=2 2ð2l0 þ k0 Þ 33 fP ð2l0 þ k0 Þ

l0f ð0Þ 3 

ie2

x

k0

j ¼ 1; 2;

ð0Þ ð0Þ @ g31 @ g32 þ @x1 @x2

ð65Þ !! ðx; x? Þ;

ð66Þ

which also gives (27) and (28). By the continuity of the fields across z ¼ L and the radiation condition at 1, we get:

ej ðx; x? ; x3 ¼ ðLÞþ Þ ¼ 0; a

j ¼ 1; 2; 3:

ð67Þ

Let us next introduce the six-dimensional vector

1 e1 a B e C B b1 C B C Ba e2 C C  e ðx; x? ; x3 Þ ¼ B X B e Cðx; x? ; x3 Þ: B b2 C B C C Ba @ e3 A e b 0

ð68Þ

3

 e satisfies in the section x3 2 ðL; 0Þ the linear system From (47)–(51), the vector X

e dX  e ðx; x? ; x3 ÞX  e; ¼A dx3

ð69Þ

 e ðx; x? ; x3 Þ A x







i 1h 3  1;1 þ exp i x3 ðkP  kS Þ A  1;1 þ exp i x3 ðkP þ kS Þ A  1;1 þ exp i x3 ðkP  kS Þ A  1;1 ¼ 2 exp i 4 ðkP þ kS Þ A 4 4 4 e e   e  e       e    1 x3 2kP x3  2;0 2kP x3  2;0 1 x3 2kS x3  0;2 2kS x3  0;2 ~ ~ B þ exp i 4 B B þ exp i 4 B þ mS x? ; 2 exp i 4 þ mP x? ; 2 exp i 4 e e e e e e e e       1 x3  0;0 1 x3  0;0  0;0 kP x3  2;0 kP x3  2;0 kS x3  0;2 þ exp 2i 4 A þ mP x? ; 2 BP þ mS x? ; 2 BS þ A þ exp 2i 4 A þ exp 2i 4 A e e e e e e e   kS x3  0;2 þ exp 2i 4 A ; ð70Þ

e

and the two-point boundary conditions

 e ðx; x? ; LÞ þ H0 X  e ðx; x? ; 0Þ ¼ Vð  x; x? Þ: HL X

ð71Þ

 j;k and B  j;k are given in Appendix A. The vector Vð  x; x? Þ and the matrices HL and H0 are defined by The matrices A

1 0  C Bb B inc;1 ðx; x? Þ C C B C B 0  x; x? Þ ¼ B C; Vð C Bb  B inc;2 ðx; x? Þ C C B A @ 0  binc;3 ðx; x? Þ 0

0

HL

1 B0 B B B0 ¼B B0 B B @0

0 0 0 0 0 0 0 1 0 0

0 0 0 0 0

0

1

0 0 0C C C 0 0 0C C; 0 0 0C C C 0 1 0A 0 0 0

0

0 B0 B B B0 0 H ¼B B0 B B @0

0 0 0

0 0

1

0 0C C C 0 0 0 0 0C C; 0 0 1 0 0C C C 0 0 0 0 0A 0 0 0 0 0 1 1 0 0

ð72Þ

246

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

the processes mS ðxÞ and mP ðxÞ are defined by (45) and (46), and

~ S ðxÞ ¼ ðml  mq ÞðxÞ; m

ð73Þ

k0 ~ P ðxÞ ¼ ðml  mq ÞðxÞ þ m ðml þ mk ÞðxÞ: 2l0 þ k0

ð74Þ

 e we have neglected terms of order e and smaller, and kept only the terms of order e2 ; e1 and 1. In the expression (70) of A The mode decomposition (i.e. the definitions of the mode velocities and impedances) has canceled the terms of order e4 that are present in the original Eqs. (47)–(51). If the source is x? -independent and if the medium is homogeneous ðml ¼ mq ¼ mk ¼ 0Þ, then the solution of (69)–(71) is a constant vector which corresponds to a collection of three downward-going modes: two shear-wave modes that propagate with the velocity cS and one pressure-wave mode that propagates with the velocity cP . If the source is x? -dependent and/or the medium is heterogeneous, then transverse spatial effects and/or random effects have to be taken into account:  j;k in (70) correspond to deterministic transverse spatial effects. For ðj; kÞ–ð0; 0Þ, these terms have large (1) The matrices A amplitudes, of order e2 , but they also have rapid phases that vary at the scale e4 , and we will see that they give rise to effective terms of order one in the limit e ! 0 through the application of an averaging theorem for highly oscillatory differential equations. These effects are purely deterministic, they consist of coupling between modes with different velocities, including coupling between upward- and downward-going modes, and it is crucial to take them into  0;0 that would be the natural candiaccount in the determination of the paraxial wave equation. Indeed, the matrix A date for the right-hand side of the paraxial wave equation does not capture the full paraxial mechanism. In other words, ‘‘brutally” suppressing the upward-going (i.e. reflected) modes leads to wrong results. As we will see, the reflected modes are vanishing in the limit e ! 0, but the coupling terms between upward- and downward-going modes give important contributions to the dynamics of the transmitted waves.  j;k in (70) correspond to downward and upward scattering and coupling between modes due to (2) The matrices B the random heterogeneities of the medium. These terms have large amplitudes, of order e1 , but they vary rapidly at the scale e2 , and the driving processes have mean zero and mixing properties. They will also give rise to effective terms of order one in the limit e ! 0 through the application of a diffusion approximation theorem.  e will be characterized by a stochastic partial differential equation driven by two correlated BrownThe limit of X ian fields.

5. The derivation of the paraxial wave equations in the homogeneous case In this section we assume that there are no random fluctuations ml ¼ mq ¼ mk ¼ 0 and we prove Proposition 1. We transform the two-point boundary value problem (69)–(71) into an initial value problem. This is done by an invariant imbedding step in which we introduce transmission and reflection matrices. The principle of the invariant imbedding approach is explained in detail in [Chapter 4][11] and was used in the context of elastic plane waves in layered media in [Chapter 5] [25]. First, we define the lateral Fourier modes for j ¼ 1; 2; 3:

^ej ðx; j? ; x3 Þ ¼ a ^ e ð x; j ; x Þ ¼ b ? 3 j

Z Z

ej ðx; x? ; x3 Þ expðij?  x? Þdx? ; a e ðx; x ; x Þ expðij  x Þdx ; b ? 3 ? ? ? j

^ej ðx; j? ; x3 ¼ 0Þ and the transmitted wave by where we have denoted j? ¼ ðj1 ; j2 Þ. The reflected wave is characterized by a ^e ðx; j? ; x3 ¼ LÞ. The parameters x and j? are frozen in the problem, so we shall not write explicitly the ðx; j? Þ-depenb j b e ðx3 Þ is defined as in (68), however, in terms of the Fourier modes a ^ej and dence of the vectors and matrices. The vector X e ^ bj and it is the solution of the two-point boundary value problem:

be dX b e ðx3 Þ X b e; ¼A dx3

b; b e ðLÞ þ H0 X b e ð0Þ ¼ V HL X

b e ðx3 Þ given by with the 6  6 matrix A

h i







b e ðx3 Þ ¼ 1 exp i x3 ðkP þ kS Þ A b 1;1 þ exp i x3 ðkP  kS Þ A b 1;1 þ exp i x3 ðkP þ kS Þ A b 1;1 þ exp i x3 ðkP  kS Þ A b 1;1 A e2 e4  e4 e4  e4       k x k x k x k x P 3 P 3 S 3 S 3 b 2;0 þ exp 2i b 2;0 þ exp 2i b 0;2 þ exp 2i b 0;2 ; b 0;0 þ exp 2i A A A A þA 4 4 4 4

e

e

e

e

ð75Þ

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

247

b j;k are given in Appendix A, HL and H0 are given by (72), and A

0

1

0

^ C Bb B inc;1 C C B C B b ¼ B 0 C; V C Bb ^ B inc;2 C C B @ 0 A ^ b

^ ð x; j Þ ¼ b ? inc;j

Z

 ðx; x? Þ expðij?  x? Þdx? ; b inc;j

j ¼ 1; 2; 3:

inc;3

By applying Proposition B.1 of Appendix B, we get that the reflected and transmitted modes are given by

b e ðx3 ¼ 0Þ V b ; ^ej ðx3 ¼ 0Þ ¼ ½ R a 2j1 ^e ðx3 ¼ LÞ ¼ ½ T b e ðx3 ¼ 0Þ V b ; b j 2j

j ¼ 1; 2; 3;

ð76Þ

j ¼ 1; 2; 3;

ð77Þ

b e and T b e are solution of the initial value problem where the reflection and transmission matrices R

be dR b e ðx3 Þ R b e H0 Þ A b e; R b e ðx3 ¼ LÞ ¼ I; ¼ ðI  R dx3 be dT b e ðx3 Þ R b e H0 A b e; T b e ðx3 ¼ LÞ ¼ I; ¼ T dx3

ð78Þ ð79Þ

and I is the 6  6 identity matrix. Note that the linear boundary value problem has been transformed into a nonlinear initial value problem, that has the form of a matrix Riccati equation. We can now apply Proposition C.1 which gives the asymptotic behavior of the reflection and transmission matrices as e ! 0:

h i h i b dR i i b 0;0 R b 1;1 R; b 1;1 R b 1;1 R; b 1;1 R bþ b ðI  RH b 0ÞA b þ b ðI  RH b 0ÞA b ; b 0ÞA b 0ÞA b 0ÞA ¼ ðI  RH ðI  RH ðI  RH dx3 kP þ kS kP  kS h i h i b dT i i b 0;0 R b 1;1 R; b 1;1 R b 1;1 R; b 1;1 R b 0A bþ b 0A b TH b 0A b þ b 0A b TH b 0A b ; ¼  TH TH TH dx3 kP þ kS kP  kS where the brackets stand for Lie brackets as defined in Proposition C.1. It turns out that these equations can be dramatically simplified. Indeed, the Lie brackets can be computed

b 1;1  A b 1;1 Þ R; b 1;1 R; b 1;1 R b 1;1 A b 1;1 A b ðI  RH b 0ÞA b ¼ ðI  RH b 0 Þð A b 0ÞA b ½ðI  RH b 1;1  A b 1;1 Þ R; b 1;1 R; b 1;1 R b 1;1 A b 1;1 A b ðI  RH b 0ÞA b ¼ ðI  RH b 0 Þð A b b 0ÞA ½ðI  RH b 1;1 R; b 1;1 R b 1;1  A b 1;1 ÞH0 R; b 1;1 A b 1;1 A b TH b 0A b ¼ Tð b A b b 0A ½ TH 0 0 0 1;1 1;1 1;1 1;1 1;1 1;1 b b b b b b b A b TH b A b ¼ Tð b A b ½ TH R; R A A A ÞH R; and we have the remarkable identity

b 0;0 þ A

i i b 1;1  A b 1;1 Þ þ b 1;1  A b 1;1 Þ ¼ D; b 1;1 A b 1;1 A b 1;1 A b 1;1 A b ðA ðA kP þ kS kP  kS

with

0 B B 2B ijj? j B b B D¼ 2 B B B @

1

1=kS

0

0

0

0

0

0

1=kS

0

0

0

0

0

1=kS

0

0

0 0

0 0

0 0

1=kS 0

0 1=kP

0 C C C 0 C C: 0 C C C 0 A

0

0

0

0

0

1=kP

b 0;0 and the Lie brackets coming from the averaging theorem combine It is absolutely remarkable that the zero-th order term A to give a diagonal matrix, while each contribution is rather complicated. This gives the following equations for the reflection and transmission matrices:

248

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

b dR b 0ÞD b R; b b 3 ¼ LÞ ¼ I; ¼ ðI  RH Rðx dx3 b dT b R; b b 3 ¼ LÞ ¼ I: b 0D Tðx ¼  TH dx3

ð80Þ ð81Þ

The solution of (80) is a diagonal matrix whose even diagonal terms are equal to one:

0 B B B B b 3Þ ¼ B Rðx B B B B @

b S ðx3 Þ 0 R

0

0

0

0

0

1

0

0

0 0

0 0

0 0

1 0

0

0

0

0

b S ðx3 Þ 0 R

0

0

1

C 0C C C 0 0C C; 0 0C C b P ðx3 Þ 0 C R A 0 1 0

where

bS dR ijj? j2 b R S ðx3 Þ; ðx3 Þ ¼  dx3 2kS bP dR ijj? j2 b R P ðx3 Þ; ðx3 Þ ¼  dx3 2kP

b S ðx3 ¼ LÞ ¼ 1; R b P ðx3 ¼ LÞ ¼ 1: R

Therefore, for j ¼ 1; 2; 3, we have

b 3 ¼ 0Þ V b a ej ðx3 ¼ 0Þ ¼ ½ Rðx lim b 2j1 ¼ 0: e!0

^ej ðx3 ¼ 0Þ This shows that the paraxial (or forward-scattering) approximation is valid in the sense that the reflected modes a are vanishing in the limit e ! 0. The solution of (81) is also a diagonal matrix:

0

1

B0 B B B0 b Tðx3 Þ ¼ B B0 B B @0 0

0

0

Tb S ðx3 Þ 0

0

0

0

0

0

0

0

0

0

0

1

0

0

0

0

0

1

0

0

0

0

0

Tb P ðx3 Þ

Tb S ðx3 Þ 0

0

1 C C C C C; C C C A

where

d Tb S ijj? j2 b T S ðx3 Þ; ðx3 Þ ¼  dx3 2kS 2 d Tb P ijj? j b T P ðx3 Þ; ðx3 Þ ¼  dx3 2kP

Tb S ðx3 ¼ LÞ ¼ 1; Tb P ðx3 ¼ LÞ ¼ 1:

By using (77) and by taking an inverse Fourier transform in j? , we obtain the expression of the transmitted wave field and the result of Proposition 1. Finally, the solutions of the Schrödinger equations in a homogeneous medium can be computed explicitly and we obtain

  i 2  ð x ; j Þ exp  jj j L þ ij  x dj? expðixsÞdx; b ? ? ? ? inc;3 2kP ð2pÞ3   Z Z f1=2  ðx; j? Þ exp  i jj? j2 L þ ij?  x? dj? expðixsÞdx; j ¼ 1; 2: vSj3 ðs; x? Þ ¼ S 3 b inc;j 2kS ð2pÞ

vP33 ðs; x? Þ ¼

f1=2 P

Z Z

6. The derivation of the paraxial wave equations in the random case In this section we prove Proposition 2. We transform the two-point boundary value problem (69)–(71) into an initial value problem. This is very important in the random case so that we can deal exclusively with quantities that are adapted to the filtration of the driving processes ml ; mq , and mk (which is necessary for the application of a diffusion approximation theorem). This transformation is done by an invariant imbedding step in which we introduce transmission and reflection operators. The algebra is more complicated than in the homogeneous case since the random medium fluctuations involve

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

249

coupling not only between the six modes (as in the homogeneous case) but also between different j? -modes. This is why b e defined as in (68) is solution of the two-point boundary value we need to introduce matrix operators. The vector X problem:

be dX b e; ¼c A e ðx3 Þ X dx3

b e ðLÞ þ H0 X b e ð0Þ ¼ V b: HL X

b ðj? Þ as Here c A e ðx3 Þ is the matrix operator acting on six-dimensional vector fields Y

b ðj? Þ ¼ ½c A e ðx3 Þ Y

Z

c b ðj0 Þdj0 ; A e ðj? ; j0? ; x3 Þ Y ? ?

with the kernel

 h

b 0;0 þ 1 exp i x3 ðkP þ kS Þ A b 1;1 c A e ðj? ; j0? ; x3 Þ ¼ dðj?  j0? Þ A e2 e4 x





io 3 b 1;1 þ exp i x3 ðkP þ kS Þ A b 1;1 þ exp i x3 ðkP  kS Þ A b 1;1 þ exp i 4 ðkP  kS Þ A 4 4 e  e   e 

 1 1 b 2kP x3 b 2;0 2kP x3 b 2;0 0 x3 ~ exp i B B þ m j  j ; þ exp i P ? ? e ð2pÞ2 e2 e4 e4      

1 1 b 2kS x3 b 0;2 2kS x3 b 0;2 0 x3 ~ B B exp i m j  j ; þ exp i þ S ? ? 2 e ð2pÞ2 e e4 e4



1 1 x b 0;0 1 1 x b 0;0 b j  j0 ; 3 B b j  j0 ; 3 B þ : m m þ e ð2pÞ2 P ? ? e2 P e ð2pÞ2 S ? ? e2 S

ð82Þ

By applying Proposition B.1, we get that the reflected and transmitted modes are given by

b ; ^ej ðx3 ¼ 0Þ ¼ ½ c a R e ðx3 ¼ 0Þ V 2j1

j ¼ 1; 2; 3;

ð83Þ

b ; bj ðx3 ¼ LÞ ¼ ½ T ðx3 ¼ 0Þ V 2j

j ¼ 1; 2; 3;

ð84Þ

^e

ce

where the reflection and transmission matrix-operator kernels are solution of the initial value problem

dc Re b c b; ¼ ðI Re c A e ðx3 Þ c Re; c R e ðx3 ¼ LÞ ¼ I H 0 Þc dx3 dc Te b; H0 c ¼ c T ec Re; c T e ðx3 ¼ LÞ ¼ I A e ðx3 Þ c dx3

ð85Þ ð86Þ

with

b ðj? ; j0 Þ ¼ dðj?  j0 ÞI; I ? ?

c H 0 ðj? ; j0? Þ ¼ dðj?  j0? ÞH0 :

Explicitly, the transmitted wave components are, for j ¼ 1; 2:

Z Z Z X 3

fS1=2

vj3 ðt; e2 x? ; x3 ¼ LÞ ¼

ð2pÞ3

^ ðx; j0 Þ expðij?  x? Þdj0 dj? Tb e2j;2k ðx; j? ; j0? ; x3 ¼ 0Þb inc;k ? ?

k¼1

  kS  exp i 4 ðL  cS tÞ dx;

e

v33 ðt; e2 x? ; x3 ¼ LÞ ¼

Z Z Z X 3

f1=2 P ð2pÞ

3

j ¼ 1; 2;

Tb e6;2k ðx; j? ; j0? ; x3 ¼ 0Þ b b inc;k ðx; j0? Þ expðij?  x? Þdj0? dj?

k¼1

  kP  exp i 4 ðL  cP tÞ dx;

e

f1 S

f1 P

vj ¼ vj3 and v3 ¼ v33 while the other components are obtained by (52). We first note that the rapid phase in x will give a localization in time of the transmitted waves (provided we show that c T e has a limit). Therefore we focus our attention on



 L þ e4 s; e2 x? ; x3 ¼ L ; j ¼ 1; 2; cS   L e v33 ðs; x? Þ ¼ v33 t ¼ þ e4 s; e2 x? ; x3 ¼ L : cP

vej3 ðs; x? Þ ¼ vj3 t ¼

250

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

The proof of the convergence follows closely the strategy adopted in [14], where the paraxial wave equation is obtained from the acoustic wave equations in the same distinguished limit. The proof is similar because the matrix operator c Ae has a diagonal block structure (with 2  2 blocks along the diagonal), which means that it has the same structure as the 2  2 operator encountered in the acoustic case. The main step of the proof consists in showing the convergence of the general moments of the transmitted wave components to the ones given by the limit system of stochastic partial differential Eqs. (34) and (35). For N 2 N; jr 2 f1; 2; 3g; mr 2 N; sr 2 R and x?r 2 R2 ; r ¼ 1; . . . ; N, the general moment of transmitted wave components

" # N Y mr e I ¼E vjr 3 ðsr ; x?r Þ e

r¼1

can be expressed as a sum of 3M multiple integrals, for M ¼

X

Ie ¼

MS

Ik ¼

MP

fS2 fP2

3M

Z



r¼1 mr :

Iek ;

k¼ðkr;h Þh¼1;...;mr ;r¼1;...;N 2f1;2;3g

e

PN

M

Z Y mr N Y

dj0?r;h dj?r;h dxr;h 

ð2pÞ r¼1 h¼1 " # Y E Tb e2jr ;2kr;h ðxr;h ; j?r;h ; j0?r;h ; 0Þ ;

Y

b b inc;kr;h ðxr;h ; j0?r;h Þ expðiðj?r;h  x?r  xr;h sr ÞÞ

r;h

r;h

for M P ¼ cardfkr;h ¼ 3; h ¼ 1; . . . ; mr ; r ¼ 1; . . . ; Ng and M S ¼ M  M P . Therefore, the convergence of the general moment of the transmitted wave components in the white-noise limit will follow from the convergence of the specific moments E½J e ð0Þ of the transmission matrix-operator kernels, where

J e ðx3 Þ ¼

M Y

Tb ejr ;kr ðxr ; j?r ; j0?r ; x3 Þ:

ð87Þ

r¼1

We call these moments ‘‘specific” because we restrict our attention to the case in which the frequencies xr ; r ¼ 1; . . . ; N, are T , but all distinct in (87). It is important to note that the transmission operator kernels c T e themselves do not converge to c only certain moments (expectations of products of components with distinct frequencies), which are those needed to ensure the convergence of the transmitted fields. We use diffusion approximation theorems in the same way as in [14] combined with the application of the averaging theorem in Appendix C as in the homogeneous case to obtain the limit of the expectation of J e ð0Þ as e ! 0. We find that this limit can be expressed as the expectation of the corresponding product of components of an ‘‘effective” transmission kernel in the following way:

" M Y

lim E½J e ð0Þ ¼ E e!0

# Tb jr ;kr ðxr ; j?r ; j0?r ; 0Þ :

r¼1

The expectation on the right-hand side is taken with respect to the following Itô-Schrödinger model for the transmission operator:

0 B B B B B 0 c T ðx; j? ; j? ; x3 Þ ¼ B B B B @

dðj?  j0? Þ 0

0 b TS

0 0 j0? Þ

0

0

0

0 0 0 dðj?  j0? Þ 0

0

0

0 0

0 0

0 0

0 b TS 0

0

0

0

0

dðj? 

0

1

C 0 C C 0 C C C; 0 C C C 0 A Tb P

where the kernels Tb S ðx; j? ; j0? ; x3 Þ and Tb P ðx; j? ; j0? ; x3 Þ are solution of 2

ijj0? j2 b k C 0;SS ð0Þ b T S ðx; j? ; j0? ; x3 Þdx3 T S ðx; j? ; j0? ; x3 Þdx3  S 8 2kS Z ikS b S ðj00  j0 ; x3 Þdj00 ; Tb S ðx; j? ; j00? ; x3 Þd B þ ? ? ? 2 2ð2pÞ

d Tb S ðx; j? ; j0? ; x3 Þ ¼ 

251

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254 2

ijj0? j2 b k C 0;PP ð0Þ b T P ðx; j? ; j0? ; x3 Þdx3 T P ðx; j? ; j0? ; x3 Þdx3  P 8 2kP Z ikP b P ðj00  j0 ; x3 Þdj00 ; þ Tb P ðx; j? ; j00? ; x3 Þd B ? ? ? 2ð2pÞ2

d Tb P ðx; j? ; j0? ; x3 Þ ¼ 

starting from Tb S ðx; j? ; j0? ; x3 ¼ LÞ ¼ dðj?  j0? Þ and Tb P ðx; j? ; j0? ; x3 ¼ LÞ ¼ dðj?  j0? Þ. Here we use the standard Itô stob P ÞT is the partial Fourier transform of the field B ¼ ðBS ; BP ÞT . It has the folb ¼ ðB bS; B chastic integral and the Brownian field B lowing operator-valued spatial covariance:

h i b 0 ðj? Þdðj? þ j0 Þ; b ? ; x3 Þ B b T ðj0 ; x0 Þ ¼ minfx3 ; x0 gð2pÞ2 C E Bðj ? 3 3 ?

ð88Þ

where

b 0 ðj? Þ ¼ C

Z R2

C0 ðx? Þ expðij?  x? Þdx? :

ð89Þ

By considering the transmission kernels in the original spatial variables:

T S;P ðx; x? ; x0? ; x3 Þ ¼

1

Z Z

ð2pÞ2

  exp iðj?  x?  j0?  x0? Þ Tb S;P ðx; j? ; j0? ; x3 Þdj? dj0? ;

we obtain that they satisfy the system (34) and (35) where we have used the Stratonovich integral instead of Itô integral.

7. Conclusion In this paper we have used invariant imbedding and limit theorems to obtain the paraxial equations for elastic waves in random media in the parabolic white-noise regime. We have shown that the paraxial systems for the pressure waves and shear waves have the form of a system of Schrödinger equations that are driven by a pair of correlated Brownian fields. We have identified the structure of the correlation between the Brownian fields and between the wave fields and we have clarified the form of the diffraction operators. These equations give a new and efficient way of modeling elastic wave propagation in a consistent manner via averaging and diffusion approximation results. Owing to the central role elastic wave propagation have in for instance geophysics we expect these result to be relevant in a number of applications, in particular for analyzing schemes for imaging as well as for numerical wave propagation.

Appendix A. Expressions of the matrices  are obtained from A b by substitution ij1 ! @ and ij2 ! @ . The matrices A @x1 @x2

0

b 1;1 A

B B B B ¼B B B B @

b 1;1

A

0

0

0

0

0

0 ij1 aa

0

0

0

0

0

0

0

0 ij2 aa

0

0

0

0

0

0 0

0 0 1 ij1 ab 0 0 0C C C ij2 ab 0 C C; 0 0C C C 0 0A

0 ij2 ac 0 0

with aa ¼

1=2

fP

1=2

fS



k0 4l0

0

0

0

0

0

0

0

0

0

0

0

0 0 0 ij2 ac

0 ij1 ac

0

0

0

ij1 ac 0 0 B0 B B B0 ¼ B B0 B B @0

0

1

0

0

0C C C 0C C; 0C C C 0A

b 1;1 A

0

0

0

fP

0

B0 0 B B B0 0 ¼B B0 0 B B @ 0 ij1 ac

1=2 1=2 fS fP l0 þk0 k0  2ffSP ; ab ¼ 1=2 þ 2ffSP , and ac ¼ 1=2 , 4l 2l þk0 fS

0

0

0

1

0

0

0 ij1 aa

0 0

0 0

0

0

0 ij2 ac

0 0 C C C 0 ij2 aa C C; 0 0 C C C 0 0 A

0

0

0

0

252

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

0

b 0;0 A

1

B B0 B 2 B0 ijj? j B B ¼ 2kS B B0 B B0 @

0 0 0 0

b 0;2 A

b 2;0 A

0

0 0

0

0

0

0

0

j21

B B 0 B B jj l0 þ k 0 B B 1 2 l0 B B 0 B B @ 0

0 0

0 0 0

j1 j2

0

j21

0

j1 j2

0

j22

0

j1 j2

0

j22

0

0

0

C 0 0C C C 0 0C C C 0 0C C C 0 0A

0

0

0

0 0

0 0

1

0

0 0

1

C 0C C C 0 0 0 0 0C C; C 0 0 0 0 0C C 0 0 0 1 0 C A 0 0 0

0

0 0 0 0

0 1

0 0 0 0 0 0

B B1 B 2 B0 ijj? j B B ¼ 2kS B B0 B B0 @

1

C 0 0C C C 1 0 0 0C C i C 0 1 0 0 C 2kS C 0 0 0 0C A

1 0

0 0 0 0 B B0 B B 2 B0 ijj? j k0 B  2kP 2l0 þ k0 B B0 B B0 @ 0

0

0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0

0 0 0 0 0 0 B B0 B B B0 ijj? j2 k0 B ¼ 2kP 2l0 þ k0 B B0 B B0 @

1 0

C 0C C C 0C Cþ i C 0 C 2kS C 0C A

0 B 2 B j1 B B 0 l0 þ k0 B B B l0 B j1 j2 B B 0 @

0 0

0

0 0 0 0 0 0 0 0 0 0 0 0 0 0 0

0 0

0 0 0

0

0 0 0

1

C j1 j2 0 0 0 C C C 0

0 0 0C C; C 0 0 0C C 0 0 0C A

2 2

0

j

0

0

0

0

0 0 0

1

C 0 0C C C 0 0C C; C 0 0C C 0 0C A

0 0 0 0 1 0 b 0;2 ¼ ð A b 0;2 ÞH ; A b 1;1 ¼ ð A b 1;1 ÞH , and A b 2;0 ¼ ð A b 2;0 ÞH , where H denotes the conjugate transpose. The matrices B b j;k are given and A by

0

b 0;2 B

0 B1 B B ikS B B0 ¼ B 2 B0 B B @0 0

b 0;0 B S

0 1

B0 B B ikS B B0 ¼ B 2 B0 B B @0 0

0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0

1

0C C C 0C C C; 0C C C 0A

0 0 0 0 0 0

0

1 0

0

0 0

0

0

0

0

0

b 2;0 B

1

0C C C 0C C C; 0C C C 0A

0

1

0

0 1 0

0

0

0

0

0

0

0

0 0

0 B0 B B ikP B B0 ¼ B 2 B0 B B @0

b 0;0 B P

0 0 0 0 0

1

0 0 0 0 0C C C 0 0 0 0 0C C C; 0 0 0 0 0C C C 0 0 0 0 0A

0 0 0 0 B0 B B ikP B B0 ¼ B 2 B0 B B @0

0 0 1 0

1

0 0 0 0

0

0 0 0 0

0 C C C 0 C C C; 0 C C C 0 A

0 0 0 0 0 0 0 0 0 0 0 1

0 0 0 0 0 1

b 0;2 ¼ ð B b 0;2 ÞH and B b 2;0 ¼ ð B b 2;0 ÞH . and B

Appendix B. An invariant imbedding theorem Let us consider the two-point boundary value problem:

dX ðzÞ ¼ AðzÞXðzÞ; dz

XðzÞ 2 Rm ;

ðB:1Þ

253

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

with the boundary condition HL XðLÞ þ H0 Xð0Þ ¼ V 0 , where AðzÞ; HL and H0 are m  m-matrices and V 0 is an m-dimensional vector. Assume that HL þ H0 is invertible. In this linear framework the invariant imbedding approach leads to the following proposition [2]. Proposition B.1. Let us define the matrices ðRðfÞÞL6f60 and ðQ ðz; fÞÞL6z6f60 as the solutions of the initial value problems:

dR ðfÞ ¼ AðfÞRðfÞ  RðfÞH0 AðfÞRðfÞ; df

L 6 f 6 0;

ðB:2Þ

starting from f ¼ L: Rðf ¼ LÞ ¼ ðHL þ H0 Þ1 , and

@Q ðz; fÞ ¼ Q ðz; fÞH0 AðfÞRðfÞ; @f

L 6 z 6 f 6 0;

ðB:3Þ

starting from f ¼ z: Q ðz; f ¼ zÞ ¼ RðzÞ. If the solution RðfÞ of the nonlinear initial value problem (B.2) exists up to f ¼ 0, then PðzÞ ¼ Q ðz; 0Þ is the solution of:

dP ðzÞ ¼ AðzÞPðzÞ; dz

0 6 z 6 L;

with HL PðLÞ þ H0 Pð0Þ ¼ I;

and consequently XðzÞ ¼ PðzÞV 0 is the solution of (B.1) with the boundary condition HL XðLÞ þ H0 Xð0Þ ¼ V 0 . Appendix C. An averaging theorem The following proposition can be found in [18,20]. We give here an elementary proof for consistency. Proposition C.1. Let us consider the solution Xe ¼ ðX ejk Þj;k¼1;...;N of the initial value problem

 ðqÞ    Q  dXe 1 X a z aðqÞ z FðqÞ ðXe Þ exp i 4 ¼ Fð0Þ ðXe Þ þ 2 þ GðqÞ ðXe Þ exp i 4 ; dz e q¼1 e e 2

2

2

ðC:1Þ

2

starting from Xe ðz0 Þ ¼ Xini . Here FðqÞ : RN ! RN and GðqÞ : RN ! RN are smooth functions and aðqÞ 2 R n f0g are distinct along with their sums and differences. Then Xe converges as e ! 0 to X the solution of the effective equation Q X dX i ½FðqÞ ðXÞ; GðqÞ ðXÞ; ¼ Fð0Þ ðXÞ þ ðqÞ dz a q¼1

Xðz0 Þ ¼ Xini ;

where ½;  stands for the Lie brackets defined by

½FðXÞ; GðXÞ ¼

N  X @FðXÞ

@X jk

j;k¼1

Gjk ðXÞ 

 @GðXÞ F jk ðXÞ : @X jk

Proof. The result can be obtained by a multi-scale expansion ansatz in powers of harmonic phases:

Xe ðzÞ ¼ Xð0Þ ðzÞ þ e2

 ðqÞ    Q  X a z aðqÞ z Xðq;1Þ ðzÞ exp i 4 þ Xðq;1Þ ðzÞ exp i 4

e

q¼1

þ e4

Q  X

Xðq;2Þ ðzÞ exp 2i

q¼1

þ e4

X 16q
þ e4



X 16q
 aðqÞ z

e4

e

  aðqÞ z þ Xðq;0Þ ðzÞ þ Xðq;2Þ ðzÞ exp 2i 4

e

  ðqÞ   ðqÞ  0 0 0 ða þ aðq Þ Þz ða  aðq Þ Þz ðq;q0 ;1;1Þ þ X Xðq;q ;1;1Þ ðzÞ exp i ðzÞ exp i 4 4

e

e

     0 0 0 ðaðqÞ þ aðq Þ Þz ðaðqÞ  aðq Þ Þz ðq;q0 ;1;1Þ Xðq;q ;1;1Þ ðzÞ exp i ðzÞ exp i þ X þ Oðe6 Þ: 4 4

e

e

ðC:2Þ

Recall that aðqÞ along with their sums and differences are distinct so that these constitute different harmonics. We now substitute the ansatz (C.2) in (C.1). By subsequently collecting the terms of order 0ðe2 Þ and 0ðe0 Þ and identifying the terms with the same harmonic (i.e. the same rapid phase), we obtain:

254

J. Garnier, K. Sølna / Wave Motion 46 (2009) 237–254

 iaðqÞ Xðq;1Þ ¼ GðqÞ ðX0 Þ;

iaðqÞ Xðq;1Þ ¼ FðqÞ ðX0 Þ;

Q X X dXð0Þ @FðqÞ ðXð0Þ Þ ðq;1Þ X @GðqÞ ðXð0Þ Þ ðq;1Þ X jk þ X jk ; ¼ Fð0Þ ðXð0Þ Þ þ @X jk @X jk dz q¼1 j;k j;k

X @GðqÞ ðXð0Þ Þ ðq;1Þ X jk ; @X jk j;k

 2iaðqÞ Xðq;2Þ ¼ 0

0

iðaðqÞ þ aðq Þ ÞXðq;q ;1;1Þ ¼

 iðaðqÞ þ aðq Þ ÞXðq;q ;1;1Þ ¼ 0

0

iðaðqÞ  aðq Þ ÞXðq;q ;1;1Þ ¼ 0

0

X @FðqÞ ðXð0Þ Þ ðq;1Þ X jk ; @X jk j;k

X @FðqÞ ðXð0Þ Þ ðq0 ;1Þ X @Fðq0 Þ ðXð0Þ Þ ðq;1Þ X jk þ X jk ; @X jk @X jk j;k j;k

0

0

2iaðqÞ Xðq;2Þ ¼

X @GðqÞ ðXð0Þ Þ ðq0 ;1Þ X @Gðq0 Þ ðXð0Þ Þ ðq;1Þ X jk þ X jk ; @X jk @X jk j;k j;k

X @FðqÞ ðXð0Þ Þ ðq0 ;1Þ X @Gðq0 Þ ðXð0Þ Þ ðq;1Þ X jk þ X jk ; @X jk @X jk j;k j;k

iðaðqÞ þ aðq Þ ÞXðq;q ;1;1Þ ¼

X @GðqÞ ðXð0Þ Þ ðq0 ;1Þ X @Fðq0 Þ ðXð0Þ Þ ðq;1Þ X jk þ X jk : @X jk @X jk j;k j;k

Substituting the expressions of Xðq;1Þ and Xðq;1Þ given by the first equation into the second one gives the desired result. The other equations give the expressions of the corrective terms. h References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30]

F. Bailly, J.-F. Clouet, J.-P. Fouque, Parabolic and white noise approximation for waves in random media, SIAM J. Appl. Math. 56 (1996) 1445–1470. R. Bellman, G.M. Wing, Introduction to Invariant Imbedding, Wiley, New York, 1975. P. Blomgren, G. Papanicolaou, H. Zhao, Super-resolution in time-reversal acoustics, J. Acoust. Soc. Am. 111 (2002) 230–248. J.F. Claerbout, Imaging the Earths Interior, Blackwell Scientific Publications, Palo Alto, 1985. J. Corones, B. DeFacio, R.J. Krueger, Parabolic approximations to the time-independent elastic wave equations, J. Math. Phys. 23 (1982) 577–586. D. Dawson, G. Papanicolaou, A random wave process, Appl. Math. Opt. 12 (1984) 97–114. A.T. de Hoop, Handbook of Radiation and Scattering of Waves, Academic Press, London, 1995 (Chapter 10) . M.V. de Hoop, A.T. de Hoop, Elastic wave up/down decomposition in inhomogeneous and anisotropic media: an operator approach and its approximations, Wave Motion 20 (1994) 57–82. A. Fannjiang, Self-averaging scaling limits for random parabolic waves, Arch. Rational Mech. Anal. 175 (2005) 343–387. A. Fannjiang, K. Sølna, Propagation and time-reversal of wave beams in atmospheric turbulence, SIAM Multiscale Model. Simul. 3 (2005) 522–558. J.-P. Fouque, J. Garnier, G. Papanicolaou, K. Sølna, Wave Propagation and Time Reversal in Randomly Layered Media, Springer, New York, 2007. J.-P. Fouque, G. Papanicolaou, Y. Samuelides, Forward and Markov approximation: the strong-intensity-fluctuations regime revisited, Waves Random Media 8 (1998) 303–314. J. Garnier, K. Sølna, Effective transport equations and enhanced backscattering in random waveguides, SIAM J. Appl. Math. 68 (2008) 1574–1599. J. Garnier, K. Sølna, Coupled paraxial wave equations in random media in the white-noise regime, Ann. Appl. Probab. 19 (2009) 318–346. J.A. Hudson, A parabolic approximation for elastic waves, Wave Motion 2 (1980) 207–214. W. Kohler, G. Papanicolaou, B. White, Localization and mode conversion for elastic waves in randomly layered media I, Wave Motion 23 (1996) 1–22. W. Kohler, G. Papanicolaou, B. White, Localization and mode conversion for elastic waves in randomly layered media II, Wave Motion 23 (1996) 181– 201. J. Kurzweil, J. Jarnı´k, Limit processes in ordinary differential equations, J. Appl. Math. Phys. 38 (1987) 241–256. T. Landers, J.F. Claerbout, Numerical calculations of elastic waves in laterally inhomogeneous media, J. Geophys. Res. 77 (1972) 1476–1482. W. Liu, Averaging theorems for highly oscillatory differential equations and iterated Lie brackets, SIAM J. Control Opt. 35 (1997) 1989–2020. J.J. McCoy, A parabolic theory of stress wave propagation through inhomogeneous linearly elastic solids, J. Appl. Mech. 44 (1977) 462–468. G. Papanicolaou, L. Ryzhik, K. Sølna, Statistical stability in time reversal, SIAM J. Appl. Math. 64 (2004) 1133–1155. G. Papanicolaou, L. Ryzhik, K. Sølna, Self-averaging from lateral diversity in the Itô-Schrödinger equation, SIAM Multiscale Model. Simul. 6 (2007) 468– 492. L. Ryzhik, J. Keller, G. Papanicolaou, Transport equations for elastic and other waves in random media, Wave Motion 24 (1996) 327–370. S.A. Shapiro, P. Hubral, Elastic waves in random media: fundamentals of seismic stratigraphic filtering, Lecture Notes in Earth Sciences, Springer, Berlin, 1999. J.W. Strohbehn (Ed.), Laser Beam Propagation in the Atmosphere, Springer, Berlin, 1978. F. Tappert, The parabolic approximation method, in: J.B. Keller, J.S. Papadakis (Eds.), Wave Propagation and Underwater Acoustics, Springer, Berlin, 1977, pp. 224–287. B. Ursin, Review of elastic and electromagnetic wave propagation in horizontally layered media, Geophysics 48 (1983) 1063–1081. K. Viswanathan, R.N. Biswas, Elastic waves from a circular disc source, Pure Appl. Geophys. 95 (1972) 27–32. S.C. Wales, J.J. McCoy, A comparison of parabolic wave theories for linearly elastic solids, Wave Motion 5 (1983) 99–113.