Wave propagation in randomly perturbed weakly coupled waveguides

TL;DR

This paper analyzes how random interface perturbations in coupled waveguides cause power leakage and blurring of power transfer, with implications for waveguide design and signal integrity.

Contribution

It provides a detailed mode amplitude analysis quantifying power transfer and leakage effects due to interface randomness in coupled waveguides.

Findings

Random interface fluctuations cause significant power leakage at large distances.

Power transfer becomes blurred and eventually reaches equipartition due to scattering.

The analysis quantifies the impact of small interface perturbations on waveguide performance.

Abstract

We present an analysis of wave propagation in a two step-index, parallel waveguide system. The goal is to quantify the effect of scattering at randomly perturbed interfaces between the guiding layers of high index of refraction and the host medium. The analysis is based on the expansion of the solution of the wave equation in a complete set of guided, radiation and evanescent modes with amplitudes that are random fields, due to scattering. We obtain a detailed characterization of these amplitudes and thus quantify the transfer of power between the two waveguides in terms of their separation distance. The results show that, no matter how small the fluctuations of the interfaces are, they have significant effect at sufficiently large distance of propagation, which manifests in two ways: The first effect is well known and consists of power leakage from the guided modes to the radiation ones. The second effect consists of blurring of the periodic transfer of power between the waveguides and the eventual equipartition of power. Its quantification is the main practical result of the paper.

Equations

$$\displaystyle\big{(}\partial_{x}^{2}+\partial_{z}^{2}\big{)}{p}(z,x)+\big{[}k{\rm n}^{(\varepsilon)}(z,x)\big{]}^{2}{p}(z,x)={f}(x)\delta^{\prime}(z),$$

