Transport of power in random waveguides with turning points

# Transport of power in random waveguides with turning points

Liliana Borcea111Department of Mathematics, University of Michigan, Ann Arbor, MI 48109. borcea@umich.edu and derekmw@umich.edu    Josselin Garnier222Centre de Mathématiques Appliquées, Ecole Polytechnique, 91128 Palaiseau Cedex, France. josselin.garnier@polytechnique.edu    Derek Wood111Department of Mathematics, University of Michigan, Ann Arbor, MI 48109. borcea@umich.edu and derekmw@umich.edu
###### Abstract

We present a mathematical theory of time-harmonic wave propagation and reflection in a two-dimensional random acoustic waveguide with sound soft boundary and turning points. The boundary has small fluctuations on the scale of the wavelength, modeled as random. The waveguide supports multiple propagating modes. The number of these modes changes due to slow variations of the waveguide cross-section. The changes occur at turning points, where waves transition from propagating to evanescent or the other way around. We consider a regime where scattering at the random boundary has significant effect on the wave traveling from one turning point to another. This effect is described by the coupling of its components, the modes. We derive the mode coupling theory from first principles, and quantify the randomization of the wave and the transport and reflection of power in the waveguide. We show in particular that scattering at the random boundary may increase or decrease the net power transmitted through the waveguide depending on the source.

Key words. mode coupling, turning waves, scattering, random waveguide.

AMS subject classifications. 60F05, 35Q99, 78M35.

## 1 Introduction

Guided waves have been studied extensively due to their numerous applications in electromagnetics [13], optics and communications [27, 29], quantum waveguides [17], ocean acoustics [32], etc. The classic waveguides, with straight boundaries and filled with media that are either homogeneous or do not vary along the direction of propagation of the waves, are well understood. The wave equation in such waveguides is separable and the solution is a superposition of independent modes, which are propagating and evanescent modes and, in the case of penetrable boundary, radiation modes. When the geometry of the waveguide varies, or the medium filling it is heterogeneous, the modes are coupled. Examples of mathematical studies that account for mode coupling and lead to numerical methods that model such waveguides can be found in [11, 12, 15, 14, 23]. We are interested in mode coupling due to small random perturbations of the waveguide, which can be quantified more explicitly using asymptotic analysis.

Random models are useful for waveguides with rough boundaries and filled with composite media that vary at small scale, comparable to the wavelength. Such variations are typically small and unknown, so they introduce uncertainty in the model of wave propagation. The random models of the boundary and the wave speed are used to quantify how this uncertainty propagates to the uncertainty of the solution of the wave equation. This is useful information in applications like imaging [10, 7, 9, 2]. The mode coupling theory in waveguides filled with random media has been developed in [24, 16, 19, 21, 20] for sound waves and in [27, 4] for electromagnetic waves. The theory has also been extended to waveguides with random perturbations of straight boundaries [5, 7, 22].

In this paper we develop a mode coupling theory in random waveguides with turning points. We consider for simplicity a two-dimensional acoustic waveguide with sound soft boundary, but the theory can be extended to three dimensions and to electromagnetics. The waveguide has a slowly bending axis, a slowly changing opening and a randomly perturbed boundary. The slow variations occur on a long scale with respect to the wavelength, whereas the random perturbations are at a scale similar to the wavelength.

In the absence of the random perturbations, waves in slowly changing waveguides can be analyzed with a local decomposition in modes that are approximately independent [3, 31]. This is known as the adiabatic approximation [29, Section 19-2], and the result differs from that in waveguides with straight boundaries in one important aspect: The change in the opening of the waveguide causes the number of propagating modes to increase or decrease by at locations called turning points. Modes transition there from propagating to evanescent or the other way around, and due to energy conservation, the impinging propagating mode is turned back i.e., is reflected. Turning waves in slowly changing waveguides are studied mathematically in [6]. A recent study of their interaction with a random boundary is given in [8], for the case of weak random fluctuations that affect only the turning modes. Here we extend the results in [8] to stronger fluctuations, that couple all the waveguide modes.

Starting from the wave equation, and using the separation of scales of variation of the waveguide, we derive an asymptotic model for the wave field that accounts for coupling of all the modes, propagating and evanescent. This coupling is described by a stochastic system of differential equations for the random mode amplitudes, endowed with initial conditions that model the source excitation and radiation conditions. The excitation is due to a point source, but due to the linearity of the equation, other source excitations can be handled by superposition. We obtain an extension of the diffusion approximation theorem proved in [28] to carry out the asymptotic analysis of the mode amplitudes. The result simplifies when the random fluctuations are smooth, because the forward and backward going components of the propagating modes become independent. This is known as the forward scattering approximation, and applies to propagation between the turning points. At the turning points there is strong coupling of the components of the turning waves, described by random reflection coefficients. With the diffusion approximation theorem, and in the forward scattering approximation regime, we quantify the net effect of the random boundary on the transmitted and reflected power in the waveguide. This is the main result of the paper, and shows that the random boundary can be useful for increasing the transmitted power.

The paper is organized as follows: We begin in section LABEL:sect:formulation with the formulation of the problem, and the derivation of the asymptotic model for the wave field. Then we give in section LABEL:sect:wavedec the mode decomposition and derive a closed system of stochastic differential equations for the random amplitudes of the propagating modes between turning points. These equations are complemented by source excitation conditions, radiation conditions, as well as continuity and reflection conditions at the turning points. The asymptotic limit of the solution of these equations and the forward scattering approximation are in section LABEL:sect:results. We use these results to quantify the transmitted and reflected power in the waveguide in section LABEL:sect:netscat. The diffusion approximation theorem used to carry the asymptotic limit is stated and proved in section LABEL:app:adif. We end with a summary in section LABEL:sect:sum.

## 2 Formulation of the problem

In this section we give the mathematical model for time-harmonic waves in a random waveguide with variable cross-section and bending axis. We begin in section LABEL:sect:form1 with the setup, and describe the scaling in section LABEL:sect:form2 in terms of a small, dimensionless parameter . We use it in section LABEL:sect:form3 to write the wave problem in a form that can be analyzed in the asymptotic limit .

### 2.1 Setup

We consider a two-dimensional acoustic waveguide with sound-soft boundary. The waveguide occupies the semi-infinite domain , bounded above and below by two curves and , as shown in Figure LABEL:fig:setup. The top boundary is perturbed by small random fluctuations about the curve shown in the figure with the dotted line. The axis of the waveguide is at half the distance between and . It is a smooth curve parametrized by the arc length , that bends slowly, meaning that its tangent and curvature vary on a scale which is large with respect to the waveguide width . The width function has bounded first two derivatives, and to avoid complications in the analysis of scattering of the waves at the random boundary, we also assume that it is monotonically increasing.

Because of the changing geometry, it is convenient to use orthogonal curvilinear coordinates with axes along and , where is the unit vector orthogonal to , pointing toward the upper boundary. For any , written henceforth as , we have

 \itbfx=\itbfx∥(z)+r\itbfn(zL), \hb@xt@.01(2.1)

where is along the waveguide axis at arc length , satisfying

 ∂z\itbfx∥(z)=τ(zL), \hb@xt@.01(2.2)

and is the coordinate in the normal direction. The domain is the set

 Ω={(r,z): z∈R, r∈(r−(z),r+(z))}, \hb@xt@.01(2.3)

where

 r−(z)=−D(z/L)2, \hb@xt@.01(2.4)

is at the bottom boundary and

 r+(z)=D(z/L)2[1+1(−ZM,ZM)(z)σν(zℓ)], \hb@xt@.01(2.5)