$$\displaystyle{\rm n}^{(0)}(x)=\left\{\begin{array}[]{ll}\it{n}&\mbox{ if }x\in(-d/2-D,-d/2)\cup(d/2,d/2+D),\\ 1&\mbox{ otherwise,}\end{array}\right.$$

$${\rm n}^{(\varepsilon)}(z,x)=\left\{\begin{array}[]{ll}n&\mbox{ if }x\in({\cal D}_{1}^{(\varepsilon)}(z),{\cal D}_{2}^{(\varepsilon)}(z))\cup({\cal D}_{3}^{\varepsilon}(z),{\cal D}_{4}^{(\varepsilon)}(z)),\\ 1&\mbox{ otherwise,}\end{array}\right.$$

$${\mathcal{D}}_1^{(\epsilon )}(z)$$

$${\mathcal{D}}_2^{(\epsilon )}(z)$$

$${\mathcal{D}}_3^{(\epsilon )}(z)$$

$${\mathcal{D}}_4^{(\epsilon )}(z)$$

$$\mathcal{R}(z)=\mathbb{E}[{\nu }_q(0){\nu }_q(z)],q=1,\dots ,4.$$

$${p}(z,x)\in{\cal C}^{0}\big{(}(0,+\infty),H_{0}^{1}(\mathbb{R})\cap H^{2}(\mathbb{R})\big{)}\cap{\cal C}^{2}\big{(}(0,+\infty),L^{2}(\mathbb{R})\big{)},\qquad z>0,$$

$$kD\sqrt{n^2-1}<\pi .$$

$$\displaystyle{{\rm n}}_{s}(x)=\left\{\begin{array}[]{ll}\it{n}&\mbox{ if }x\in(d/2,d/2+D),\\ 1&\mbox{ otherwise,}\end{array}\right.$$

$$\frac{\sqrt{n^{2}k^{2}-\beta^{2}}}{\sqrt{\beta^{2}-k^{2}}}\tan\Bigg{(}\frac{D}{2}\sqrt{k^{2}n^{2}-\beta^{2}}\Bigg{)}=1,\qquad\beta\in(k,nk).$$

$$\displaystyle\phi(x)=\left\{\begin{array}[]{ll}\Big{(}\frac{2}{\eta}+D\Big{)}^{-1/2}\cos(\xi\frac{D}{2})\exp[\eta(x-\frac{d}{2})],&x\leq\frac{d}{2},\\ \Big{(}\frac{2}{\eta}+D\Big{)}^{-1/2}\cos\big{[}\xi(x-\frac{d}{2}-\frac{D}{2})\big{]},&\frac{d}{2}\leq x\leq\frac{d}{2}+D,\\ \Big{(}\frac{2}{\eta}+D\Big{)}^{-1/2}\cos(\xi\frac{D}{2})\exp\big{[}-\eta(x-\frac{d}{2}-D)\big{]},&x\geq\frac{d}{2}+D.\end{array}\right.$$

$$\xi =\sqrt{n^2{k}^2-{\beta }^2},\eta =\sqrt{{\beta }^2-k^2}.$$

$$\displaystyle\big{(}\partial_{x}^{2}+\partial_{z}^{2})p^{(0)}(z,x)+\big{[}k{\rm n}^{(0)}(x)\big{]}^{2}{p}^{(0)}(z,x)={f}(x)\delta^{\prime}(z),$$

$$p^{(0)}(z,x)=\sum _{t\in \{e,o\}}\frac{a_t^{(0)}}{\sqrt{{\beta }_t}}{e}^{i{\beta }_{t}z}{\varphi }_t(x)+O(z^{-2}),z>0.$$

$$\mathcal{H}(x)={\partial }_x^2+{\mathrm{n}}^{(0)}(x),$$

$${\beta }_t\in (k,nk),t\in \{e,o\}.$$

$$a_t^{(0)}=\frac{\sqrt{{\beta }_t}}{2}{\int }_{\mathbb{R}}dx\overline{{\varphi }_t(x)}f(x),t\in \{e,o\},$$

$$\exp (-\eta d)\ll 1,$$

$${\varphi }_e(x)$$

$$\displaystyle\beta_{e}=\beta+\beta^{\prime}e^{-\eta d}+o\big{(}e^{-\eta d}\big{)},\quad\beta_{o}=\beta-\beta^{\prime}e^{-\eta d}+o\big{(}e^{-\eta d}\big{)},$$

$${\beta }^{\prime }$$

$$p^{(0)}(z,x)=\frac{\phi(|x|)}{\sqrt{\beta}}e^{i\beta z}\left[1_{(0,\infty)}(x)u_{+}^{(0)}(z)+1_{(-\infty,0)}(x)u_{-}^{(0)}(z)\right]+O\big{(}e^{-\eta d}\big{)}+O(z^{-2}).$$

$$u_{+}^{(0)}(z)$$

$$u_{-}^{(0)}(z)$$

$$u_{\pm }^{(0)}(z)$$

$$z=\exp (\eta d)Z,Z>0,$$

$$a_e^{(0)}=a_o^{(0)}=\frac{a^{(0)}}{2},$$

$$\displaystyle u_{+}^{(0)}\big{(}e^{\eta d}Z\big{)}=a^{(0)}\cos(\beta^{\prime}Z),\quad\quad u_{-}^{(0)}\big{(}e^{\eta d}Z\big{)}=ia^{(0)}\sin(\beta^{\prime}Z).$$

$$1\gg \exp (-\eta d)\gg {\epsilon }^2.$$

$$\exp (-\eta d)=O({\epsilon }^2).$$

Full text

1 Introduction

Guided waves have applications in electromagnetics [9], optics and communications [22, 27], imaging underwater [19, 5, 14, 15], imaging of and in tunnels [3] and so on. The classical theory of guided waves is for ideal waveguides with perfectly reflecting straight walls and filled with homogeneous media, where the wave equation can be solved using separation of variables. The wave field is represented as a superposition of finitely many guided modes, which are waves that propagate along the axis of the waveguide, and infinitely many evanescent modes which decay away from the source. These modes do not interact with each other and thus have constant amplitude determined by the wave source [9, 27].

Motivated by applications in imaging and communications, the classical theory has been extended to waveguides filled with random media [22, 12, 13, 1], with randomly perturbed boundaries [2, 1] and with slowly changing cross-section [6, 7, 22]. Weakly guiding waveguides with confining graded-index profile affected by small random perturbations have also been analyzed in [11, 25]. The resulting mode coupling theory quantifies the interaction between the modes induced by scattering, and the consequent randomization of the wave field.

We consider waveguides with penetrable boundaries, where the guiding effect is due to a medium of high index of refraction embedded in a homogeneous background. Such waveguides are analyzed in [19, 15, 17] in the context of underwater acoustics [28] and are of great importance in optics and communications [22, 26, 27]. Motivated by the latter applications, we consider a waveguide system made of two parallel step-index waveguides, which is known as a directional coupler in integrated optics [4], [27, Chapter 10]. The classical analysis of this system, described in [18], [27, Chapter 10] and [22, Chapter 10], is based on the observation that the transverse profiles of the guided waves are essentially supported in the step-index waveguides and decay outside. When the step-index waveguides are nearby, the decaying tails penetrate the neighboring waveguide and transfer of power can occur. This is the sole coupling mechanism in ideal directional couplers and for synchronous waveguides (with identical guided mode phase velocity) there is a complete, periodic transfer of power. That is to say, if the source emits power in one step-index waveguide, this is transferred to the other waveguide and back in a periodic manner, at regular distance intervals. These intervals depend on the separation between the step-index waveguides. The larger this distance is, the weaker the coupling and the farther the waves must travel for the transfer of power to occur.

We introduce a mathematical analysis of a randomly perturbed directional coupler, where the interfaces that separate the medium with high index of refraction from the background have small amplitude random fluctuations on a scale similar to the wavelength. The classic approach in [27, Chapter 10], which is based solely on the guided modes, is inadequate in this case, because scattering at the random interfaces induces mode coupling. We take into account all the modes, the guided, radiation and evanescent ones, and quantify how their interaction affects the performance of the directional coupler. The analysis is focused on the case of well separated waveguides, where the deterministic coupling is weak. It applies to an arbitrary number of guided modes, but we describe in depth the results for the case of single guided mode step-index waveguides. We show that mode coupling induced by the random fluctuations is present no matter how far apart the waveguides are and it has two effects: The first effect is well known [22, Chapter 10] and consists of power leakage from the guided modes to the radiation modes. Our analysis captures it and shows that the leaked power is self-averaging i.e., it is independent of the realization of the random processes that model the fluctuations of the interfaces. The other effect consists of the blurring of the periodic transfer of power between the waveguides and the eventual equipartition of power between the guided waves. Its quantification in terms of the waveguide separation and amplitude of the random fluctuations is the main practical result of the paper.

The paper is organized as follows: We begin in section 2 with the mathematical formulation of the problem and then state the results in section 3. Their derivation is in sections 4–6. We end with a summary in section 7.

2 Formulation of the problem

We study the propagation of a time harmonic wave in a medium with index of refraction ${\rm n}^{(\varepsilon)}(z,x)$ defined below, that models two step-index waveguides of width $D$, separated in the transverse direction $x$ by the distance $d$, as illustrated in Figure 1. The wave field is denoted by $p(z,x)$ and it solves the two-dimensional Helmholtz equation

[TABLE]

with radiation condition at infinity, where $k$ is the wavenumber and $f(x)$ models a source supported at the origin of the range coordinate $z$.

In the case of ideal (unperturbed) waveguides, the index of refraction is range independent and equal to

[TABLE]

with $n>1$. We are interested in perturbed waveguides, where the index of refraction

[TABLE]

jumps across four randomly fluctuating interfaces

[TABLE]

The fluctuations are modeled by the zero-mean, bounded, independent and identically distributed stationary random processes $\{\nu_{q}(z)\}_{1\leq q\leq 4}$ with smooth covariance function

[TABLE]

These satisfy strong mixing conditions as defined for example in [23, section 2]. The typical amplitude of the fluctuations is much smaller than $D$ and it is modeled in (6) by the small and positive dimensionless parameter $\varepsilon$.

We study the wavefield at $z>0$, satisfying

[TABLE]

and to set radiation conditions, we suppose that the random fluctuations are supported in the range interval $(0,L^{(\varepsilon)})$. We will see that net scattering effect of these fluctuations becomes of order one at range scales of order $\varepsilon^{-2}$, so we let $L^{(\varepsilon)}=L/\varepsilon^{2}$. We will also see that for the assumed smooth covariance ${\cal R}(z)$, the guided waves propagate mostly in the forward direction and do not interact with any fluctuations at $z<0$, which is why we neglect them.

The goal of the paper is to quantify how scattering at the random interfaces (6) affects the coupling of the two step-index waveguides centered at $x=\pm(d+D)/2$. We consider in particular the case of a sufficiently large separation distance $d$ between the waveguides, where the deterministic coupling is very weak but the coupling induced by the random fluctuations is still present.

3 Statement of results

We state here the main results of the paper, derived in sections 4–6 by decomposing the wavefield into guided, radiation and evanescent modes of the waveguide system made of two parallel step-index waveguides illustrated in Figure 1. While our analysis applies to an arbitrary number of guided modes, in this section we consider only the case where an isolated step-index waveguide has only one guided mode. This case captures all the essential aspects of the problem and arises when

[TABLE]

The two ideal step-index waveguides centered at $\pm(d/2+D/2)$ do not interact in the waveguide system when $d\to\infty$. Thus, we can write the wavefield for large $d$ in terms of the unique guided mode $\phi(x)e^{i\beta z}$ of the step-index waveguide centered at $x=d/2+D/2$, modeled by the index of refraction

[TABLE]

and in terms of the unique guided mode $\phi(-x)e^{i\beta z}$ of the step-index waveguide centered at $x=-d/2-D/2$, modeled by the index of refraction ${{\rm n}}_{s}(-x)$. Here $\phi(x)$ is the eigenfunction of the Helmholtz operator $\partial_{x}^{2}+{\rm n}_{s}(x)$ for the eigenvalue $\beta^{2}$. This $\beta$ is defined in Lemma 1 as the unique solution in the interval $(k,nk)$ of

[TABLE]

The expression of the eigenfunction $\phi(x)$ is [20, Section 2], [22, Chapter 1]

[TABLE]

It has a peak centered at $x=d/2+D/2$ of width $1/\xi$, and an exponentially decaying tail, at the rate $\eta$, where

[TABLE]

3.1 Coupling of ideal waveguides

We show in section 4 that under the assumption (9) and for sufficiently large separation distance $d$, the solution of

[TABLE]

with radiation condition at infinity, takes the form***We use consistently the index ${(0)}$ for the mode amplitudes and wave field in the ideal (unperturbed) two step-index wave guide system.

[TABLE]

The terms in the sum model the guided modes of the waveguide system, whereas the $O(z^{-2})$ remainder accounts for the radiation and evanescent modes [20, Section 3]. The guided modes $\{\phi_{t}(x)e^{i\beta_{t}z}\}_{t=e,o}$ are waves that propagate along the range coordinate $z$ and have transverse profiles essentially supported in the two step-index waveguides. They are defined by the even and odd eigenfunctions $\phi_{e}(x)$ and $\phi_{o}(x)$ of the operator

[TABLE]

with ${\rm n}^{(0)}(x)$ given in (4), for the eigenvalues $\beta_{e}^{2}$ and $\beta_{o}^{2}$ satisfying

[TABLE]

These eigenfunctions are given explicitly in section 4.1. They have peaks close to the center axes $x=\pm(d/2+D/2)$ of the step-index waveguides and decay outside their transverse support, as illustrated in Figure 2. Each guided mode is multiplied in (20) by the constants

[TABLE]

called the mode amplitudes, that are determined by the source. Throughout the paper the bar denotes complex conjugate.

We assume that the waveguides separation distance $d$ is sufficiently large, so that

[TABLE]

with $\eta$ defined in (18), and obtain that the eigenfunctions can be approximated by

[TABLE]

where “sgn” is the sign function. The accuracy of this approximation is illustrated in the right plot of Figure 2. Consequently, the transverse profile of the even guided mode presents two positive peaks centered at $x=\pm(d/2+D/2)$, with exponentially decaying tails, while the transverse profile of the odd guided mode presents one positive peak centered at $x=d/2+D/2$, one negative peak centered at $x=-d/2-D/2$, and exponentially decaying tails. The form of these peaks is proportional to the unique eigenfunction (17) of the single-mode step-index waveguide.

We also obtain that the wave numbers are

[TABLE]

with

[TABLE]

so the wavefield (20) in the ideal waveguide system takes the form

[TABLE]

Here we introduced the range-dependent amplitudes of the waves propagating in the two step-index waveguides

[TABLE]

with indexes “$\pm$” corresponding to the waveguide centered at $x=\pm(d+D)/2$.

Equation (28) shows that the wavefield consists of two components, the first one is centered at $x=d/2+D/2$ with the form $\phi(|x|)e^{i\beta z}{\bf 1}_{(0,\infty)}(x)$ and the second one is centered at $x=-d/2-D/2$ with the form $\phi(|x|)e^{i\beta z}{\bf 1}_{(-\infty,0)}(x)$. This is similar to the case of two independent, single-mode step-index waveguides, except that in (28) the amplitudes $u_{\pm}(z)$ vary in $z$, due to coupling. We obtain from (24) and (29–30) that for $z$ of the order of the wavelength,

[TABLE]

However, at large $z$, satisfying

[TABLE]

$u_{\pm}^{(0)}(z)$ oscillate periodically in $z$. For example, if the source gives the amplitudes

[TABLE]

according to (23), so that $u_{+}^{(0)}(0)=a^{(0)}$ and $u_{-}^{(0)}(0)=0$, at range (32) we have

[TABLE]

In conclusion, the total wave power in the ideal waveguide system at large range $z$ is essentially supported in the two step-index waveguides. The wave power in the waveguide centered at $x=\pm(d/2+D/2)$ is proportional to $|u_{\pm}^{(0)}(z)|^{2}$, and equation (34) shows that it oscillates slowly and periodically. At scaled distance $Z=m\pi/\beta^{\prime}$, $m\in\mathbb{N}$, the wave power is concentrated in the waveguide centered at $x=d/2+D/2$, whereas at $Z=(1/2+m)\pi/\beta^{\prime}$, $m\in\mathbb{N}$, the wave power is concentrated in the waveguide centered at $x=-d/2-D/2$. These periodic oscillations have been reported in the literature [27, Chapter 10]. The standard method to analyze them is not to start from the analysis of the modes of the waveguide system, as we do in section 4, but to simplify by assuming that the modes can be represented as a weighted sum of the guided modes of the two waveguides. This simplified approach does not allow to take into account the role of evanescent and radiation modes, which are critical to the study of random waveguides in sections 5–6, with results described next.

3.2 Coupling of random waveguides

The analysis of the solution $p(z,x)$ of the Helmholtz equation (1) with index of refraction (5) is carried out in sections 5–6. It shows that under the assumptions (9), (24) and at large range $z/\varepsilon^{2}$, where scattering at the random interfaces (6) becomes significant, there are three distinguished regimes that determine the coupling between the waveguides:

The “moderate coupling” regime, where the separation distance $d$ is moderately large, satisfying

[TABLE] 2. 2.

The “weak coupling” regime, where $d$ is large enough so that

[TABLE] 3. 3.

The “very weak coupling” regime, where $d$ is so large that

[TABLE]

We now describe the results in each of these three regimes.

3.2.1 Moderate coupling

At large range $z/\varepsilon^{2}$ and in the regime defined by (35), the solution of (1) with radiation condition at infinity and with index of refraction (5) is

[TABLE]

where $u_{\pm}(z)$ are random processes and $o(1)$ denotes a residual that tends to zero as $\varepsilon\to 0$. This residual accounts for the radiation and evanescent components of $p$.

Similar to (28), we have a propagating wave $\phi(|x|)e^{i\beta z/\varepsilon^{2}}{\bf 1}_{(0,\infty)}(x)$ centered at $x=d/2+D/2$ and another wave $\phi(|x|)e^{i\beta z/\varepsilon^{2}}{\bf 1}_{(-\infty,0)}(x)$ centered at $x=-d/2-D/2$. The coupling between the two waveguides is described by the random variations of the complex amplitudes $u_{\pm}(z)$, the analogues of (29–30). To write their expressions, we introduce the notation

[TABLE]

which takes into account that the residual in (26) is not negligible under the assumption (35) at range $O(\varepsilon^{-2})$. We have

[TABLE]

where \big{(}a_{e}(z),a_{o}(z)\big{)} is a random, Markovian process defined at $z\geq 0$, with initial condition $a_{t}(0)=a_{t}^{(0)}$ given in (23) for $t\in\{e,o\}$, and with the infinitesimal generator in Theorem 2. To describe the results, it suffices to state from there that

[TABLE]

with probability one, and that

[TABLE]

with positive $\Lambda$ and $\Gamma$ defined by

[TABLE]

in terms of

[TABLE]

and the power spectral density $\widehat{\cal R}\geq 0$, the Fourier transform of the covariance (7).

We have therefore from (40) and (41) that the total power of the guided waves decays exponentially, at the rate $\Lambda$,

[TABLE]

This decay models the transfer of power from the guided modes to the radiation modes, induced by scattering at the random interfaces (6).

The imbalance of power between the two waveguides is quantified by

[TABLE]

and its expectation is

[TABLE]

To explain what this gives, consider the source excitation (33) with $a_{e}^{(0)}=a_{o}^{(0)}$, so that the initial wavefield is supported in the waveguide centered at $x=d/2+D/2$. Then, equation (48) becomes

[TABLE]

and it describes the competition between the deterministic and random coupling of the waveguides. The cosine in (49) models the deterministic coupling which induces periodic oscillations of the power, as in section 3.1. The random coupling is modeled by the exponential decay in $z$, at the rate $\Gamma$. It shows that as the range increases, the power tends to become equally distributed among the two waveguides. This decay is present in (48) as well, so the equipartition of power at large $z$ is independent of the initial condition generated by the source.

3.2.2 Weak coupling

When the separation distance $d$ between the waveguides satisfies (36), the wavefield $p(z/\varepsilon^{2},x)$ has the same expression as in (38), (40), but the random processes $(u_{+}(z),u_{-}(z))$ have different statistics.

The total power of the guided waves is still given by (46), and decays at the same rate $\Lambda$ defined in (43). However, the expectation of the imbalance of power between the two waveguides satisfies the damped harmonic oscillator equation

[TABLE]

with $\beta^{\prime}$ defined in (27) and

[TABLE]

Based on the value of $\theta$, which is independent of $\varepsilon$ by assumption (36), we distinguish three regimes, which we describe for the source excitation (33), with $u_{-}^{(\varepsilon)}(0)=0$ and therefore $\mathbb{E}[\mathcal{P}(0)]=1$:

1. When $2\theta\beta^{\prime}<\Gamma$, the solution of (50) is

[TABLE]

This tends to $1$ as $\theta\to 0$, meaning that when the two waveguides are very far apart, there is no transfer of power between them. This is just as in the ideal (deterministic) waveguide system. However, unlike in ideal waveguides, the power is transferred from the guided modes to the radiation ones, as described by the exponential decay in (46).

For a finite $\theta$, we have $\mathbb{E}[\mathcal{P}(z)]\to 0$ as $z\to\infty$, so the random coupling distributes the power evenly among the two waveguides.

2. In the critical case $2\theta\beta^{\prime}=\Gamma$ we have

[TABLE]

As in the previous case, the random coupling equidistributes the power among the two waveguides, in the limit $z\to\infty$.

3. When $2\theta\beta^{\prime}>\Gamma$, the solution of (50) is

[TABLE]

It displays periodic oscillations induced by the deterministic coupling of the waveguides, but these oscillations are damped due to the random coupling. In particular, if $2\theta\beta^{\prime}\gg\Gamma$, we get

[TABLE]

in agreement with (49).

We plot in Figure 3 the imbalance ratio $\mathbb{E}[\mathcal{P}(z)]$ as a function of $z$ normalized by the deterministic coupling distance†††At the distance $\pi z_{\theta}$ the power is fully transferred from one step-index waveguide to the other one, when there are no random perturbations. $z_{\theta}=1/(2\theta\beta^{\prime})$. From this plot and from equation (46) we conclude that wave scattering at the random interfaces (6) has two net effects:

It induces a self-averaging loss (or leakage) of total power, due to the coupling of the guided modes with the radiation modes. 2. 2.

It causes a blurring of the periodic (deterministic) power transfer from one waveguide to the other.

The blurring effect at item 2 is the main practical result of the paper. It depends on the effective parameter $\Gamma$ and it is significant as soon as $\Gamma z_{\theta}$ becomes of order one. Therefore, the deterministic transfer of power is very sensitive to the random fluctuations of the interfaces (6).

3.2.3 Very weak coupling

When the distance $d$ between the waveguides is so large that (37) holds, the wavefield has the same expression as (38), but the wave amplitudes have different statistics. They model the only coupling in this regime, between the guided and radiation modes, which generates effective wave power leakage. Explicitly, we show in section 6 that the wave amplitudes converge in probability, as $\varepsilon\to 0$, to the deterministic function $(|u_{+}(z)|^{2},|u_{-}(z)|^{2})$ satisfying

[TABLE]

Although deterministic, this function is not as in the ideal waveguide system, where $(|u_{+}^{(0)}(z)|^{2},|u_{-}^{(0)}(z)|^{2})$ is constant in $z$ (because $\beta^{\prime}e^{-\eta d}z\ll 1$ in (29-30) in this regime). Instead, it decays exponentially at rate $\Lambda$, due to the power leakage. The imbalance of power between the waveguides is constant

[TABLE]

so in the case of the source excitation (33), ${\cal P}(z)\simeq{\cal P}(0)=1.$

4 Analysis in ideal waveguides

The analysis in this section is classical and follows the lines of [20]. It derives the results stated in section 3.1 by expanding the wave field on a complete set of eigenmodes. The proof of the completeness of this set is the most delicate part and it is carried out in [20] by the Levitan-Levinson method [8, Chapter 9].

Recall the Helmholtz operator (21) in the transverse coordinate, with index of refraction ${\rm n}^{(0)}(x)$ given in (4), and note that it is self-adjoint with respect to the scalar product associated to the $L^{2}$-norm,

[TABLE]

Its spectrum is (-\infty,k^{2})\cup\big{\{}\big{(}\beta_{t,j}^{2}\big{)}_{1\leq j\leq N_{t}},~{}~{}t\in\{e,o\}\big{\}}, where $\beta_{t,j}$ are called the guided mode wavenumbers. They are positive and satisfy the order relation

[TABLE]

where the index $t\in\{e,o\}$ stands for the even and odd eigenfunctions in the transverse coordinate $x$.

We describe next the eigenfunctions for the discrete spectrum and the improper eigenfunctions for the continuum spectrum, and explain that they form a complete set. We use them to decompose the wavefield into guided, radiation and evanescent modes with amplitudes determined by the source.

4.1 Discrete spectrum

There are two sets of discrete eigenvalues and eigenfunctions: the first one associated with the even modes in $x$ and the second one associated with the odd modes.

The $j$-th even eigenfunction $\phi_{e,j}(x)$ for the eigenvalue $\beta_{e,j}^{2}$ is defined by

[TABLE]

and $\phi_{e,j}(-x)=\phi_{e,j}(x)$ for $x\geq 0$. Here

[TABLE]

and $A_{e,j}>0$ is the normalization constant

[TABLE]

calculated so that $\phi_{e,j}$ has unit $L^{2}$-norm. Moreover, $\beta_{e,j}\in(k,nk)$ satisfies the dispersion relation

[TABLE]

with $\xi_{e,j},\eta_{e,j}$ given by (62), which ensures the continuity of (61) and its derivative. The number of solutions $\beta_{e,j}\in(k,nk)$ of (64) is denoted by $N_{e}$ and when $d$ is large, it depends on the value of $kD\sqrt{n^{2}-1}$.

The $j$-th odd eigenfunction for the eigenvalue $\beta_{o,j}^{2}$ is defined similarly,

[TABLE]

with $\phi_{o,j}(-x)=-\phi_{o,j}(x)$ for $x\geq 0$ and

[TABLE]

The normalization constant $A_{o,j}>0$ is given by

[TABLE]

so that $\phi_{o,j}$ has unit $L^{2}$ norm and $\beta_{o,j}\in(k,nk)$ satisfies the dispersion relation

[TABLE]

with $\xi_{o,j},\eta_{o,j}$ given by (69), for $j=1,\ldots,N_{o}$, which ensures that (68) and its derivative are continuous.

4.2 Continuous spectrum

For $\gamma\in(-\infty,k^{2})$ there are two improper eigenfunctions, even and odd, denoted by $\phi_{t,\gamma}(x)$, for $t\in\{e,o\}$. We write their expression below in terms of the parameters

[TABLE]

The even eigenfunctions satisfy $\phi_{e,\gamma}(-x)=\phi_{e,\gamma}(x)$ for $x\geq 0$, and are defined by

[TABLE]

and by

[TABLE]

and by

[TABLE]

with normalization constant $A_{e,\gamma}>0$ given by

[TABLE]

The odd eigenfunctions satisfy $\phi_{o,\gamma}(-x)=-\phi_{o,\gamma}(x)$ for $x\geq 0$, and are defined by

[TABLE]

and by

[TABLE]

and by

[TABLE]

with normalization constant $A_{o,\gamma}>0$ given by

[TABLE]

Remark 4.1**.**

Note that the improper eigenfunctions $\phi_{t,\gamma}(x)$, for $t\in\{e,o\}$, are not in $L^{2}(\mathbb{R})$. However, we can define for any $\varphi\in L^{2}(\mathbb{R})$

[TABLE]

where the limit holds in $L^{2}(-\infty,k^{2})$. The normalizing constants $A_{t,\gamma}$ given in (76) and (80) are such that, for any test function $a\in L^{2}(-\infty,k^{2})$,

[TABLE]

They depend continuously on $\gamma$ and are bounded on $(-\infty,k^{2})$ with A_{t,\gamma}\simeq\big{(}2\pi\sqrt{|\gamma|}\big{)}^{-1/2} as $\gamma\to-\infty$.

4.3 Completness

The proof of completness is based on the Levitan-Levinson method [8, Chapter 9] and is the same as that in [20, Section 2.3]. We state directly the result:

For any $\varphi\in L^{2}(\mathbb{R})$ we have the Parseval identity

[TABLE]

and the map which assigns to $\varphi$ the coefficients of its spectral decomposition

[TABLE]

is an isometry from $L^{2}(\mathbb{R})$ onto $\mathbb{C}^{N_{e}+N_{o}}\times L^{2}(-\infty,k^{2})^{2}$. Moreover, there exists a resolution of the identity $\Pi$ of the operator (21) such that: for any $\varphi\in L^{2}(\mathbb{R})$ and any $-\infty\leq r\leq r^{\prime}\leq+\infty$,

[TABLE]

and for any $\varphi$ in the domain of ${\cal H}$,

[TABLE]

4.4 Modal decomposition

Using the results in the previous section, we can write the solution $p^{(0)}(z,x)$ of the Helmholtz equation (19) as a superposition of modes

[TABLE]

The first sum contains the guided modes, where $p_{t,j}^{(0)}(z)$ are one-dimensional waves that propagate along $z$, with wavenumber $\beta_{t,j}$,

[TABLE]

The other two sums contain the radiation modes and evanescent modes, where

[TABLE]

The radiation modes correspond to $\gamma\in(0,k^{2})$, so that ${p}_{t,\gamma}^{(0)}(z)$ are one-dimensional waves that propagate along $z$ at wavenumber $\sqrt{\gamma}$. The evanescent modes are decaying exponentially‡‡‡ We neglect the case $\gamma=0$ because there is no coupling of the modes in ideal waveguides, and therefore it has a negligible contribution in (81). Mode coupling occurs in randomly perturbed waveguides, so the case $\gamma=0$ must be dealt with carefully, as explained in section 5.6. in $z$ on the range scale $1/\sqrt{|\gamma|}$, for $\gamma<0$.

The solutions of (82–83) that satisfy the radiation conditions are

[TABLE]

for $t\in\{e,o\}$, with constant mode amplitudes determined by the source in (1),

[TABLE]

Substituting in (81), we obtain the modal expansion of the wavefield $p^{(0)}(z,x)$.

4.5 Single guided component waveguides

We describe here the special situation considered in section 3.1, where the assumptions (9) and (24) hold and there is a single guided mode per step-index waveguide, as we now explain.

Lemma 1**.**

Under the assumption $kD\sqrt{n^{2}-1}<\pi$, equation

[TABLE]

has a unique solution $\beta$ in $(k,nk)$, whereas equation

[TABLE]

has no solution in $(k,nk)$. Furthermore, when the distance $d$ between the waveguides is large enough, so that $\exp(-\eta d)\ll 1$, with $\eta=\sqrt{\beta^{2}-k^{2}}$, equations (64) and (71) have a single solution, meaning that $N_{e}=N_{o}=1$.

The proof of the lemma is in appendix A. We let henceforth $\beta$ be the unique solution of (87) in the interval $(k,nk)$ and recall the definition (18) of $\xi$ and $\eta$ in terms of $\beta$. Since $N_{e}=N_{o}=1$, we simplify notation as $\beta_{t,1}\leadsto\beta_{t},$ for $t\in\{e,o\}$, and obtain from (64) and (71) that

[TABLE]

with

[TABLE]

Similarly, we let $\phi_{t,1}(x)\leadsto\phi_{t}(x)$, for $t\in\{e,o\}$, and obtain from (61)–(70) that

[TABLE]

Substituting in (81) and obtaining from [20, Section 3.2] that the radiation and evanescent modes have an $O(z^{-2})$ contribution§§§Note that in [10] the decay is $O(z^{-1})$ because the source is $f(x)\delta(z)$. In our case the decay is $O(z^{-2})$, because we have the $z$ derivative of the solution in [10], due to the source $f(x)\delta^{\prime}(z)$., we obtain equations (20) and (28–30) used in section 3.1 to describe the deterministic coupling between the ideal step-index waveguides.

5 Analysis in randomly perturbed waveguides

Let us introduce the notation

[TABLE]

and rewrite the index of refraction defined in (5–6) as

[TABLE]

with

[TABLE]

To explain this reformulation, we illustrate in Figure 4 the fluctuations $\varepsilon D\nu_{1}(z)$ of the interface at $x=-d/2-D$. At range $z_{1}$ we have $\nu_{1}(z_{1})<0$, so $(z_{1},x)$ with $x\in(-d/2-D+\varepsilon D\nu_{1}(z_{1}),-d/2-D)$ is in the waveguide and therefore,

[TABLE]

At range $z_{2}$ we have $\nu_{1}(z_{2})>0$, so $(z_{2},x)$ with $x\in(-d/2-D,-d/2-D+\varepsilon D\nu_{1}(z))$ is outside the waveguide and

[TABLE]

The same reasoning applies to all the interfaces.

Our goal in this section is to analyze the solution $p(z,x)$ of (1), with index of refraction (92) and radiation condition at infinity, and in particular, to derive the results stated in section 3.2.

5.1 Modal decomposition

The completness result of section 4.3 allows the expansion of $p(z,x)$, $z$ by $z$, in terms of the eigenfunctions defined in sections 4.1–4.2,

[TABLE]

Here $p_{t,j}(z)$ are complex-valued amplitudes of guided modes with wavenumber $\beta_{t,j}$ satisfying (64) or (71). They propagate along $z$ and satisfy the following one-dimensional Helmholtz equations at $z>0$,

[TABLE]

Similarly, $p_{t,\gamma}(z)$ are complex-valued amplitudes of modes that are propagating for $\gamma\in(0,k^{2})$ (radiation modes) and decaying for $\gamma<0$ (evanescent modes), and satisfy the one-dimensional Helmholtz equations

[TABLE]

The source terms in these equations are due to the random fluctuations, which are supported at $z\in(0,L/\varepsilon^{2})$. These induce mode coupling, modeled by the random coefficients

[TABLE]

Recalling the definition (92) of $V^{(\varepsilon)}(z,x)$ and using Taylor expansions of the eigenfunctions $\phi_{t,j}(x)$ and $\phi_{t,\gamma}(x)$ around $x=\pm d/2$ and $x=\pm(d/2+D)$, we obtain a power series (in $\varepsilon$) expression of these coefficients

[TABLE]

and similarly for $C_{t,j,t^{\prime},\gamma^{\prime}}^{(\varepsilon)},C_{t,\gamma,t^{\prime},l^{\prime}}^{(\varepsilon)},C_{t,\gamma,t^{\prime},\gamma^{\prime}}^{(\varepsilon)}$.

Since we consider propagation distances of the order of $\varepsilon^{-2}$, we neglect from now on the $o(\varepsilon^{2})$ terms in these expansions, because they give no contribution in the limit $\varepsilon\to 0$. To simplify the presentation, we will also not write the $\varepsilon^{2}$ terms, even though they are not negligible. They contribute to the expression of the infinitesimal generators of the limit processes via terms that are proportional to $\frac{\Delta k^{2}}{2\beta_{t,j}}\mathbb{E}[c_{t,j,t,j}(z)]$. Because $\partial_{x}(\phi_{t,j}^{2})$ is odd, for $t\in\{e,o\}$, we obtain from (102) that

[TABLE]

Moreover, in the case of two weakly coupled single guided component waveguides,

[TABLE]

We will incorporate these terms in the statements of the results without giving additional details of their derivation. They correspond to effective deterministic phases of the mode amplitudes, in the limit $\varepsilon\to 0$, as explained in Remark 5.2.

5.2 Forward and backward going waves

The guided modes can be decomposed further in forward and backward guided modes. This decomposition is basically the method of variation of parameters for the perturbed Helmholtz equations (94–95), where we define the complex valued amplitudes

[TABLE]

for $t\in\{e,o\}$, such that

[TABLE]

and

[TABLE]

Substituting (105–106) in (94–95), we obtain that the amplitudes (104) satisfy the following first-order system of stochastic differential equations

[TABLE]

[TABLE]

[TABLE]

[TABLE]

The right-hand sides in these equations are supported at $z\in(0,L/\varepsilon^{2})$ and model the mode coupling induced by scattering at the random interfaces (6). This coupling causes the randomization of the amplitudes (104) and therefore of the wavefield.

5.3 Role of the evanescent modes

Equations (107–110) show that the forward and backward going amplitudes of both the guided and radiation modes are coupled to each other and to the evanescent modes $p_{t,\gamma}(z)$, for $\gamma\in(-\infty,0)$ and $t\in\{e,o\}$. These mode amplitudes satisfy

[TABLE]

where $g_{t,\gamma}(z)$ and $g^{\rm ev}_{t,\gamma}(z)$ are supported at $z\in(0,L/\varepsilon^{2})$ and are defined by

[TABLE]

We now explain that the solution of (111) can be expressed in terms of the amplitudes (104) when $\varepsilon$ is small enough. This result has already been obtained in similar frameworks [15] and it allows us to obtain a closed system of equations for the guided and radiation mode amplitudes.

Let us invert the operator $\partial_{z}^{2}+\gamma$ in (111) using the Green’s function that satisfies the radiation condition (i.e., it decays away from the source). This gives

[TABLE]

at $z>0$, for $\gamma<0$ and $t\in\{e,o\}$. The first term in this expression comes from the random mode coupling induced by scattering at the random interfaces (6). The second term comes from the source. Since we will consider propagation distances $z$ of the order of $\varepsilon^{-2}$, we can neglect this second term.

Using (113) and the notation ${\itbf p}^{\rm ev}_{t}=({p}_{t,\gamma}(z))_{z>0,\gamma<0}$ for $t\in\{e,o\}$, we obtain that

[TABLE]

where $I$ is the identity operator, ${\cal J}_{t}$ is the integral operator

[TABLE]

coming from the $g_{t,\gamma}^{\rm ev}$ term in (114), and the right-hand side is

[TABLE]

for $z>0,\gamma<0$. The operator ${\cal J}_{t}$ is bounded from ${\cal C}^{0}((0,\infty),L^{1}(-\infty,0))$ to itself and as in [15], it can be shown that

[TABLE]

and therefore

[TABLE]

Furthermore, equations (107–110) give that ${a}_{t^{\prime},l^{\prime}}(z^{\prime})$, ${b}_{t^{\prime},l^{\prime}}(z^{\prime})$, ${a}_{t^{\prime},\gamma^{\prime}}(z^{\prime})$, ${b}_{t^{\prime},\gamma^{\prime}}(z^{\prime})$ are equal to ${a}_{t^{\prime},l^{\prime}}(z)$, ${b}_{t^{\prime},l^{\prime}}(z)$, ${a}_{t^{\prime},\gamma^{\prime}}(z)$, ${b}_{t^{\prime},\gamma^{\prime}}(z)$ up to terms of order $\varepsilon$, as long as $|z-z^{\prime}|=O(1)$, where the exponential in (117) contributes.

Gathering the results we obtain

[TABLE]

for $z>0,\gamma<0$ and $t\in\{e,o\}$. Substituting this expression into (107-110) we get a closed system for the guided and radiation mode amplitudes (104).

5.4 Long range scaling

The resulting stochastic system for the guided and radiation mode amplitudes has a right-hand side which is proportional to $\varepsilon$. Therefore, the amplitudes are not affected by scattering at the random interfaces (6) until the waves travel at long enough range, dependent of $\varepsilon$. Note that the $O(\varepsilon)$ coupling terms in (107–110) have zero expectation, so we need a central limit theorem type scaling, with range $z\leadsto{z}/{\varepsilon^{2}},$ in order to observe net scattering effects.

We rename the mode amplitudes in this long range scaling as

[TABLE]

for $t\in\{e,o\}$. Recalling that the random fluctuations are supported in the range interval $(0,L/\varepsilon^{2})$ and using the radiation conditions (i.e., the waves are outgoing at $z>L$), we obtain

[TABLE]

where ${a}_{t,j}^{(0)}$ and ${a}_{t,\gamma}^{(0)}$ are given in (84–85), for $t\in\{e,o\}$ . These boundary conditions and the closed system of stochastic differential equations described at the end of the previous section define the random amplitudes (120).

5.5 Forward scattering approximation

We now introduce the forward scattering approximation, where the coupling between forward and backward going modes can be neglected. It follows from the following facts, under the assumption that the power spectral density $\widehat{\cal R}(\kappa)$ has compact support or fast decay so that

[TABLE]

1. In the limit $\varepsilon\to 0$, the coupling between the forward and backward guided modes depends on the coefficients

[TABLE]

while the coupling among the forward-guided modes depends on

[TABLE]

Recalling definitions (92) and (96–99) and that $\beta_{t,j}\in(k,nk)$, we conclude that the coefficients (124), which are proportional to $\widehat{\cal R}(\beta_{t^{\prime},l^{\prime}}+\beta_{t,j})$, are negligible under the assumption (123). Therefore, we can neglect the coupling between the forward and backward guided modes. Nevertheless, the forward guided modes are coupled among themselves, because the coefficients (125), which are proportional to $\widehat{\cal R}(\beta_{t^{\prime},l^{\prime}}-\beta_{t,j})$, are not negligible.

2. The coupling between forward guided and backward radiation modes depends in the limit $\varepsilon\to 0$ on the coefficients

[TABLE]

for $\gamma^{\prime}\in(0,k^{2})$. These are proportional to $\widehat{\cal R}(\beta_{t,j}+\sqrt{\gamma^{\prime}})$ and since $\beta_{t,j}\in(k,nk)$, this type of coupling is negligible under the assumption (123).

3. The coupling between any radiation modes is negligible in our regime $\varepsilon\to 0$, as we will see in the following. This is because these modes are essentially spatially supported outside the region where the medium is fluctuating.

We suppose henceforth that assumption (123) holds. Because there is no coupling between the forward and backward going modes, and the backward mode amplitudes satisfy the homogeneous boundary conditions (121–122) at $z=L$, we can set

[TABLE]

The expression (93) of the wavefield simplifies to

[TABLE]

and the forward guided mode amplitudes ${a}_{t,j}^{(\varepsilon)}$ and ${a}_{t,\gamma}^{(\varepsilon)}$ satisfy

[TABLE]

and

[TABLE]

where

[TABLE]

5.6 Diffusion limit theorem

The next theorem gives the $\varepsilon\to 0$ limit of the random forward-going mode amplitudes, the solutions of (128–129) with the initial conditions given in (121–122).

Theorem 2**.**

Suppose that the wavenumbers $\beta_{t,j}$ are distinct, for $1\leq j\leq N_{t}$ and $t\in\{e,o\}$. Then, the process \big{(}({a}_{t,j}^{(\varepsilon)}(z))_{j=1}^{N_{t}},({a}_{t,\gamma}^{(\varepsilon)}(z))_{\gamma\in(0,k^{2})},t\in\{e,o\}\big{)} converges in distribution in {\cal C}^{0}\big{(}[0,{L}],\mathbb{C}^{N}\times L^{2}(0,k^{2})^{2}\big{)}, where $\mathbb{C}^{N}\times L^{2}(0,k^{2})^{2}$ is equipped with the weak topology, to the Markov process \big{(}({a}_{t,j}(z))_{j=1}^{N_{t}},({a}_{t,\gamma}(z))_{\gamma\in(0,k^{2})},t\in\{e,o\}\big{)}. The infinitesimal generator of this limit process is ${\cal L}={\cal L}^{1}+{\cal L}^{2}+{\cal L}^{3},$ where $\{{\cal L}^{j}\}_{1\leq j\leq 3}$ are the differential operators:

[TABLE]

In these definitions we use the classical complex derivative:

[TABLE]

and the coefficients of the operators (131)–(133) are defined for indexes $j=1,\ldots,N_{t},$ $l^{\prime}=1,\ldots,N_{t^{\prime}}$ and $t,t^{\prime}\in\{e,o\}$, as follows:

For all $(t,j)\neq(t^{\prime},l^{\prime})$:

[TABLE]

For all $(t,j),(t^{\prime},l^{\prime})$:

[TABLE]

For all $(t,j)$:

[TABLE]

The proof of this diffusion limit theorem is based on a martingale approach using the perturbed test function method. It is an extension of the diffusion approximation theorem in [23, 10] that is carried out in [15] and [16, Sections 3, 4.1]. Here we reproduce the result but we simplify its statement. The rigorous statement [16, Theorem 4.1] has two steps: The first step deals with the convergence of truncated processes, obtained using the resolution of the identity $\Pi_{\mathbb{R}\setminus(-\gamma_{\star},\gamma_{\star})}$ defined in section (4.3). This removes the vicinity of $\gamma=0$, in order to take the diffusion limit $\varepsilon\to 0$. The second step deals with the convergence of the truncated processes themselves, in the limit $\gamma_{\star}\to 0$.

Remark 5.1**.**

Note that:

1. The convergence result in Theorem 2 holds in the weak topology.

*2. The infinitesimal generator ${\cal L}$ does not involve the derivatives $\partial_{a_{t,\gamma}}$. Therefore \Big{(}({a}_{t,j}^{(\varepsilon)}(z))_{j=1}^{N_{t}},t\in\{e,o\}\Big{)} converges in distribution in ${\cal C}^{0}([0,{L}],\mathbb{C}^{N})$ to the Markov process \Big{(}({a}_{t,j}(z))_{j=1}^{N_{t}},t\in\{e,o\}\Big{)} with generator ${\cal L}$. Of course, the weak and strong topologies are the same in $\mathbb{C}^{N}$.

3. If the generator ${\cal L}$ is applied to a test function that depends only on the mode powers \Big{(}(P_{t,j}=|{a}_{t,j}|^{2})_{j=1}^{N_{t}},t\in\{e,o\}\Big{)}, then the result is a function that depends only on \Big{(}(P_{t,j})_{j=1}^{N_{t}},t\in\{e,o\}\Big{)}. Thus, the mode powers \Big{(}(P_{t,j}(z))_{j=1}^{N_{t}},t\in\{e,o\}\Big{)} define a Markov process, with infinitesimal generator

[TABLE]

4. The radiation mode amplitudes remain constant on $L^{2}(0,k^{2})^{2}$, equipped with the weak topology, as $\varepsilon\to 0$. However, this does not describe the power transported by the radiation modes $\displaystyle\sum_{t\in\{e,o\}}\int_{0}^{k^{2}}|a_{t,\gamma}|^{2}d\gamma$, because the convergence does not hold in the strong topology of $L^{2}(0,k^{2})^{2}$.

The third point in Remark 5.1 implies that the mean mode powers satisfy the closed system of equations

[TABLE]

with initial conditions

[TABLE]

These equations (135) can be found in the literature [22, 24]. Here we derived them from first principles, taking into account the coupling between all types of modes.

Remark 5.2**.**

The three terms $\{{\cal L}^{j}\}_{1\leq j\leq 3}$ of the the limit infinitesimal generator ${\cal L}$ account for the mode coupling as follows:

1. The operator ${\cal L}^{1}$ accounts for coupling of the guided modes, modeled by the coefficients $\{\Gamma_{t,j,t^{\prime},j^{\prime}}^{c}\}_{1\leq j\leq N_{t},1\leq j^{\prime}\leq N_{t^{\prime}},t,t^{\prime}\in\{e,o\}}.$ This coupling results in the power exchange between these modes.

2. The operator ${\cal L}^{2}$ is due to coupling between the guided and radiation modes, which causes power leakage from the guided modes to the radiation ones (effective diffusion), modeled by the coefficients $\{\Lambda_{t,j}^{c}\}_{1\leq j\leq N_{t},t\in\{e,o\}}$.

3. The operator ${\cal L}^{3}$ accounts for coupling between the guided and evanescent modes, which gives the term proportional to $\{\kappa_{t,j}^{\rm ev}\}_{1\leq j\leq N_{t},t\in\{e,o\}}.$ The other term, proportional to $\{\kappa_{t,j}\}_{1\leq j\leq N_{t},t\in\{e,o\}}$, is due to the neglected terms discussed at the end of section 5.1. Both contributions give only additional phase terms on the guide modes i.e., an effective dispersion.

6 Two weakly coupled single mode random waveguides

Using the results of the previous section and the decay in range (like $({z}/{\varepsilon^{2}})^{-2}=\varepsilon^{4}/z^{2}$) of the contribution of the radiation and evanescent modes, we obtain that

[TABLE]

where the amplitudes $\{a_{t,j}(z)\}_{1\leq j\leq N_{t},t\in\{e,o\}}$ are described in Theorem 2 and the residual $o(1)$ tends to zero as $\varepsilon\to 0$.

We now specialize this result to the case $N_{e}=N_{o}=1$, where we can apply Lemma 1, under the assumptions (9) and (24), and we can use the expansion (90) of the eigenfunctions. The goal is to derive the results stated in section 3.2, in the three coupling regimes defined by the scaling relations (35–37).

6.1 Moderate coupling regime

In the regime defined by the scaling relation (35), the expression (137) becomes

[TABLE]

with $\phi$ defined in (17) and $\beta$ defined in Lemma 1. It corresponds to a single guided wave in each waveguide, with transverse profile defined by $\phi(|x|)$ and range dependent random amplitude

[TABLE]

This is the expression of the wave field given in section 3.2, where we simplified the notation as $a_{t,1}(z)\leadsto a_{t}(z)$ and introduced

[TABLE]

satisfying $1\gg|\Delta\beta_{t}|\gg\varepsilon^{2},$ by equations (35) and (89).

The statistics of the amplitudes (139) follows from Theorem 2, in the case $N_{e}=N_{o}=1$. The coefficients in the infinitesimal generator ${\cal L}$ are described in the next corollary, proved in Appendix B.

Corollary 3**.**

The effective coefficients in the operator ${\cal L}^{1}$ are given by

[TABLE]

and

[TABLE]

where $\Gamma$ has the expression (44). The effective coefficients in the operator ${\cal L}^{2}$ are

[TABLE]

where $\Lambda$ is given by (43) and

[TABLE]

The effective coeffients in the operator ${\cal L}^{3}$ are

[TABLE]

where

[TABLE]

and

[TABLE]

Using Theorem 2 and Corollary 3, we obtain the joint moments of the mode amplitudes $a_{e}(z)$ and $a_{o}(z)$: The mean amplitudes are given by

[TABLE]

with $a_{t}^{(0)}$ defined in (23) for $t\in\{e,o\}$. The mean powers and field covariance are

[TABLE]

and the fourth-order moments are

[TABLE]

We conclude from equations (149–150) and (152) that

[TABLE]

This implies the results (41) and (46), which model the effective leakage of power from the guided modes to the radiation modes. Equations (149–151) also give the expression (48) of the imbalance of power between the waveguides.

Remark 6.1**.**

These results show that of all the coefficients in Corollary 3, $\Gamma$ and $\Lambda$ are the important ones. The coefficient $\Gamma$ determines the coupling between the two guided modes, whereas $\Lambda$ determines the power leakage to the radiation modes. Note from (43)–(44) that $\Gamma$ depends on power spectral density $\widehat{\cal R}$, the Fourier transform of the covariance ${\cal R}$, evaluated at zero wavenumber, while $\Lambda$ depends on $\widehat{\cal R}$ evaluated at wavenumbers larger than $\beta-k$, where we recall that $\beta\in(k,nk)$. This result is consistent with the observation in [21] that radiation losses should only depend on the power spectral density evaluated at wavenumbers between $\beta-k$ and $\beta+k$.

It is possible to encounter situations where $\Gamma>0$ and $\Lambda=0$, or the opposite. Indeed, if the power spectral density vanishes for wavenumbers larger than $\beta-k$, then $\Lambda=0$. Otherwise, if $\widehat{\cal R}(0)=0$ (this happens for instance when the random processes $\nu_{j}$ are derivatives of stationary processes), then $\Gamma=0$. However, in general, both $\Gamma$ and $\Lambda$ are positive.

6.2 Weak coupling regime

In the regime defined by the scaling relation (36), the wave field has the same expression (138), but the amplitudes $u_{\pm}(z)$ have different statistics. They are obtained as the $\varepsilon\to 0$ limit of

[TABLE]

where we used the expression of $\beta_{e,1}$ and $\beta_{o,1}$ given in equations (89) and set $\exp(-\eta d)=\varepsilon^{2}\theta$ with $\theta\in[0,\infty)$. To obtain this limit, we let

[TABLE]

and then use equations (128)–(129) to derive a closed system of stochastic differential equations satisfied by the process \big{(}\alpha_{t}^{(\varepsilon)}(z),\big{(}a_{t,\gamma}^{(\varepsilon)}(z)\big{)}_{\gamma\in(0,k^{2})},t\in\{e,o\}\big{)}. The limit of this process is given in the next theorem, with proof outlined in Appendix C. As was the case in Theorem 2, the infinitesimal generator does not involve the radiation mode amplitudes, so we describe directly the limit of \big{(}(\alpha_{t}^{(\varepsilon)}(z),t\in\{e,o\}\big{)}.

Theorem 4**.**

In the limit $\varepsilon\to 0$, the process \big{(}{\alpha}_{e}^{(\varepsilon)}(z),{\alpha}_{o}^{(\varepsilon)}(z)\big{)} converges in distribution in {\cal C}^{0}\big{(}[0,{L}],\mathbb{C}^{2}\big{)} to the Markov process ${\bm{\alpha}}(z)=(\alpha_{e}(z),\alpha_{o}(z))$ with infinitesimal generator

[TABLE]

where $\theta$ is defined in (51), and the operators are defined by

[TABLE]

and

[TABLE]

with coefficients given in Corollary 3.

The result of Theorem 4 can be reformulated as stated in the following corollary, by the change of variables:

[TABLE]

Corollary 5**.**

In the limit $\varepsilon\to 0$, the process $(u_{+}^{(\varepsilon)}(z),u_{-}^{(\varepsilon)}(z))$ converges in distribution in ${\cal C}^{0}([0,L],\mathbb{C}^{2})$ to the Markov process ${\itbf u}(z)=(u_{+}(z),u_{-}(z))$ with infinitesimal generator ${\cal L}_{{\itbf u},\theta}={\cal L}_{u_{+}}+{\cal L}_{u_{-}}+\theta{\cal L}_{u_{+},u_{-}},$ defined by

[TABLE]

Note that for any test function $F$, the generator gives

[TABLE]

This shows that $|u_{+}(z)|^{2}+|u_{-}(z)|^{2}=2(|\alpha_{o}(z)|^{2}+|\alpha_{e}(z)|^{2})$ is deterministic and satisfies (46), as stated in section 3.2.2. The total power transported by the two guided modes decays as $\exp(-\Lambda z)$ with probability one, due to an effective leakage towards the radiation modes:

[TABLE]

To determine the imbalance of power between the two guided modes

[TABLE]

we apply the infinitesimal generator ${\cal L}_{{\itbf u},\theta}$ to the two functions

[TABLE]

and we get

[TABLE]

Denoting ${\cal I}(z)=2{\rm Im}(\overline{u_{+}}u_{-}(z))/(|u_{+}(z)|^{2}+|u_{-}(z)|^{2})$, we find that

[TABLE]

which gives the harmonic oscillator equation (50) satisfied by $\mathbb{E}[\mathcal{P}(z)]$.

6.3 Very weak coupling regime

When the separation $d$ between the waveguides is so large that $\exp(-\eta d)\ll\varepsilon^{2}$, we obtain from (89) that

[TABLE]

The $\varepsilon\to 0$ limit in this case is similar to that in the previous section, and the result is the same as setting $\beta^{\prime}=0$ in Theorem 4. The expression (138) of the wave field still holds, but the mode amplitudes have different statistics, as described in the next corollary:

Corollary 6**.**

In the limit $\varepsilon\to 0$, the process $(u_{+}^{(\varepsilon)}(z),u_{-}^{(\varepsilon)}(z))$ converges in distribution in ${\cal C}^{0}([0,{L}],\mathbb{C}^{2})$ to the process $(u_{+}(z),u_{-}(z))$ where $u_{+}(z)$ and $u_{-}(z)$ are independent and identically distributed Markov processes with infinitesimal generator ${\cal L}_{u}$ defined by (162).

Note that, for any test function $F$, we have

[TABLE]

This shows that the processes $(|u_{+}(z)|^{2},|u_{-}(z)|^{2})$ are deterministic and exponentially decaying with the rate $\Lambda$. We conclude that $(|u_{+}^{(\varepsilon)}|^{2}(z),|u_{-}^{(\varepsilon)}|^{2}(z))$ converges in probability to the deterministic function $(|u_{+}(z)|^{2},|u_{-}(z)|^{2})$ that decays exponentially as $\exp(-\Lambda z)$. Therefore, there is no mode coupling in this regime, except between guided and radiation modes, which results in effective leakage. In particular, the imbalance of power between the two waveguides is constant, as stated in section 3.2.3.

7 Summary

In this paper we introduced an analysis of wave propagation in a directional coupler consisting of two parallel step-index waveguides. The waveguide effect is due to a medium of high index of refraction separated from a uniform background by randomly fluctuating interfaces. The fluctuations occur on a length scale (correlation length) that is similar to the wavelength and have small amplitude modeled by a dimensionless parameter $\varepsilon$ satisfying $0<\varepsilon\ll 1$. The analysis is based on the decomposition of the wave field in a complete set of guided, radiation and evanescent modes. It accounts for the interaction of all these modes and derives a closed system of stochastic differential equations for the guided and radiation mode amplitudes. These equations are driven by the random fluctuations of the interfaces and model the net scattering effect that becomes significant at distances of propagation of the order of $\varepsilon^{-2}$. We analyze them in the asymptotic limit $\varepsilon\to 0$, under the assumption that the covariance function of the fluctuations is smooth. This allows us to use the forward scattering approximation and characterize the limit mode amplitudes as a Markovian process with infinitesimal generator that we calculate explicitly. The analysis applies to waveguides that support an arbitrary number of guided modes. Since many directional couplers use single guided mode waveguides, we studied in detail this case and obtained a detailed quantification of the transfer of power for three different regimes, where the coupling between the guided waves is stronger or weaker, depending on how far apart the waveguides are. In all regimes there is a self-averaging power leakage from the guided to the radiation modes. This is the only coupling if the waveguides are very far apart. Otherwise, the guided modes are coupled and the random boundary fluctuations induce a blurring of the periodic transfer of power that would occur in the absence of the random fluctuations. We quantified this blurring and showed that at sufficiently long distances of propagation the power becomes evenly distributed among the waveguides, independent of the initial condition defined by the wave source.

Acknowledgements

This research is supported in part by the Air Force Office of Scientific Research under award number FA9550-18-1-0131 and in part by the NSF grant DMS1510429.

Appendix A Proof of Lemma 1

The solution $\beta$ of (87) can be obtained as follows: First, note that $cos(q)/q$ is a monotone, strictly decreasing function in the interval $q\in(0,\pi/2)$, with

[TABLE]

Thus, there is a unique solution $q^{\star}$ of

[TABLE]

which gives

[TABLE]

Substituting $q^{\star}=\frac{D}{2}\sqrt{n^{2}k^{2}-\beta^{2}}$ into this equation we obtain that (87) is satisfied, with \beta=\sqrt{k^{2}n^{2}-\Big{(}\frac{2q^{\star}}{D}\Big{)}^{2}}. Obviously, $\beta<kn$ and the condition $\beta>k$ is consistent with

[TABLE]

where the last inequality is by the assumption (9).

A similar construction shows that (88) has no solution in the interval $(k,nk)$. Therefore, the right-hand sides of equations (64) and (71) vanish at the single root $\beta\in(k,nk)$ of (87).

If $d$ is large, then the left-hand sides of equations (64) and (71) are small, and by the implicit function theorem, they have unique solutions, close to $\beta$ calculated above. Consequently, $N_{e}=N_{o}=1$.

Appendix B Proof of Corollary 3

Theorem 2 states that

[TABLE]

where $\beta_{e}\approx\beta_{o}\approx\beta$ and the expectation follows by definitions (101), (17) and equation (90),

[TABLE]

Here the approximation is up to an error of the order of $\exp(-\eta d)\ll 1$, which we neglect in this appendix. Substituting in (170) we get the results (141) and (44). Proceeding similarly, we find $\Gamma^{1}_{t,1,t^{\prime},1}=\Gamma,$ for all $t,t^{\prime}\in\{e,o\}.$

The expression of $\Lambda^{c}_{t,1}$ in Theorem 2 is

[TABLE]

where the expectation is calculated using the definition (97),

[TABLE]

As explained above, $\phi_{t}^{2}\approx\phi^{2}$ for $t\in\{e,o\}$, and definition (17) gives that

[TABLE]

We also have from the expressions of $\phi_{e,\gamma}$ and $\phi_{o,\gamma}$ in section 4.2 that

[TABLE]

with $A_{e,\gamma}$ given by (76). Substituting in (172), we obtain

[TABLE]

As $d$ is large, we take the weak limit $d\to+\infty$ in this expression, seen as a function of $\gamma$, and we get

[TABLE]

by the following calculation: for any $\alpha_{1},\alpha_{2},\alpha_{3},\alpha_{4},\alpha_{5},\alpha_{6}$, with $\alpha_{4}>\sqrt{\alpha_{5}^{2}+\alpha_{6}^{2}}>0$, we have

[TABLE]

This gives $\Lambda_{e,1}^{c}=\Lambda$, as stated in Corollary 3. We can deal with

[TABLE]

and \mathbb{E}\big{[}C_{t,1,o,\gamma}(0)C_{t,1,o,\gamma}(z)\big{]} in the same way, which gives $\Lambda_{o,1}^{c}=\Lambda$. The coefficients $\Lambda_{t,1}^{s}$, $\kappa_{t,1}$, $t\in\{e,o\}$ are obtained in a similar way.

Appendix C Proof of Theorem 4

The result follows from an extended version of Theorem 2, when all $\beta_{t,j}$ are equal, as was considered in [19]. We obtain

[TABLE]

where $c.c.$ stands for complexe conjugate. The statement of Theorem 4 follows from this expression, by using the definition (101) of the random coefficients $C_{t_{1},1,t_{2},1}(z)$ and $C_{t^{\prime},\gamma^{\prime},t_{2},1}(z)$ and the symmetry of the eigenfunctions i.e., that $\phi_{e,1}$ is even and $\phi_{o,1}$ is odd.