is at the randomly perturbed top boundary . The perturbation is modeled by the random process and it extends over the interval , the support of the indicator function , where is a long scale needed to impose outgoing boundary conditions on the waves.***In practice may be chosen based on the duration of the observation time of the wave, using the hyperbolicity of the wave equation in the time domain. The single frequency wave analyzed in this paper is the Fourier transform of the time dependent wave field. We let the boundaries of the waveguide be straight and parallel for .

The random process is stationary with zero mean

 E[ν(ζ)]=0, \hb@xt@.01(2.6)

and auto-correlation function

 R(ζ)=E[ν(0)ν(ζ)]. \hb@xt@.01(2.7)

We assume that is mixing, with rapidly decaying mixing rate, as defined for example in [28, section 2], and it is bounded, with bounded first two derivatives, almost surely. This implies in particular that is integrable and has at least four bounded derivatives. We normalize by

 R(0)=1,∫∞−∞dζR(ζ)=O(1), \hb@xt@.01(2.8)

so that in (LABEL:eq:coord3) is the standard deviation of the fluctuations of , and quantifies their correlation length.

The waves are generated by a point source at , which emits a complex signal at frequency . We take the origin of at the source, so that . The waveguide is filled with a homogeneous medium with wave speed , and the wave field satisfies the Helmholtz equation

 Δp(ω,\itbfx)+k2p(ω,\itbfx)=f(ω)δ(\itbfx−\itbfx⋆),\itbfx=(r,z)∈Ω, \hb@xt@.01(2.9)

where is the wavenumber. In curvilinear coordinates this takes the form

 ⎡⎢ ⎢ ⎢⎣∂2r−1Lκ(zL)∂r1−rLκ(zL)+∂2z[1−rLκ(zL)]2+rL2κ′(zL)∂z[1−rLκ(zL)]3+k2⎤⎥ ⎥ ⎥⎦p(ω,r,z) =∣∣∣1−r⋆Lκ(0)∣∣∣−1f(ω)δ(z)δ(r−r⋆), \hb@xt@.01(2.10)

as shown in appendix LABEL:ap:CurvCoord, where is the derivative of the curvature . The sound soft boundary is modeled by the homogeneous Dirichlet boundary conditions

 p(ω,r+(z),z)=p(ω,r−(z),z)=0, \hb@xt@.01(2.11)

and at points with we have radiation conditions that state that is outgoing and bounded.

### 2.2 Scaling

There are four length scales in the problem: The wavelength , the width of the waveguide , the scale of the slow variations of the waveguide, and the correlation length of the random fluctuations of the boundary . They satisfy

 L≫D∼λ∼ℓ, \hb@xt@.01(2.12)

where denotes “of the same order as”, and we model the separation of scales using the dimensionless parameter

 ε=ℓL,0<ε≪1. \hb@xt@.01(2.13)

Our analysis of the wave field is in the asymptotic limit .

As shown in section LABEL:sect:wavedec, the ratio of and defines , the number of propagating components of the wave, called modes, where denotes the integer part. The assumption in (LABEL:eq:sc1) means that

 Nmin≤N(z)≤Nmax, \hb@xt@.01(2.14)

for all , where and are natural numbers, independent of .

The scales and are of the same order in (LABEL:eq:sc1) so that the waves interact efficiently with the random fluctuations of the boundary. This interaction, called cumulative scattering, randomizes the wave field as it propagates in the waveguide. The distance from the source at which the randomization occurs depends on the standard deviation of the fluctuations. We scale as

 σ=√ε~σ,~σ=O(1), \hb@xt@.01(2.15)

so that we observe the randomization at distances .

The scaled variables are defined as follows: The arc length is scaled by ,

 ~z=zL, \hb@xt@.01(2.16)

and the similar lengths , and are scaled by , to obtain

 ~D(~z)=D(z/L)ℓ,  ~r=rℓ,  ~k=kℓ. \hb@xt@.01(2.17)

We also introduce the scaled bound of the support of the random fluctuations, which is a large number, independent of .

### 2.3 Asymptotic model

Let us multiply equation (LABEL:eq:we3) by and use the scaling relations (LABEL:eq:sc4)-(LABEL:eq:sc5p). Dropping the tilde to simplify notation, because all variables are scaled henceforth, we obtain

 [∂2z+(1−εrκ(z))2ε2(∂2r+k2)−κ(z)(1−εrκ(z))ε∂r+εrκ′(z)(1−εrκ(z))∂z]p(ω,r,z) =f(ω)[1−εr⋆κ(0)]εδ(r−r⋆)δ(z), \hb@xt@.01(2.18)

with homogeneous Dirichlet boundary conditions (LABEL:eq:we2) at

 r−(z)=−D(z)2,  r+(z)=D(z)2[1+1(−ZM,ZM)(z)√εσν(zε)], \hb@xt@.01(2.19)

and appropriate radiation conditions for . These equations define the asymptotic model for the wave field, and we wish to analyze it in the limit .

The boundary has dependent fluctuations, so to ensure that the boundary conditions are satisfied at all orders of , we change variables to

 r=ρ+[2ρ+D(z)]4√εσν(zε), \hb@xt@.01(2.20)

for , and denote the transformed wave field by

 pε(ω,ρ,z)=p(ω,ρ+(2ρ+D(z))4√εσν(zε),z). \hb@xt@.01(2.21)

At there are no fluctuations so the transformation is the identity . We use the same notation for the wave field at all , and analyze it separately in the regions with the random fluctuations and without. The results are connected by continuity at .

The change of variables (LABEL:eq:AM3) makes the boundary conditions independent of ,

 pε(ω,±D(z)2,z)=0, \hb@xt@.01(2.22)

and maps the random fluctuations to the differential operator in the wave equation. Explicitly, we show in appendix LABEL:ap:DerAs that the wave equation becomes

 ∞∑j=0εj/2−1Lεjpε(ω,ρ,z)=ˆf(ω)[1+O(√ε)]δ(ρ−ρ⋆)δ(z), \hb@xt@.01(2.23)

for and , with the leading term in the expansion of the operator

 Lε0=(ε∂z)2+∂2ρ+k2. \hb@xt@.01(2.24)

This is the Helmholtz operator in a perfect waveguide, with straight and parallel boundaries. The random fluctuations appear in the first perturbation operator,

 \hb@xt@.01(2.25)

The second perturbation operator has a deterministic part, due to the curvature of the axis of the waveguide, and a random part, quadratic in the random fluctuations,

 Lε2= −κ(z)[2ρ(∂2ρ+k2)+∂ρ]+σ24{3ν2(zε)+[2ρ+D(z)]24ν′2(zε)}∂2ρ +[2ρ+D(z)]σ24{ν(zε)ν′(zε)ε∂2ρz+[ν′2(zε)+12ν′′(zε)ν(zε)]∂ρ}. \hb@xt@.01(2.26)

The remaining operators in the asymptotic series in (LABEL:eq:AM5) depend on higher powers of the fluctuations , but play no role in the limit .

By assumption, there are no variations of the waveguide at , so the operator in the left hand side of (LABEL:eq:AM5) reduces to in this region.

## 3 Mode decomposition and coupling

To analyze the solution of the wave equation (LABEL:eq:AM5) with boundary conditions (LABEL:eq:AM6) in the limit , we begin in section LABEL:sect:dec1 with the mode decomposition of the wave field. The modes are special solutions of the wave equation with operator (LABEL:eq:AM7). They represent propagating and evanescent waves which are coupled by the perturbation operators (LABEL:eq:AM8)-(LABEL:eq:AM9), as explained in section LABEL:sect:dec2. We are interested in the propagating modes, which are left and right going waves with random amplitudes satisfying a stochastic system of equations derived in section LABEL:sect:dec3. It is this system that we analyze in the asymptotic limit to quantify the cumulative scattering effects in the waveguide.

### 3.1 Mode decomposition

The second term in (LABEL:eq:AM7) is the Sturm-Liouville operator acting on functions that vanish at , for any given . Its eigenvalues are real and distinct

 λj(z)=k2−μ2j(z),μj(z)=πjD(z),  j=1,2,… \hb@xt@.01(3.1)

and the eigenfunctions

 yj(ρ,z)=√2D(z)sin[(2ρ+D(z))2μj(z)], \hb@xt@.01(3.2)

form an orthonormal basis in . We use this basis to decompose the solution of (LABEL:eq:AM7) in one dimensional waves called modes, for each ,

 pε(ω,ρ,z)=∞∑j=1pεj(ω,z)yj(ρ,z). \hb@xt@.01(3.3)

As shown in appendix LABEL:ap:ModeCoup, the modes can be written as

 pεj(ω,z)=uεj(ω,z)[1+O(√ε)], \hb@xt@.01(3.4)

with satisfying a coupled system of one dimensional wave equations that we now describe:

In the perturbed section of the waveguide, satisfies

 1ε[(ε∂z)2+k2−μ2j(z)]uεj(ω,z)+σ√εμ2j(z)ν(zε)uεj(ω,z)+σ2gεj(ω,z)uεj(ω,z) ≈Cεj(ω,z)+f(ω)yj(ρ⋆,0)δ(z), \hb@xt@.01(3.5)

where the approximation indicates that we dropped lower order terms that have no contribution in the limit . The coefficient in the left hand side is

 gεj(ω,z)=−34μ2j(z)ν2(zε)−[(πj)212+116]ν′2(zε), \hb@xt@.01(3.6)

and

 Cεj(ω,z) =∞∑q=1,q≠j[σΓjq√εν′′(zε)+σ2γjq(zε)+γojq(z)]uεq(ω,z) \hb@xt@.01(3.7)

models the coupling between the modes. The leading coupling coefficients and are constants, given in equation (LABEL:eq:DC6) in appendix LABEL:ap:ModeCoup. The second order coefficients and are quadratic in , as described in equations (LABEL:eq:DC7)-(LABEL:eq:DC8), and the coefficients and , given in equations (LABEL:eq:DC9)-(LABEL:eq:DC10), are due to the slow changes in the waveguide.

In the region , where the waveguide has straight and parallel boundaries, the wave equation simplifies to

 1ε[(ε∂z)2+k2−μ2j(z)]uεj(ω,z)=0. \hb@xt@.01(3.8)

Depending on the index , its solution is either an outgoing propagating wave or a decaying evanescent wave. This wave is connected to the solution of (LABEL:eq:DC3) by the continuity of and at .

### 3.2 Random mode amplitudes

Equations (LABEL:eq:DC3) are perturbations of the wave equation with operator , where the perturbation term models the coupling of the modes. This coupling is similar to that in waveguides with randomly perturbed parallel boundaries, studied in [5, 10], but the slow variation of the waveguide introduces two differences: The first is the presence of the extra terms and in (LABEL:eq:DC5), given by (LABEL:eq:DC9)-(LABEL:eq:DC10), which turn out to play no role in the limit . The second difference is important, as it gives a dependent number

 N(z)=⌊kD(z)π⌋ \hb@xt@.01(3.9)

of mode indexes for which . These modes are oscillatory functions in , and represent left and right going waves. For indexes the modes are decaying evanescent waves.

#### 3.2.1 Turning points

The function (LABEL:eq:DC12) that defines the number of propagating modes is piecewise constant. Starting from the origin, where we denote the number of propagating modes by , the function (LABEL:eq:DC12) increases by at arc lengths , for , and decreases by at , for . The jump locations , ordered as

 −ZM<…

are the zeroes of the eigenfunctions (LABEL:eq:eigval), and are called turning points [26, 6]. We assume henceforth that the monotonically increasing satisfies

 D′(z(t)±)>0,∀t≥1, \hb@xt@.01(3.10)

so that the turning points are simple and isolated. Consistent with our scaling, they are spaced at order one scaled distances.

Between any two consecutive turning points and , where we set by convention the number of propagating modes equals the constant

 N(t−1)±=N(0)±(t−1). \hb@xt@.01(3.11)

This number is bounded above and below as in (LABEL:eq:ASN), with and so the bounds and on the indexes are

 \hb@xt@.01(3.12)

Beginning from the source location , which we assume is not a turning point, is defined as the unique, negative arc-length satisfying

 k=πN(t−1)−D(z(t)−),t=1,…,t−M, \hb@xt@.01(3.13)

where the uniqueness is due to the monotonicity of . Similarly, the jump location is defined as the unique, positive arc length satisfying

 k=π(N(t−1)++1)D(z(t)+),t=1,…,t+M. \hb@xt@.01(3.14)

The analysis of the modes is similar on the left and right of the source, so we focus attention in this section on a sector of the waveguide, for some , and simplify the notation for the number (LABEL:eq:DC13pp) of propagating modes

 N=N(t−1)−. \hb@xt@.01(3.15)

These modes are a superposition of right and left going waves, with random amplitudes that model cumulative scattering in the waveguide, as we explain in the next section.

#### 3.2.2 The left and right going waves

We decompose the propagating modes in left and right going waves, using a flow of smooth and invertible matrices ,

 (aεj(ω,z)bεj(ω,z))=Mε,−1j(ω,z)(uεj(ω,z)vεj(ω,z)), \hb@xt@.01(3.16)

where denotes the inverse of and

 vεj(ω,z)=−iε∂zuεj(ω,z),j=1,…,N. \hb@xt@.01(3.17)

We obtain from (LABEL:eq:DC3) that

 ∂z(aεj(ω,z)bεj(ω,z)) =Mε,−1j(ω,z){iε(01k2(ω)−μ2j(z)0)Mεj(ω,z)−∂zMεj(ω,z) −iCεj(ω,z)Mε,−1j(ω,z)(01), \hb@xt@.01(3.18)

and the decomposition is achieved by a flow that removes to leading order the large deterministic term in (LABEL:eq:FM6), the first line in the right hand side.

The matrix must have the structure

 Mεj(ω,z)=⎛⎝Mεj,11(ω,z)−¯¯¯¯¯¯¯¯¯¯¯¯Mεj,11(ω,z)Mεj,21(ω,z)¯¯¯¯¯¯¯¯¯¯¯¯Mεj,21(ω,z)⎞⎠, \hb@xt@.01(3.19)

where the bar denotes complex conjugate, so that the decomposition (LABEL:eq:FM5) conserves energy. The expression of the components in (LABEL:eq:res6) depends on the mode index, more precisely on the mode wave number denoted by

 βj(ω,z)=√k2−μ2j(z). \hb@xt@.01(3.20)

Note that is bounded away from zero for all , and it approaches zero as , for . This last mode is a turning wave which transitions from a propagating wave at to an evanescent wave at , as described in section LABEL:sect:TP2. Here we give the decomposition of the modes indexed by .

The entries of (LABEL:eq:res6) are defined by

 Mεj,11(ω,z) =1√βj(ω,z)exp[iε∫z0dz′βj(ω,z′)], Mεj,21(ω,z) =βj(ω,z)Mεj,11(ω,z), \hb@xt@.01(3.21)

for . This definition looks the same as in perfect waveguides with straight and parallel boundary, except that the mode wave number varies with . We obtain from (LABEL:eq:res6)-(LABEL:eq:FM8) that the determinant of is constant

 \hb@xt@.01(3.22)

so the matrix is invertible, and the decomposition (LABEL:eq:FM5) can be rewritten as

 uεj(ω,z)=1√βj(ω,z)[aεj(ω,z)eiε∫z0dz′βj(ω,z′)−bεj(ω,z)e−iε∫z0dz′βj(ω,z′)], \hb@xt@.01(3.23)

and

 ε∂zuεj(ω,z)=i√βj(ω,z)[aεj(ω,z)eiε∫z0dz′βj(ω,z′)+bεj(ω,z)e−iε∫z0dz′βj(ω,z′)]. \hb@xt@.01(3.24)

Note that equations (LABEL:eq:FM10)-(LABEL:eq:FM11) are just the method of variation of parameters for the perturbed wave equation satisfied by the -th mode. They decompose the mode in a right going wave with amplitude and a left going wave with amplitude . In perfect waveguides these amplitudes would be constant, meaning physically that the waves are independent. In our case the amplitudes are random fields, satisfying the system of stochastic differential equations

 ∂z(aεj(ω,z)bεj(ω,z))=Hεj(ω,z)(aεj(ω,z)bεj(ω,z))−i2Cεj(ω,z)(¯¯¯¯¯¯¯¯¯¯¯¯Mεj,11(ω,z)Mεj,11(ω,z)), \hb@xt@.01(3.25)

obtained by substituting (LABEL:eq:res6) and (LABEL:eq:FM8) in (LABEL:eq:FM6). Here is the matrix valued random process

 \hb@xt@.01(3.26)

with entries satisfying the relations

 Hε(ba)j(ω,z)=¯¯¯¯¯¯¯¯¯¯¯¯¯¯Hε(ab)j(ω,z),Hε(bb)j(ω,z)=¯¯¯¯¯¯¯¯¯¯¯¯¯¯Hε(aa)j(ω,z), \hb@xt@.01(3.27)

and taking the values

 Hε(aa)j(ω,z)≈i2βj(ω,z)[σ√εμ2j(z)ν(zε)+σ2gεj(ω,z)], \hb@xt@.01(3.28)

and

 Hε(ab)j(ω,z)≈[¯¯¯¯¯¯¯¯¯¯¯¯¯¯Hε(aa)j(ω,z)−∂zβj(ω,z)2βj(ω,z)]exp[−2iε∫z0dz′βj(ω,z′)]. \hb@xt@.01(3.29)

As before, the approximation means up to negligible terms in the limit .

Equations (LABEL:eq:FM12) show that the amplitudes of the -th mode are coupled to each other by the process , and to the other modes by , defined by the series (LABEL:eq:DC5). The first terms in this series involve the propagating waves , for , decomposed as in (LABEL:eq:FM10)-(LABEL:eq:FM10). We describe in the next two sections the turning and the evanescent waves.

#### 3.2.3 The turning waves

The mode indexed by transitions at from propagating to evanescent. This transition is captured by the matrix , which has the same structure as in (LABEL:eq:res6), but its entries are defined in terms of Airy functions [1, chapter 10]. This is because near the simple turning point , equation (LABEL:eq:DC3) for is a perturbation of Airy’s equation. We refer to [6, 29] for classic studies of turning waves in waveguides, and to [8] for an analysis of their interaction with the random boundary. The setup in [8] is the same as here, with the exception that we consider a larger standard deviation of the random fluctuations, to observe mode coupling in the waveguide.

We use the same as in [8], with entries

 \hb@xt@.01(3.30)

and

 MεN,21(ω,z)=−iε∂zMεN,11(ω,z), \hb@xt@.01(3.31)

for where is a small, positive number, independent of . We go slightly beyond the turning point to capture the transition of the wave to an evanescent one. The phase in definition (LABEL:eq:TW1) is given by the function

 ϕN(ω,z)=∫zz(t)−dz′√|k2−μ2N(z′)|, \hb@xt@.01(3.32)

evaluated at the source location , and the amplitude factor

 QN(ω,z)=∣∣3ϕN(ω,z)/2∣∣1/6∣∣k2−μ2N(z)∣∣1/4, \hb@xt@.01(3.33)

is shown in [8, Section 3.1.1] to be bounded, and at least twice continuously differentiable. The Airy functions and are evaluated at

 \hb@xt@.01(3.34)

where is of order one in the vicinity of the turning point, and it is much larger than one in the rest of the domain .

We recall from [8, Lemma 3.1] that the matrix is invertible, with constant determinant

 detMεN(ω,z)=2,∀z∈(z(t)−−δ,z(t−1)−), \hb@xt@.01(3.35)

so the decomposition (LABEL:eq:FM5) is well defined. Moreover, [8, Lemma 3.2] shows that at the expressions (LABEL:eq:TW1)-(LABEL:eq:TW2) become like (LABEL:eq:FM8),

 MεN,11(ω,z) =1√βN(ω,z)exp[iε∫z0dz′βN(ω,z′)]+O(ε), MεN,21(ω,z) =βN(ω,z)MεN,11(ω,z)+O(ε), \hb@xt@.01(3.36)

so the turning wave behaves just like any other propagating wave until it reaches the vicinity of the turning point from the right. On the left side of the turning point, at , the entries of grow exponentially, as given in [8, Lemma 3.3]. The wave is evanescent in this region, and must be decaying in order to have energy conservation. This is ensured by the radiation condition

 \hb@xt@.01(3.37)

which sets to zero the coefficients of the growing Airy function and its derivative in the expression of and at the end of the domain. We refer to [8, Section 3.1] for more details, and for the proof that the result does not depend on the particular choice of which is small, but bounded away from in the limit .

The evolution equation of the turning mode amplitudes is of the same form as in (LABEL:eq:FM12), with the following entries of the matrix (LABEL:eq:res24)-(LABEL:eq:res24p) indexed by ,

 Hε(aa)N(ω,z)≈i∣∣MεN,11(ω,z)∣∣22[σ√εμ2j(z)ν(zε)+σ2gεj(ω,z)], \hb@xt@.01(3.38)

and

 Hε(ab)N(ω,z)≈−i[MεN,11(ω,z)]22[σ√εμ2j(z)ν(zε)+σ2gεj(ω,z)]. \hb@xt@.01(3.39)

These expressions reduce to those in (LABEL:eq:FM14)-(LABEL:eq:FM14p) at , with the extra term involving in (LABEL:eq:FM14p) coming from an correction of the amplitudes, induced by the residual in (LABEL:eq:FM8E).

#### 3.2.4 Coupling with the evanescent waves

The modes indexed by in equations (LABEL:eq:DC3) are evanescent waves, with wavenumber

 βj(ω,z)=√μ2j(z)−k2. \hb@xt@.01(3.40)

We analyze these waves in appendix LABEL:ap:evanesc, and show that they can be expressed in terms of the propagating ones. Explicitly, at arc length , satisfying , we obtain that

 uεj(ω,z)≈ −σ√ε2βj(ω,z)N∑q=1∫∞−∞dξ[γ(e)jq(ω,zε+ξ)aεq(ω,z)√βq(ω,z)e1ε∫z0dz′βq(ω,z′)+iξβq(ω,z) −¯¯¯¯¯¯¯¯γ(e)jq(ω,zε+ξ)bεq(ω,z)√βq(ω,z)e−1ε∫z0dz′βq(ω,z′)−iξβq(ω,z)]e−|ξ|βj(ω,z), \hb@xt@.01(3.41)

where the approximation means that we neglect the terms that are negligible in the limit . Here we introduced the notation

 γ(e)jq(ω,zε+ξ)=Γjqν′′(zε+ξ)+iβq(ω,z)Θjqν′(zε+ξ)), \hb@xt@.01(3.42)

with coefficients and defined in (LABEL:eq:DC6), and recall that the bar denotes complex conjugate. For the derivative we have

 ε∂zuεj(ω,z)≈ −σ√ε2βj(ω,z)N∑q=1∫∞−∞dξ[θ(e)jq(ω,zε+ξ)aεq(ω,z)√βq(ω,z)e1ε∫z0dz′βq(ω,z′)+iξβq(ω,z) −¯¯¯¯¯¯¯¯θ(e)jq(ω,zε<