Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

Analytical study on the evolutionary behavior of transfer matrix element moments in strongly coupled multimode systems

Open Access Open Access

Abstract

In this paper, the evolutionary behavior of the transfer matrix element moments in strongly coupled multimode systems is studied analytically. The randomly coupled multimode system is modeled by a set of coupled stochastic differential equations (SDEs), which can be used to find the coupled ordinary differential equations (ODEs) for the arbitrary order moments of the matrix elements. Since the ODEs are with the constant coefficients, it is possible to obtain the analytical solutions. The asymptotic behaviors of the solutions are investigated by comparing with the existing results derived from the property of the Haar matrix, and a perfect agreement is observed. The evolutionary behaviors of the transfer matrix element moments computed by the analytical formulas have excellent match with the Monte Carlo simulation results. The analytical method can be highly beneficiary for the multimode system design and analysis.

© 2021 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

1. Introduction

Multimode waveguides, such as the multimode fibers, have received particular attentions recently due to their promising applications in various areas, such as communications and sensing [15]. Inside the multimode waveguides, the modes couple to each other due to the random index fluctuations, which can be either weak or strong, depending on the strength of the fluctuations and the phase matching condition for the modal propagation constant.

The coupling effect can be characterized by a transfer matrix which links the input modal amplitudes to the output modal amplitudes [6]. The matrix elements behave randomly and their moments are important for the average modal power analysis, which have been studied quite extensively during the fiber nonlinear effect analysis [610]. During the numerical integration of the nonlinear Schrödinger equations (NLSE), the nonlinear term should be averaged [610]. Specifically, the fourth order moments [6] directly affect the Kerr nonlinear interaction term. The important propagation equation for the nonlinear coupling inside multimode fibers, i.e. the generalized Manakov equation, has been derived using the statistical property of the 4th order matrix element moments [610]. The multimode Raman amplification process also requires the 4th order transfer matrix element moments to conduct accurate analysis [11,12]. For example, the nonlinear coupling coefficients Qplmn in Eq. (6) of Ref. [11] is vital for the modal amplification characterization, whose averaged value is related to the fourth order moments of the transfer matrix elements. The third example to implement the statistical information of the transfer matrix moments is the multimode Brillouin sensing systems [13]. The random coupling plays an important role in re-allocating the Brillouin modal dependent gain. As a simple case, the Brillouin gain in a fiber with two degenerated linearly polarized modes is highly dependent on the transfer matrix element moments as shown in [14]. The study on the transfer matrix element moments is also essential to characterize the modal dispersion [1517], which is an important parameter in the mode division multiplexing systems and the multimode sensing systems [1822].

In the multimode fibers, the coupling between the different modal groups is weak while the coupling within the same modal group is strong [1822]. The strongly coupled mode system gives rise to the random nature of the transfer matrix and has been paid special attention to during the analysis [710,15,16,18,21]. Up to date, most of the studies on the strongly coupled multimode system have assumed a pure random transfer matrix which is uniformly distributed in the matrix space [610]. In this case, the transfer matrix can be regarded as the Haar matrix and the moments are analytically derived [6,23]. The treatment is valid after the modes have propagated a relatively long distance in the strongly coupled mode system.

While the outstanding precedent works have been carried out, there remains one important problem unsolved: before the matrix element moments reach their eventual asymptotic values, how do they evolve in the strongly coupled mode systems? If the related multimode fiber is relatively short, such evolutionary statistical information becomes vital. This is the case for the multimode parametric amplifiers with short high nonlinear multimode fibers as well as the sensing systems with short fibers. For the system with relatively long fiber, it is also interesting to study the modal power evolutionary behavior at the beginning of the fiber. Although the Monte Carlo simulations can be conducted to obtain the evolutionary behaviors of the transfer matrix element moments [6], it is quite time consuming and lacks physical implications.

In our recent study [24], a stochastic differential equation (SDE) based method was proposed to study the Mueller matrix element moments in the polarization state analysis. In this paper, we propose to study the evolutionary behaviors of the transfer matrix element moments in multimode waveguides/fibers through the SDE approach. With the fundamental SDE for the transfer matrix, the SDE for the arbitrary order moments of the matrix elements can be obtained in a systematic manner. The ordinary differential equations (ODEs) can be obtained afterwards by taking the average on the both sides of the equation. Such ODEs for the moments of the matrix elements fully avoid the repeated computation of many random realizations in the Monte Carlo simulations. In case the random coupling strength is stable along the propagation direction, the coefficients of the ODEs will be constant and the analytical solutions are available. Henceforth, it is possible to study the evolutionary behavior of the random mode coupling effect analytically.

The rest of the paper is organized as follows. Section 2 will provide the basic theory of the SDEs and present the main results, especially the related ODEs for the 2nd order and the 4th order moments of the transfer matrix elements. Section 3 will discuss the asymptotic behaviors of the moments and it is found that the predicted values agree well with the published literatures [6,23]. Section 4 presents the Monte Carlo simulation results and perfect agreement is observed compared with the analytical solutions. Section 5 summarizes the paper.

2. Theory

Here, we adopt the theoretical framework established in [15] for the randomly coupled modes in multimode fibers. Assuming that there exist totally N linearly polarized (LP) modes in the multimode fiber with each of them divided into 2 polarized sub-modes (totally 2N modes), the modal propagation inside the randomly coupled multimode fibers can be described by [15]

$$|{\psi (z )} \rangle = {\textbf U}|{\psi (0 )} \rangle ,$$
where $|{\psi (z )} \rangle $ is the vector standing for the amplitudes of the different modes, U the transfer matrix, whose element moments are to be characterized. It is shown in [15,17] that
$$\frac{{d{\textbf U}}}{{dz}} = j\frac{1}{{2N}}({{{{\beta}}_0}{\textbf I} + ({\vec{{\boldsymbol{\mathrm{\beta}} }}\cdot \vec{{\boldsymbol{\mathrm{\Lambda}} }}} )} ){\textbf U},$$
where β0I stands for the common phase in propagation and vector $\vec{{\boldsymbol{\mathrm{\beta}} }}$ with 4N2−1 elements is a generalized birefringence vector containing the birefringence vectors of the N individual propagation modes as well as the coupling between spatial modes, $\vec{{\boldsymbol{\mathrm{\Lambda}} }}$ the matrix vector whose 4N2−1 elements are the 2N×2N matrices characterizing the coupling between polarized modes within one LP mode, and the coupling between the different LP modes. The detailed description can be found in [15] and it is also detailed in Appendix A. One may choose the Gell-Mann matrix [16] as the matrix basis. However, a slightly different formulism should be adopted in Eq. (2).

The 4N2-1 elements of $\vec{{\boldsymbol{\mathrm{\beta}} }}$ are non-correlated Gaussian random processes with identical distribution [15]. The uncorrelated white noise assumption for the generalized birefringence vector β is a common practice, which has been widely adopted by the research community during the study of strongly coupled multimode systems, e.g. [1517,25,26]. Similar treatment has also been implemented during the analysis of the PMD (polarization mode dispersion) impact [27], which has been verified experimentally. The term of β0I may be dropped, because the common phase term does not play a role in the matrix element moment evaluation [15]. In this work, it is not dropped and is combined with $\vec{{\boldsymbol{\mathrm{\beta}} }}$. Similarly, the identity matrix I is combined with the matrix vector $\vec{{\boldsymbol{\mathrm{\Lambda}} }}$

$$\begin{array}{l} \vec{{\boldsymbol{\mathrm{\beta}} ^{\prime}}} = \left( {\begin{array}{{c}} {{{{\beta}}_0}}\\ {\vec{{\boldsymbol{\mathrm{\beta}} }}} \end{array}} \right)\\ \vec{{\boldsymbol{\mathrm{\Lambda}} ^{\prime}}} = \left( {\begin{array}{{c}} {\textbf I}\\ {\vec{{\boldsymbol{\mathrm{\Lambda}} }}} \end{array}} \right), \end{array}$$

Therefore, Eq. (2) becomes

$$\begin{array}{l} \frac{{d{\textbf U}}}{{dz}} = j{\textbf TU}\\ {\textbf T} = j\frac{1}{{2N}}\left( {{{{\beta}} _0}{\textbf I} + \left( {{{\vec{\boldsymbol{\mathrm{\beta}}} }}\cdot {{\vec{\boldsymbol{\mathrm{\Lambda}}} }}} \right)} \right) = \frac{1}{{2N}}{{\vec{\boldsymbol{\mathrm{\beta}}} }^\prime}\cdot {{\vec{\boldsymbol{\mathrm{\Lambda}}} }^\prime}, \end{array}$$

Since the elements of $\vec{{\boldsymbol{\mathrm{\beta}} }}$ are non-correlated Gaussian random processes with identical distribution and the value of β0 is not relevant (cancelled due to the presence of Ui and Uj* in the moments), one may assume β0 to be an identically distributed Gaussian process as the elements of $\vec{{\boldsymbol{\mathrm{\beta}} }}$ and is not correlated with any elements of $\vec{{\boldsymbol{\mathrm{\beta}} }}$ without loss of generality. Henceforth, we have

$$\begin{array}{l} \left\langle {{{{{\beta}}^{\prime}}_m}(z ){{{{\beta}}^{\prime}}_n}({z^{\prime}} )} \right\rangle = \frac{{{\sigma ^2}}}{{4{N^2}}}{\delta _{mn}}\delta ({z - z^{\prime}} )\\ \left\langle {\vec{{\boldsymbol{\mathrm{\beta}}^{\prime}}}\cdot \vec{{\boldsymbol{\mathrm{\beta}}^{\prime}}}} \right\rangle = {\sigma ^2}\delta ({z - z^{\prime}} ), \end{array}$$
where σ2 is the variance of the random vector, which is related to the mean-square magnitude of the generalized birefringence parameter γ defined in [15]. It is shown in [15] that $\left\langle {\vec{{\boldsymbol{\mathrm{\beta}} }}\cdot \vec{{\boldsymbol{\mathrm{\beta}} }}} \right\rangle = {\omega ^2}{\gamma ^2}\delta ({z - z^{\prime}} )$, where $\omega $ is the angular frequency away from the reference signal center wavelength, and γ is the mean-square magnitude of the generalized birefringence parameter which is similar to the PMD value in the birefringent fiber. It is worth noting that when $\omega = 0$, it does not simply indicate that the randomness disappears. It just suggests that the signal propagates in the same spatial rotational reference frame as the center frequency signal does. One may also choose the reference signal frequency to be 0, like the case in [28,29] to analyze the PMD induced randomness and $\omega $ represents the optical angular frequency in this case. With the above equation and Eq. (5), one has,
$${\sigma ^2} = \frac{{{\omega ^2}{\gamma ^2}4{N^2}}}{{4{N^2} - 1}},$$

According to Eqs. (4) and (5), one may conclude that the elements of T are non-correlated identically distributed Gaussian processes,

$$\left\langle {{t_{mn}}(z ){t_{m^{\prime}n^{\prime}}}^\ast ({z^{\prime}} )} \right\rangle = {\delta _{mm^{\prime}}}{\delta _{nn^{\prime}}}\frac{{{\sigma ^2}}}{{8{N^3}}}\delta ({z - z^{\prime}} ),$$
$$\left\langle {{\textbf T}dz{\textbf T}dz} \right\rangle = \left\langle {{\textbf T}dz{{\textbf T}^H}dz} \right\rangle = \frac{{{\sigma ^2}}}{{4{N^2}}}{\textbf I}dz,$$

The detailed proof of the above conclusion can be found in Appendix A.

Equation (4) is a SDE in the Stratonovich sense, which evaluates the right-hand side (RHS) of the equation in the middle of the integration interval [z, z + dz], i.e. the point of z + dz/2,

$$d{\textbf U} = j{\textbf T}dz{\textbf U}\left( {z + \frac{{dz}}{2}} \right),$$

It should be converted to the Ito sense SDE, which evaluates the RHS of the equation on the left side of the integration interval, i.e. z, to obtain the average value or to derive the SDE for the higher order moments,

$$\begin{aligned} d{\textbf U} &= j{\textbf T}dz{\textbf U}(z )- \frac{1}{2}\left\langle {{\textbf T}dz{\textbf T}dz} \right\rangle {\textbf U}(z )\\ &= j{\textbf T}dz{\textbf U} - \frac{{{\sigma ^2}}}{{8{N^2}}}{\textbf I}dz{\textbf U}, \end{aligned}$$

The above derivation resembles the treatment in [25,26]. The details for the conversion process can be found in Appendix B.

Equation (9) is the SDE which serves as the foundation for the later contents of this work. Based on Eq. (9), it is possible to have the SDE for the 2nd order moments. First of all, we have the following Eq. (10) in the mean square limit sense,

$$d{\textbf U} \otimes d{{\textbf U}^\ast } = \left\langle {{\textbf T}dz \otimes {{\textbf T}^\ast }dz} \right\rangle ({{\textbf U} \otimes {{\textbf U}^\ast }} ),$$

Hence, one has

$$\begin{aligned} & d({{\textbf U} \otimes {{\textbf U}^\ast }} )= d{\textbf U} \otimes {{\textbf U}^\ast } + {\textbf U} \otimes d{{\textbf U}^\ast } + d{\textbf U} \otimes d{{\textbf U}^\ast }\\ &\quad = \left( {j{\textbf T}dz{\textbf U} - \frac{{{\sigma^2}}}{{8{N^2}}}{\textbf I}dz{\textbf U}} \right) \otimes {{\textbf U}^\ast } + {\textbf U} \otimes \left( { - j{{\textbf T}^\ast }dz{{\textbf U}^\ast } - \frac{{{\sigma^2}}}{{8{N^2}}}{\textbf I}dz{{\textbf U}^\ast }} \right) + \left\langle {{\textbf T}dz \otimes {{\textbf T}^\ast }dz} \right\rangle ({{\textbf U} \otimes {{\textbf U}^\ast }} )\\ &\quad = ({j{\textbf T}dz \otimes {\textbf I} - {\textbf I} \otimes j{{\textbf T}^\ast }dz} )({{\textbf U} \otimes {{\textbf U}^\ast }} )+ \left( { - \frac{{{\sigma^2}}}{{4{N^2}}}{\textbf I}dz + \left\langle {{\textbf T}dz \otimes {{\textbf T}^\ast }dz} \right\rangle } \right)({{\textbf U} \otimes {{\textbf U}^\ast }} ), \end{aligned}$$
where ${\otimes}$ denotes matrix Kronecker product. It is worth noting that Eqs. (10) and (11) has been derived based on the important matrix property $({{\textbf A} \otimes {\textbf B}} )({{\textbf C} \otimes {\textbf D}} )= ({{\textbf {AC}}} )\otimes ({{\textbf {BD}}} )$[30]. The average of the first term of the last expression in (11) is zero. Because the term with single Tdz vanishes due to the fact that the average of Tdz is zero i.e. $\left\langle {{\textbf T}dz} \right\rangle = 0$. Furthermore, in the Ito SDE, Tdz is the incremental part in the Markovian process, which is independent with respect to $({{\textbf U} \otimes {{\textbf U}^\ast }} )$. Hence, taking the average on both sides of the equation, one has the ODE for the 2nd order moments as,
$$\begin{array}{l} \frac{{d\left\langle {{\textbf U} \otimes {{\textbf U}^\ast }} \right\rangle }}{{dz}} = {\textbf S}\left\langle {{\textbf U} \otimes {{\textbf U}^\ast }} \right\rangle \\ {\textbf S} ={-} \frac{{{\sigma ^2}}}{{4{N^2}}}{\textbf I} + \frac{{\left\langle {{\textbf T}dz \otimes {{\textbf T}^\ast }dz} \right\rangle }}{{dz}}, \end{array}$$

The first term in S comes from the contribution of the first order moments derivative, while the third order term originates the property of Brownian motion. Physically, the first term stands for the variation caused by the variance of the related element itself, while the second term denotes the variation arises from the coupling effect among the second order moments. All the terms have the weights of σ2. Although dz exists in the expression of S, it will be cancelled due to Eq. (7a). The elements of the second term in S can be calculated in a straight forward way. The matrix has the elements on the ((m-1)N + m)th row and ((n-1)N + n)th column to be σ2, while other elements to be 0.

Based on Eq. (11) and using the similar procedures, one may obtain the SDE for the 4th order moments as follows,

$$\begin{array}{l} d({({{\textbf U} \otimes {{\textbf U}^\ast }} )\otimes ({{\textbf U} \otimes {{\textbf U}^\ast }} )} )\\ = d({{\textbf U} \otimes {{\textbf U}^\ast }} )\otimes ({{\textbf U} \otimes {{\textbf U}^\ast }} )+ ({{\textbf U} \otimes {{\textbf U}^\ast }} )\otimes d({{\textbf U} \otimes {{\textbf U}^\ast }} )+ d({{\textbf U} \otimes {{\textbf U}^\ast }} )\otimes d({{\textbf U} \otimes {{\textbf U}^\ast }} )\\ = ({({j{\textbf T}dz \otimes {\textbf I} - {\textbf I} \otimes j{{\textbf T}^\ast }dz + {\textbf S}dz} )({{\textbf U} \otimes {{\textbf U}^\ast }} )} )\otimes ({{\textbf U} \otimes {{\textbf U}^\ast }} )\\ + ({{\textbf U} \otimes {{\textbf U}^\ast }} )\otimes ({({j{\textbf T}dz \otimes {\textbf I} - {\textbf I} \otimes j{{\textbf T}^\ast }dz + {\textbf S}dz} )({{\textbf U} \otimes {{\textbf U}^\ast }} )} )\\ - ({({{\textbf T}dz \otimes {\textbf I} - {\textbf I} \otimes {{\textbf T}^\ast }dz} )\otimes ({{\textbf T}dz \otimes {\textbf I} - {\textbf I} \otimes {{\textbf T}^\ast }dz} )} )({({{\textbf U} \otimes {{\textbf U}^\ast }} )\otimes ({{\textbf U} \otimes {{\textbf U}^\ast }} )} ), \end{array}$$
where we have used
$$d({{\textbf U} \otimes {{\textbf U}^\ast }} )\otimes d({{\textbf U} \otimes {{\textbf U}^\ast }} )={-} ({({{\textbf T}dz \otimes {\textbf I} - {\textbf I} \otimes {{\textbf T}^\ast }dz} )\otimes ({{\textbf T}dz \otimes {\textbf I} - {\textbf I} \otimes {{\textbf T}^\ast }dz} )} )({({{\textbf U} \otimes {{\textbf U}^\ast }} )\otimes ({{\textbf U} \otimes {{\textbf U}^\ast }} )} ),$$

Taking the average, the term with single Tdz vanishes due to Eq. (7a), and one has

$$\begin{array}{l} \frac{{d\left\langle {({{\textbf U} \otimes {{\textbf U}^\ast }} )\otimes ({{\textbf U} \otimes {{\textbf U}^\ast }} )} \right\rangle }}{{dz}} = {\textbf F}\left\langle {({{\textbf U} \otimes {{\textbf U}^\ast }} )\otimes ({{\textbf U} \otimes {{\textbf U}^\ast }} )} \right\rangle \\ {\textbf F} = {\textbf S} \otimes {\textbf I} + {\textbf I} \otimes {\textbf S} - \frac{{\left\langle {({{\textbf T}dz \otimes {\textbf I} - {\textbf I} \otimes {{\textbf T}^\ast }dz} )\otimes ({{\textbf T}dz \otimes {\textbf I} - {\textbf I} \otimes {{\textbf T}^\ast }dz} )} \right\rangle }}{{dz}}, \end{array}$$

The first term and the second term in F come from the contribution of the second order moments derivative, while the third order term originates from the property of Brownian motion. Physically, the first term and the second term stand for the variation caused by the variance of the related element itself, while the third term denotes the variation arises from the coupling effect among the fourth order moments. All the terms have the weights of σ2. The third term of matrix F can be split into four matrices A, B, C, and D.

$$\begin{array}{l} \frac{{\left\langle {({{\textbf T}dz \otimes {\textbf I} - {\textbf I} \otimes {{\textbf T}^\ast }dz} )\otimes ({{\textbf T}dz \otimes {\textbf I} - {\textbf I} \otimes {{\textbf T}^\ast }dz} )} \right\rangle }}{{dz}}\\ = {\textbf A} + {\textbf B} - {\textbf C} - {\textbf D}\\ {\textbf A} = \frac{{\left\langle {({{\textbf T}dz \otimes {\textbf I}} )\otimes ({{\textbf T}dz \otimes {\textbf I}} )} \right\rangle }}{{dz}}\\ {\textbf B} = \frac{{\left\langle {({{\textbf I} \otimes {{\textbf T}^\ast }dz} )\otimes ({{\textbf I} \otimes {{\textbf T}^\ast }dz} )} \right\rangle }}{{dz}}\\ {\textbf C} = \frac{{\left\langle {({{\textbf T}dz \otimes {\textbf I}} )\otimes ({{\textbf I} \otimes {{\textbf T}^\ast }dz} )} \right\rangle }}{{dz}}\\ {\textbf D} = \frac{{\left\langle {({{\textbf I} \otimes {{\textbf T}^\ast }dz} )\otimes ({{\textbf T}dz \otimes {\textbf I}} )} \right\rangle }}{{dz}}, \end{array}$$

The elements of these matrices A(i,j), B(i,j), C(i,j), D(i,j) can be analytically computed as follows:

$$\begin{array}{l} {\textbf A}({({m - 1} ){{({2N} )}^3} + ({n - 1} ){{({2N} )}^2} + ({l - 1} )2N + k,({m^{\prime} - 1} ){{({2N} )}^3} + ({n^{\prime} - 1} ){{({2N} )}^2} + ({l^{\prime} - 1} )2N + k^{\prime}} )\\ = \left\{ {\begin{array}{{c}} {{\sigma^2} \cdots \cdots \cdots m = l^{\prime},l = m^{\prime},n = n^{\prime},k = k^{\prime}}\\ {0 \cdots \cdots \cdots \cdots \cdots \cdots otherwise} \end{array}} \right.\\ {\textbf B}({({m - 1} ){{({2N} )}^3} + ({n - 1} ){{({2N} )}^2} + ({l - 1} )2N + k,({m^{\prime} - 1} ){{({2N} )}^3} + ({n^{\prime} - 1} ){{({2N} )}^2} + ({l^{\prime} - 1} )2N + k^{\prime}} )\\ = \left\{ {\begin{array}{{c}} {{\sigma^2} \cdots \cdots \cdots m = m^{\prime},l = l^{\prime},n = k^{\prime},k = n^{\prime}}\\ {0 \cdots \cdots \cdots \cdots \cdots \cdots otherwise} \end{array}} \right.\\ {\textbf C}({({m - 1} ){{({2N} )}^3} + ({n - 1} ){{({2N} )}^2} + ({l - 1} )2N + k,({m^{\prime} - 1} ){{({2N} )}^3} + ({n^{\prime} - 1} ){{({2N} )}^2} + ({l^{\prime} - 1} )2N + k^{\prime}} )\\ = \left\{ {\begin{array}{{c}} {{\sigma^2} \cdots \cdots \cdots n = n^{\prime},l = l^{\prime},m = k,m^{\prime} = k^{\prime}}\\ {0 \cdots \cdots \cdots \cdots \cdots \cdots otherwise} \end{array}} \right.\\ {\textbf D}({({m - 1} ){{({2N} )}^3} + ({n - 1} ){{({2N} )}^2} + ({l - 1} )2N + k,({m^{\prime} - 1} ){{({2N} )}^3} + ({n^{\prime} - 1} ){{({2N} )}^2} + ({l^{\prime} - 1} )2N + k^{\prime}} )\\ = \left\{ {\begin{array}{{c}} {{\sigma^2} \cdots \cdots \cdots m = m^{\prime},k = k^{\prime},n = l,n^{\prime} = l^{\prime}}\\ {0 \cdots \cdots \cdots \cdots \cdots \cdots otherwise} \end{array},} \right. \end{array}$$
Where m, n, l, k and m’, n’, l’, k’ are integer numbers starting from 1 to 2N. The derivation of the expressions for the elements of matrices A, B, C, and D is based on Eq. (7a) and the property of the matrix Kronecker product. Since Eqs. (12) and (15) are ODEs with constant coefficients, they can be solved analytically as
$$\begin{array}{l} \left\langle {{\textbf U}(z )\otimes {\textbf U}{{(z )}^\ast }} \right\rangle = \textrm{exp} ({{\textbf S}z} )\left\langle {{\textbf U}(0 )\otimes {\textbf U}{{(0 )}^\ast }} \right\rangle \\ \left\langle {({{\textbf U}(z )\otimes {\textbf U}{{(z )}^\ast }} )\otimes ({{\textbf U}(z )\otimes {\textbf U}{{(z )}^\ast }} )} \right\rangle = \textrm{exp} ({{\textbf F}z} )\left\langle {({{\textbf U}(0 )\otimes {\textbf U}{{(0 )}^\ast }} )\otimes ({{\textbf U}(0 )\otimes {\textbf U}{{(0 )}^\ast }} )} \right\rangle , \end{array}$$
where exp stands for the matrix exponential function, which can be evaluated analytically using the property of Jordan canonical form. At z=0, the transfer matrix T is assumed to be the identity matrix, and therefore the matrices in Eq. (18) should be the identity matrices at z=0 as well.

3. Asymptotic behaviors of the matrix element moments

3.1 2nd order moments

According to Eq. (12), one may find that the matrix S has the following eigen values

$$\frac{{{\sigma ^2}}}{{8{N^3}}}({ - 2N, - 2N \cdots , - 2N,0} ),$$

Matrix S can be decomposed into the summation of its eigen values multiplied by its eigen vectors and their Hermitian transpose. It can be seen from Eq. (19) that one of the eigen values is 0 and the related eigen vector will be associated with the asymptotic behavior during the evaluation of the moments. This is because Eq. (18) suggests the second order moments evolve with respect to the matrix exponential exp(Sz). The eigen values of exp(Sz) is actually the exponential of Eq. (19). Therefore, as the propagation distance increases, the matrix exp(Sz) has most of its eigen values to approach zero, but the last one to approach 1 (exp(0z)). Henceforth, we have

$$\begin{array}{l} \mathop {\lim }\limits_{z \to \infty } \left\langle {{\textbf U}(z )\otimes {\textbf U}{{(z )}^\ast }} \right\rangle \\ = \mathop {\lim }\limits_{z \to \infty } \textrm{exp} ({{\textbf S}z} )\\ = \mathop {\lim }\limits_{z \to \infty } \left( {\sum\limits_{l = 1}^{4{N^2} - 1} {{{\textbf v}_l}\textrm{exp} \left( { - \frac{{{\sigma^2}}}{{4{N^2}}}z} \right){{\textbf v}_l}^H} + {\textbf v}\textrm{exp} ({0z} ){{\textbf v}^H}} \right)\\ = {\textbf v}{{\textbf v}^H}, \end{array}$$
where v is the eigen vector associated with the eigen value 0 of S, vl the eigen vectors associated with the rest negative eigen values of S. Here we have used the property of the matrix exponential that the eigen vector of matrix S maintains as the eigen vector of the matrix exponential [31]. The expression of v is
$${\textbf v}(k )= \left\{ {\begin{array}{{c}} {\frac{1}{{\sqrt {2N} }} \cdots \cdots k = ({m - 1} )2N + m}\\ {0 \cdots \cdots \cdots \cdots \cdots \cdots otherwise} \end{array}} \right.,$$
where m is an integer, m=1,2,…2N. According to Eqs. (20) and (21), it is straight forward to have
$$\mathop {\lim }\limits_{z \to \infty } \left\langle {{{|{{U_{ij}}} |}^2}} \right\rangle = \frac{1}{{2N}},$$
which is in alignment with the Haar matrix property outlined in [19].

3.2 4th order moments

Similar to the procedures above, the matrix F in Eq. (15) has many negative eigenvalues and two zero eigenvalues. Hence, we have

$$\begin{array}{l} \mathop {\lim }\limits_{z \to \infty } \left\langle {({{\textbf U}(z )\otimes {\textbf U}{{(z )}^\ast }} )\otimes ({{\textbf U}(z )\otimes {\textbf U}{{(z )}^\ast }} )} \right\rangle \\ = \mathop {\lim }\limits_{z \to \infty } \textrm{exp} ({{\textbf F}z} )\\ = {{\textbf v}_1}\textrm{exp} ({0z} ){{\textbf v}_1}^H + {{\textbf v}_2}\textrm{exp} ({0z} ){{\textbf v}_2}^H\\ = {{\textbf v}_1}{{\textbf v}_1}^H + {{\textbf v}_2}{{\textbf v}_2}^H, \end{array}$$
where v1 and v2 are the two eigen vectors associated with the two zero eigenvalues of F. Their elements can be analytically expressed as
$$\begin{array}{l} {{\textbf v}_1}({({m - 1} ){{({2N} )}^3} + ({n - 1} ){{({2N} )}^2} + ({l - 1} )2N + k} )= \left\{ {\begin{array}{{l}} {\frac{1}{{\sqrt {2({2N - 1} )2N} }} \cdots m = n,l = k,m \ne l}\\ { - \frac{1}{{\sqrt {2({2N - 1} )2N} }} \cdots m = k,n = l,m \ne n}\\ {0 \cdots \cdots \cdots \cdots \cdots \cdots \cdots \cdots \cdots \cdots \cdots otherwise} \end{array}} \right.\\ {{\textbf v}_2}({({m - 1} ){{({2N} )}^3} + ({n - 1} ){{({2N} )}^2} + ({l - 1} )2N + k} )= \left\{ {\begin{array}{{l}} {\frac{1}{{\sqrt {({2N + 1} )N} }} \cdots \cdots \cdots m = n = l = k}\\ {\frac{1}{{\sqrt {({2N + 1} )4N} }} \cdots \begin{array}{{l}} {m = n,l = k,m \ne l;}\\ {m = k,n = l,m \ne n} \end{array}}\\ {0 \cdots \cdots \cdots \cdots \cdots \cdots \cdots \cdots \cdots \cdots \cdots otherwise} \end{array}} \right. \end{array}$$
where m, n, l, k are integer numbers starting from 1 to 2N. It is concluded from Eq. (23,24) that the asymptotic behaviors of the 4th order matrix element moments can be obtained as follows.
$$\begin{array}{l} \mathop {\lim }\limits_{z \to \infty } \left\langle {{{|{{U_{ij}}} |}^4}} \right\rangle = \frac{2}{{2N({2N + 1} )}}\\ \mathop {\lim }\limits_{z \to \infty } \left\langle {{{|{{U_{ik}}} |}^2}{{|{{U_{il}}} |}^2}} \right\rangle = \frac{1}{{2N({2N + 1} )}}\quad \quad \quad \quad \quad \quad ({k \ne l} )\\ \mathop {\lim }\limits_{z \to \infty } \left\langle {{{|{{U_{ik}}} |}^2}{{|{{U_{jl}}} |}^2}} \right\rangle = \frac{1}{{({2N - 1} )({2N + 1} )}}\quad \quad \quad \quad ({i \ne j,k \ne l} )\\ \mathop {\lim }\limits_{z \to \infty } \left\langle {{U_{ik}}{U_{jl}}{U_{il}}^\ast {U_{jk}}^\ast } \right\rangle = \frac{{ - 1}}{{({2N - 1} )2N({2N + 1} )}}\quad \quad ({i \ne j,k \ne l} ), \end{array}$$
which is in exact agreement with the Haar matrix property outlined in [6,19].

As shown in [6], the nonlinear term on the kth mode $\left\langle {{{\hat{N}}_k}} \right\rangle$ should be expressed as:

$$\left\langle {{{\hat{N}}_k}} \right\rangle = j\gamma \sum\limits_{l = 1}^N {\sum\limits_{m = 1}^N {\sum\limits_{n = 1}^N {\frac{{2\left( {\left\langle {{{|{{U_{ml}}} |}^2}{{|{{U_{nk}}} |}^2}} \right\rangle + \left\langle {{U_{ml}}{U_{nk}}{U_{mk}}^\ast {U_{nl}}^\ast } \right\rangle } \right)}}{{({1 + {d_{kl}}} )({1 + {d_{mn}}} )}}} } } {f_{mmnn}}{|{{A_l}} |^2}{A_k},$$
where N is the number of the modes, Al the amplitude of the lth mode, fmmnn the nonlinear overlap factor, γ the fiber nonlinear coefficient, dkl the number defined as 1 if mode k and mode l belong to the same mode group and 0 otherwise.

Equation (26) suggests our proposed model is readily to be directly implemented for the Kerr effect analysis in an evolutionary sense.

4. Results and discussions

Monte Carlo simulations have been conducted to verify the proposed analytical formulas to compute the evolution of the transfer matrix element moments. The Monte Carlo simulations are conducted numerically based on Eq. (4). In the equation, the vector $\vec{{{\boldsymbol{\mathrm{\beta}}} ^{\prime}}}$ has been modeled as the 4N2 dimensional Brownian motion, which has equal variance for each component and those components are not correlated to each other [15]. 3000 realizations have been computed to obtain the average values for the 2nd order and the 4th order moments of the matrix elements. The propagation distance has been normalized to σ2z/(4N2), with z being the fiber length. The analytical results are computed directly based on Eq. (18). The numerical simulation splits the normalized distance into 0.01 per step. Regarding the required computational time, the analytical method has comparable computational cost with respect to the cost of one realization in the numerical simulations. Therefore, the overall computational cost has been reduced by three orders of magnitude by implementing the analytical expressions. Two multimode systems have been studied. One is a system with one LP mode, which is split into two polarized sub-modes. The other is a system with three LP modes which are split into six polarized sub-modes. They are referred to as the 2-mode system and the 6-mode system in the rest of the paper.

First of all, the results for the 2nd order moments for the 2-mode system are illustrated in Fig. 1. Perfect agreement between the analytical predictions and the Monte Carlo simulations demonstrates the validity of the proposed analytical formulas. The relative differences between the two results for $\left\langle {{{|{{U_{11}}} |}^2}} \right\rangle$, $\left\langle {{{|{{U_{12}}} |}^2}} \right\rangle$, $\left\langle {{{|{{U_{21}}} |}^2}} \right\rangle$, and $\left\langle {{{|{{U_{22}}} |}^2}} \right\rangle$ during the evolution are shown in Fig. 2, which demonstrate that the discrepancies are less than 3%. Since the matrix is unitary, we have ${|{{U_{11}}} |^2} = {|{{U_{22}}} |^2},{|{{U_{12}}} |^2} = {|{{U_{21}}} |^2}$. As the number of realizations in Monte Carlo simulations further increases, the discrepancy can be reduced accordingly. In the rest of the work, the discrepancy between the numerical and analytical results is comparable to the results in Fig. 2. In Fig. 1, It can be seen that the 2nd order moments of the transfer matrix element gradually converge to the value of 1/2, which is in agreement with the asymptotic predictions. The evolutionary behaviors of the curves suggest that the convergence completes when σ2z/(4N2) is around 4. Such value is relatively small. As shown in [15], the mean square of the modal dispersion vector $\vec{\tau }$ is proportional to the fiber length [15,18] and is defined as

$$\left\langle {{\tau^2}} \right\rangle = {\gamma ^2}z = \frac{{{\sigma ^2}z}}{{{\omega ^2}4{N^2}}}({4{N^2} - 1} ),$$

To roughly estimate the fiber length required to reach σ2z/(4N2) =4, a zero reference frequency is chosen. Assuming the signal wavelength to be located at 1550 nm, the mean-square magnitude of the generalized birefringence parameter γ has the value of 1 ps/$\sqrt {km} $, and the number of modes in the system N=3, it can be calculated σ2/(4N2) =42.2/m, σ2=1521/m. σ2z/(4N2) =4 suggests that z is about 0.1 m, when $\sqrt {\left\langle {{\tau^2}} \right\rangle } $ is about 0.01ps. If γ has the value of 0.1 ps/$\sqrt {km} $, σ2/(4N2) =0.422/m, σ2=15.21/m and the corresponding z will be about 10 m when $\sqrt {\left\langle {{\tau^2}} \right\rangle } $ is about 0.01ps.

 figure: Fig. 1.

Fig. 1. The evolution of the 2nd order transfer matrix element moments for the 2-mode system.

Download Full Size | PDF

 figure: Fig. 2.

Fig. 2. The relative differences between the analytical and the numerical results for $\left\langle {{{|{{U_{11}}} |}^2}} \right\rangle$, $\left\langle {{{|{{U_{12}}} |}^2}} \right\rangle$, $\left\langle {{{|{{U_{21}}} |}^2}} \right\rangle$, and $\left\langle {{{|{{U_{22}}} |}^2}} \right\rangle$ during the evolution.

Download Full Size | PDF

While studying the inter-channel nonlinearity effect, i.e., inter-channel cross-phase modulation effect, the neighboring channel center frequency should be chosen as the reference frequency, and the channel angular frequency spacing should be assigned for ω. For example, the channel frequency spacing can be chosen as 100 GHz. In this case, if the mean-square magnitude of the generalized birefringence parameter γ has the value of 1 ps/$\sqrt {km} $, and the number of modes in the system N=3, it can be calculated σ2/(4N2) =1.128×10−5/m, σ2=4.06×10−4/m. σ2z/(4N2) =4 suggests that z is about 10 km, when $\sqrt {\left\langle {{\tau^2}} \right\rangle } $ is about 3.2 ps. If γ has the value of 0.1 ps/$\sqrt {km} $, σ2/(4N2) =1.128×10−7/m, σ2=4.06×10−6/m and the corresponding z will be about 1000 km when $\sqrt {\left\langle {{\tau^2}} \right\rangle } $ is about 3.2 ps.

For the same reference signal center frequency, a larger $\sqrt {\left\langle {{\tau^2}} \right\rangle } $ could correspond to a shorter L. In parameter amplifiers with high nonlinear fibers (HNLFs) or multimode sensing systems, a fiber with the length of tens/hundreds of meters is usually implemented. In the data center transmission links, the fiber length is in the order of several meters. Therefore, the evolutionary behaviors of the 2nd order moments should be negligible in the fiber transmission systems (with the fiber length of hundreds of kms), but be remarkable in the above-mentioned real applications. In those systems, the related statistical power/phase fluctuations induced by the transfer matrix elements variations suggest the designers should consider a more stringent dynamic range and a more robust fault tolerance.

Afterwards, the 4th order matrix element moments of the previous 2-mode system are shown in Figs. 3 and 4. The four important representative 4th order moments mentioned in Eq. (25) are chosen and shown in the figures. Again, perfect agreement between the analytical results and the Monte Carlo simulations is observed. Similar to the trends in Fig. 1, the convergence of the curves completes within the range of [0, 4] for the normalized propagation distance σ2z/(4N2), which further illustrates that the moment evolutionary behavior is significant in short fiber systems. The imaginary part of the dash-dot line in Fig. 4, which is very small, has been ignored for illustration convenience. The asymptotic values of the 4th order moments agree with Eq. (25) and the published literatures [6,23].

 figure: Fig. 3.

Fig. 3. The evolution of the 4th order transfer matrix element moments for the 2-mode system.

Download Full Size | PDF

 figure: Fig. 4.

Fig. 4. The evolution of the 4th order transfer matrix element moments for the 2-mode system.

Download Full Size | PDF

After demonstrating the validity of the formulas for the 2-mode system, the results for the 6-mode system are presented in Figs. (57). Figure 5 shows the 2nd order moments evolutionary behavior for the 6-mode system. Perfect agreement between the analytical results and the Monte Carlo simulation results are demonstrated as expected. The convergence process is accomplished within the range of [0, 2.5] for σ2z/(4N2). The converging values of the 2nd moments meet with the predicted asymptotic values of 1/6. Figures 6 and 7 demonstrate the results of the 4th order matrix element moments for the 6-mode system. Perfect agreement between the analytical predictions and the Monte Carlo averaged values is observed. Like the 2-mode system, the four important representative 4th order moments are shown separately in the two figures. The convergence behavior is similar to the 2nd order case, which is accomplished within the range of [0, 2.5] for σ2z/(4N2). The asymptotic values of the moments agree with the published results. It is worth mentioning that part of the 4th (dashed) line in Fig. 7 has been magnified to see the detailed value in the convergent regime so that the asymptotic behavior can be better observed.

 figure: Fig. 5.

Fig. 5. The evolution of the 2nd order transfer matrix element moments for the 6-mode system.

Download Full Size | PDF

 figure: Fig. 6.

Fig. 6. The evolution of the 4th order transfer matrix element moments for the 6-mode system.

Download Full Size | PDF

 figure: Fig. 7.

Fig. 7. The evolution of the 4th order transfer matrix element moments for the 6-mode system.

Download Full Size | PDF

5. Summary

A novel approach is proposed to analytically study the evolutionary behaviors of the transfer matrix element in multimode systems. The fundamental SDE is derived, which can serve as the basis for the further derivation of the higher order moment equations. ODEs are derived by averaging the SDEs and the analytical solutions to the ODEs are found. The research of the ODEs demonstrates that the asymptotic behaviors of the solutions meet well with the published results, which were derived through a completely different approach. Further verification of the analytical solutions is conducted by comparing with the Monte Carlo simulation results, and a perfect agreement is observed. The proposed method can be useful in the design and analysis for the future multimode fiber links or the multimode fiber sensing systems.

Appendix A

According to [15], the elements of the matrix vector $\vec{{\boldsymbol{\mathrm{\Lambda}} }}$ are matrices generated using the rules described in Appendix A of [15]. It is re-summarized in the following using a different description language.

Denoting the three Pauli matrices as

$${\sigma _1} = \left( {\begin{array}{{cc}} 1&0\\ 0&{ - 1} \end{array}} \right),{\sigma _2} = \left( {\begin{array}{{cc}} 0&1\\ 1&0 \end{array}} \right),{\sigma _3} = \left( {\begin{array}{{cc}} 0&j\\ { - j}&0 \end{array}} \right),$$

Three types of matrices are generated for the elements of the matrix vector $\vec{{\boldsymbol{\mathrm{\Lambda}} }}$:

1, 4N2-4N matrices ${{\boldsymbol{\mathrm{\Lambda}}}_{3N + 1}}$ to ${{\boldsymbol{\mathrm{\Lambda}}}_{4{N^2} - N}}$ are generated by setting the elements which do not belong to the 2X2 block matrix in main diagonal, with 1 and 1 or j and -j in the symmetrical or anti-symmetrical manner. The rest of the elements are set to be zero. The normalization factor is $\sqrt N $

2, 3N matrices ${{\boldsymbol{\mathrm{\Lambda}}}_1}$ to ${{\boldsymbol{\mathrm{\Lambda}}}_{3N}}$ are generated by setting the three Pauli matrices in the main diagonal from the top left to the bottom right. The normalization factor is $\sqrt N $.

Hence, 2N of them associated with the two Pauli matrices σ2 and σ3 can be attributed to the class 1 matrices which has the off-diagonal elements be arranged with 1 and 1 or j and -j. Together, one has 4N2-2N matrices.

N of them are in the form of

$$\sqrt N \left( {\begin{array}{{cccccc}} \ddots &{}&{}&{}&{}&{}\\ {}&0&{}&{}&{}&{}\\ {}&{}&1&{}&{}&{}\\ {}&{}&{}&{ - 1}&{}&{}\\ {}&{}&{}&{}&0&{}\\ {}&{}&{}&{}&{}& \ddots \end{array}} \right),$$
3, N-1 matrices are generated in the form of
$$\begin{array}{l} {{\boldsymbol{\mathrm{\Lambda}}}_{4{N^2} - N + n}} = \sqrt {\frac{N}{{n({n + 1} )}}} \left( {\begin{array}{{cccccccc}} 1&{}&{}&{}&{}&{}&{}&{}\\ {}&1&{}&{}&{}&{}&{}&{}\\ {}&{}& \ddots &{}&{}&{}&{}&{}\\ {}&{}&{}&{ - n}&{}&{}&{}&{}\\ {}&{}&{}&{}&{ - n}&{}&{}&{}\\ {}&{}&{}&{}&{}&0&{}&{}\\ {}&{}&{}&{}&{}&{}&0&{}\\ {}&{}&{}&{}&{}&{}&{}& \ddots \end{array}} \right),\\ n = 1 \cdots N - 1 \end{array}$$
4, Combining the matrix vector $\vec{{\boldsymbol{\mathrm{\Lambda}} }}$ with the last matrix ${{\boldsymbol{\mathrm{\Lambda}}}_0}$, one gets ${\vec{\boldsymbol{\mathrm{\Lambda}}}^{\prime}}$ with 4N2 elements. ${{\boldsymbol{\mathrm{\Lambda}}}_0}$ is the identity matrix and it is with the normalization factor of 1.

The elements of the matrix vector $\vec{{\boldsymbol{\mathrm{\Lambda}} }}$ when N=3 are listed in Appendix C.

Now considering the matrix

$${\textbf T} = \frac{1}{{2N}} \vec{{\boldsymbol{\mathrm{\beta}} }}^{\prime}\cdot \vec{\boldsymbol{\mathrm{\Lambda }}}^{\prime} = \frac{1}{{2N}}\sum\limits_{n = 0}^{4{N^2} - 1} {{{{\beta}} _n}{{\boldsymbol{\mathrm{\Lambda}}}_n}} ,$$

It is evident that the 4N2-2N off-diagonal elements of the Hermitian symmetrical matrix T are generated by the 4N2-4N matrices created following rule 1, and the 2N matrices created following rule 2 with the Pauli matrices σ2 and σ3. Therefore, each off-diagonal element on the upper triangle of the matrix is formed by two independent and identically distributed Gaussian random processes, which are multiplied by 1 and j, i.e. $\frac{1}{{2N}}\left( {\sqrt N {{{\beta}}_l} + j\sqrt N {{{\beta}}_k}} \right)$. The lower triangle of the matrix is Hermitian symmetrical with respect to the upper triangle. Considering $\left\langle {{{{{\beta}}^{\prime}}_m}(z ){{{{\beta}}^{\prime}}_n}({z^{\prime}} )} \right\rangle = \frac{{{\sigma ^2}}}{{4{N^2}}}{\delta _{mn}}\delta ({z - z^{\prime}} )$, it is straight forward to have

$$\begin{array}{l} \left\langle {{t_{mn}}(z ){t_{m^{\prime}n^{\prime}}}^\ast ({z^{\prime}} )} \right\rangle = {\delta _{mm^{\prime}}}{\delta _{nn^{\prime}}}\frac{{{\sigma ^2}}}{{8{N^3}}}\delta ({z - z^{\prime}} ),\\ m \ne n,m^{\prime} \ne n^{\prime} \end{array}$$

These off-diagonal elements are also not correlated with the diagonal elements, which means

$$\begin{array}{l} \left\langle {{t_{mn}}(z ){t_{m^{\prime}m^{\prime}}}^\ast ({z^{\prime}} )} \right\rangle = 0,\\ m \ne n \end{array}$$

The diagonal terms are calculated by

$$\frac{1}{{2N}}\left( {\sum\limits_{n = 0}^N {{{{\beta}}_n}{{\boldsymbol{\mathrm{\Lambda}}}_n}} + \sum\limits_{n = 4{N^2} - N + 1}^{4{N^2} - 1} {{{{\beta}}_n}{{\boldsymbol{\mathrm{\Lambda}}}_n}} } \right),$$

Considering the (2n+1)th and the (2n+2)th diagonal elements, one has

$$\begin{array}{l} {t_{({2n + 1} )({2n + 1} )}}(z )= \frac{1}{{2N}}\left( {{{{\beta}}_0} + \sqrt N {{{\beta}}_n} + \sum\limits_{k = n + 1}^{N - 1} {\sqrt {\frac{N}{{k({k + 1} )}}} {{{\beta}}_{4{N^2} - N + k}}} + \sqrt {\frac{N}{{n({n + 1} )}}} n{{{\beta}}_{4{N^2} - N + n}}} \right)\\ {t_{({2n + 2} )({2n + 2} )}}(z )= \frac{1}{{2N}}\left( {{{{\beta}}_0} - \sqrt N {{{\beta}}_n} + \sum\limits_{k = n + 1}^{N - 1} {\sqrt {\frac{N}{{k({k + 1} )}}} {{{\beta}}_{4{N^2} - N + k}}} + \sqrt {\frac{N}{{n({n + 1} )}}} n{{{\beta}}_{4{N^2} - N + n}}} \right), \end{array}$$
and hence
$$\begin{array}{l} \left\langle {{t_{({2n + 1} )({2n + 1} )}}(z ){t_{({2n + 1} )({2n + 1} )}}^\ast ({z^{\prime}} )} \right\rangle = \left\langle {{t_{({2n + 2} )({2n + 2} )}}(z ){t_{({2n + 2} )({2n + 2} )}}^\ast ({z^{\prime}} )} \right\rangle \\ = \frac{{{\sigma ^2}}}{{4{N^2}}}\frac{1}{{4{N^2}}}\left( {1 + N + \frac{{nN}}{{({n + 1} )}} + \sum\limits_{k = n + 1}^{N - 1} {\frac{N}{{k({k + 1} )}}} } \right)\delta ({z - z^{\prime}} )\\ = \frac{{{\sigma ^2}}}{{4{N^2}}}\frac{1}{{4{N^2}}}2N\delta ({z - z^{\prime}} )\\ = \frac{{{\sigma ^2}}}{{8{N^3}}}\delta ({z - z^{\prime}} )\\ \left\langle {{t_{({2n + 1} )({2n + 1} )}}(z ){t_{({2n + 2} )({2n + 2} )}}^\ast ({z^{\prime}} )} \right\rangle \\ = \frac{{{\sigma ^2}}}{{4{N^2}}}\frac{1}{{4{N^2}}}\left( {1 - N + \frac{{nN}}{{({n + 1} )}} + \sum\limits_{k = n + 1}^{N - 1} {\frac{N}{{k({k + 1} )}}} } \right)\delta ({z - z^{\prime}} )\\ = 0, \end{array}$$

It is also straight forward to conclude that the correlation function between the two non-adjacent main diagonal elements is zero.

Therefore, we may conclude that

$$\left\langle {{t_{mn}}(z ){t_{m^{\prime}n^{\prime}}}^\ast ({z^{\prime}} )} \right\rangle = {\delta _{mm^{\prime}}}{\delta _{nn^{\prime}}}\frac{{{\sigma ^2}}}{{8{N^3}}}\delta ({z - z^{\prime}} ),$$
which is Eq. (7a). Henceforth, one has
$$\left\langle {{\textbf T}(z ){{\textbf T}^H}({z^{\prime}} )} \right\rangle = \frac{{{\sigma ^2}}}{{4{N^2}}}{\textbf I}\delta ({z - z^{\prime}} ),$$

We assume that a small step $\Delta z$ is multiplied on T. From Eq. (38) and property of the Brownian motion, one has:

$$\begin{array}{l} \left\langle {{\textbf T}\Delta z{{\textbf T}^H}\Delta z} \right\rangle \approx \left\langle {\int\limits_0^{\Delta z} {{\textbf T}(z )dz} \int\limits_0^{\Delta z} {{{\textbf T}^H}({z^{\prime}} )dz^{\prime}} } \right\rangle \\ = \int\limits_0^{\Delta z} {dz} \int\limits_0^{\Delta z} {\left\langle {{\textbf T}(z ){{\textbf T}^H}({z^{\prime}} )} \right\rangle dz^{\prime}} \\ = \frac{{{\sigma ^2}}}{{4{N^2}}}{\textbf I}\int\limits_0^{\Delta z} {dz\int\limits_0^{\Delta z} {\delta ({z - z^{\prime}} )dz^{\prime}} } \\ = \frac{{{\sigma ^2}}}{{4{N^2}}}{\textbf I}\Delta z, \end{array}$$
one can deduce that
$$\left\langle {{\textbf T}dz{\textbf T}dz} \right\rangle = \left\langle {{\textbf T}dz{{\textbf T}^H}dz} \right\rangle = \frac{{{\sigma ^2}}}{{4{N^2}}}{\textbf I}dz,$$
which is Eq. (7b).

Appendix B

In this appendix, the conversion of the fundamental SDE from the Stratonovich sense to the Ito sense is presented. The key point for the conversion is that the Stratonovich sense the right-hand side (RHS) of the SDE is evaluated in the middle of the integration interval [z z + dz], i.e. z + dz/2, while the Ito sense SDE evaluates the RHS of the SDE on the left side of the integration interval, i.e. z.

To start with, we expand the Stratonovich sense SDE by the Taylor series

$$\begin{aligned} d{\textbf U} &= j{\textbf {TU}}\left( {z + \frac{{dz}}{2}} \right)dz\\ &= j{\textbf T}dz{\textbf U}(z )+ j{\textbf T}dz\frac{{d{\textbf U}}}{2} + O({dz} ), \end{aligned}$$

Since

$$\begin{aligned} j{\textbf T}dz\frac{{d{\textbf U}}}{2} &= j\frac{1}{2}{\textbf T}dz\left( {j{\textbf T}dz{\textbf U}(z )+ j{\textbf T}dz\frac{{d{\textbf U}}}{2} + O({dz} )} \right)\\ &={-} \frac{1}{2}{\textbf T}dz{\textbf T}dz{\textbf U}(z )+ O({dz} ), \end{aligned}$$

We have

$$d{\textbf U} = j{\textbf T}dz{\textbf U}(z )- \frac{1}{2}{\textbf T}dz{\textbf T}dz{\textbf U}(z )+ O({dz} ),$$

By omitting the higher order term O(dz), the above equation can be converted to the Ito sense in the mean square limit sense, i.e.

$$d{\textbf U} = j{\textbf T}dz{\textbf U}(z )- \frac{1}{2}\left\langle {{\textbf T}dz{\textbf T}dz} \right\rangle {\textbf U}(z ),$$
which is Eq. (9)

Appendix C

Here we list all the elements of the matrix vector $\vec{{\boldsymbol{\mathrm{\Lambda}} }}$ when N=3.

$$\begin{array}{@{}l@{}} {{\Lambda }_1} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 1&0&0&0&0&0\\ 0&{ - 1}&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 \end{array}} \right),{{\Lambda }_2} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&1&0&0&0\\ 0&0&0&{ - 1}&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_3} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 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&0\\ 0&0&0&0&0&{ - 1} \end{array}} \right)\\[4pt] {{\Lambda }_4} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&1&0&0&0&0\\ 1&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&0 \end{array}} \right),{{\Lambda }_5} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&1&0&0\\ 0&0&1&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_6} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 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&0&1\\ 0&0&0&0&1&0 \end{array}} \right)\\[4pt] {{\Lambda }_7} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&j&0&0&0&0\\ { - j}&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&0 \end{array}} \right),{{\Lambda }_8} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&j&0&0\\ 0&0&{ - j}&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_9} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 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&0&j\\ 0&0&0&0&{ - j}&0 \end{array}} \right)\\[4pt] {{\Lambda }_{10}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&1&0&0&0\\ 0&0&0&0&0&0\\ 1&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{11}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&j&0&0&0\\ 0&0&0&0&0&0\\ { - j}&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{12}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&1&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 1&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right)\\[4pt] {{\Lambda }_{13}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&j&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ { - j}&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{14}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&1&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 1&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{15}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&j&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ { - j}&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right)\\[4pt] {{\Lambda }_{16}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&1\\ 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&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{17}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&j\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ { - j}&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{18}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&1&0&0&0\\ 0&1&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right), \end{array}$$
$$\begin{array}{@{}l@{}} {{\Lambda }_{19}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&j&0&0&0\\ 0&{ - j}&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{20}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&1&0&0\\ 0&0&0&0&0&0\\ 0&1&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{21}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&j&0&0\\ 0&0&0&0&0&0\\ 0&{ - j}&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right)\\[4pt] {{\Lambda }_{22}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&1&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&1&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{23}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&j&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&{ - j}&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{24}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&1\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&1&0&0&0&0 \end{array}} \right)\\[4pt] {{\Lambda }_{25}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&j\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&{ - j}&0&0&0&0 \end{array}} \right),{{\Lambda }_{26}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&1&0\\ 0&0&0&0&0&0\\ 0&0&1&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{27}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&j&0\\ 0&0&0&0&0&0\\ 0&0&{ - j}&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right)\\[4pt] {{\Lambda }_{26}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&1\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&1&0&0&0 \end{array}} \right),{{\Lambda }_{27}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&j\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&{ - j}&0&0&0 \end{array}} \right),{{\Lambda }_{28}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&1&0\\ 0&0&0&1&0&0\\ 0&0&0&0&0&0 \end{array}} \right)\\[4pt] {{\Lambda }_{29}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&j&0\\ 0&0&0&{ - j}&0&0\\ 0&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{30}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&1\\ 0&0&0&0&0&0\\ 0&0&0&1&0&0 \end{array}} \right),{{\Lambda }_{31}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&j\\ 0&0&0&0&0&0\\ 0&0&0&{ - j}&0&0 \end{array}} \right)\\[4pt] {{\Lambda }_{32}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 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&0&0&0&0\\ 1&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{33}} = \sqrt 3 \left( {\begin{array}{{@{}cccccc@{}}} 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&j&0&0&0&0\\ { - j}&0&0&0&0&0 \end{array}} \right),{{\Lambda }_{34}} = \sqrt {\frac{3}{2}} \left( {\begin{array}{{@{}cccccc@{}}} 1&0&0&0&0&0\\ 0&1&0&0&0&0\\ 0&0&{ - 1}&0&0&0\\ 0&0&0&{ - 1}&0&0\\ 0&0&0&0&0&0\\ 0&0&0&0&0&0 \end{array}} \right)\\[4pt] {{\Lambda }_{35}} = \sqrt {\frac{1}{2}} \left( {\begin{array}{{@{}cccccc@{}}} 1&0&0&0&0&0\\ 0&1&0&0&0&0\\ 0&0&1&0&0&0\\ 0&0&0&1&0&0\\ 0&0&0&0&{ - 2}&0\\ 0&0&0&0&0&{ - 2} \end{array}} \right),{{\Lambda }_0} = \left( {\begin{array}{{@{}cccccc@{}}} 1&0&0&0&0&0\\ 0&1&0&0&0&0\\ 0&0&1&0&0&0\\ 0&0&0&1&0&0\\ 0&0&0&0&1&0\\ 0&0&0&0&0&1 \end{array}} \right). \end{array}$$

Funding

National Natural Science Foundation of China (61775168).

Acknowledgement

The authors would like to thank the two anonymous reviewers for their helpful suggestions during the review process.

Disclosures

The authors declare no conflicts of interest.

References

1. D. Richardson, J. Fini, and L. Nelson, “Space-division multiplexing in optical fibres,” Nat. Photonics 7(5), 354–362 (2013). [CrossRef]  

2. Z. He, X. Li, M. Luo, R. Hu, C. Li, Y. Qiu, S. Fu, Q. Yang, and S. Yu, “Independent component analysis based channel equalization for 6 × 6 MIMO-OFDM transmission over few-mode fiber,” Opt. Express 24(9), 9209–9217 (2016). [CrossRef]  

3. R. Ryf, M. A. Mestre, S. Randel, X. Palou, A. H. Gnauck, R. Delbue, P. Pupalaikis, A. Sureka, Y. Sun, X. Jiang, and R. Lingle, “Combined SDM and WDM transmission over 700-km Few-Mode Fiber,” in Optical Fiber Communication Conference 2013, OSA Technical Digest (online) (Optical Society of America, 2013), paper OW1I.2.

4. N. Bai, C. Xia, and G. Li, “Adaptive frequency-domain equalization for the transmission of the fundamental mode in a few-mode fiber,” Opt. Express 20(21), 24010–24017 (2012). [CrossRef]  

5. M. J. Murray and B. Redding, “Quantitative strain sensing in a multimode fiber using dual frequency speckle pattern tracking,” Opt. Lett. 45(6), 1309–1312 (2020). [CrossRef]  

6. S. Buch, S. Mumtaz, R.-J. Essiambre, A. M. Tulino, and G. P. Agrawal, “Averaged nonlinear equations for multimode fibers valid in all regimes of random linear coupling,” Opt. Fiber Technol. 48, 123–127 (2019). [CrossRef]  

7. A. Mecozzi, C. Antonelli, and M. Shtaif, “Nonlinear propagation in multi-mode fiber in the strong coupling regime,” Opt. Express 20(11), 11673–11678 (2012). [CrossRef]  

8. A. Mecozzi, C. Antonelli, and M. Shtaif, “Coupled Manakov equations in multimode fibers with strongly coupled groups of modes,” Opt. Express 20(21), 23436–23441 (2012). [CrossRef]  

9. S. Mumtaz, R.-J. Essiambre, and G. P. Agrawal, “Nonlinear propagation in multimode and multicore fibers: generalization of the Manakov equations,” J. Lightwave Technol. 31(3), 398–406 (2012). [CrossRef]  

10. C. Antonelli, M. Shtaif, and A. Mecozzi, “Modeling of nonlinear propagation in space division multiplexed fiber-optic transmission,” J. Lightwave Technol. 34(1), 36–54 (2016). [CrossRef]  

11. F. Poletti and P. Horak, “Description of ultra-short pulse propagation in multi-mode optical fibers,” J. Opt. Soc. Am. B 25(10), 1645–1654 (2008). [CrossRef]  

12. C. Antonelli, A. Mecozzi, and M. Shtaif, “Raman amplification in multimode fibers with random mode coupling,” Opt. Lett. 38(8), 1188–1190 (2013). [CrossRef]  

13. Z. Zhang, Y. Lu, Y. Pan, X. Bao, and L. Chen, “Trench-assisted multimode fiber used in Brillouin optical time domain sensors,” Opt. Express 27(8), 11396–11405 (2019). [CrossRef]  

14. A. Zadok, E. Zilka, A. Eyal, L. Thévenaz, and M. Tur, “Vector analysis of stimulated Brillouin scattering amplification in standard single-mode fibers,” Opt. Express 16(26), 21692–21707 (2008). [CrossRef]  

15. C. Antonelli, A. Mecozzi, M. Shtaif, and P. J. Winzer, “Stokes-space analysis of modal dispersion in fibers with multiple mode transmission,” Opt. Express 20(11), 11718–11733 (2012). [CrossRef]  

16. Q. Hu and W. Shieh, “Autocorrelation function of channel matrix in few-mode fibers with strong mode coupling,” Opt. Express 21(19), 22153–22165 (2013). [CrossRef]  

17. A. Andrusier, M. Shtaif, C. Antonelli, and A. Mecozzi, “Assessing the Effects of Mode-Dependent Loss in Space-Division Multiplexed Systems,” J. Lightwave Technol. 32(7), 1317–1322 (2014). [CrossRef]  

18. K.-P. Ho and J. M. Kahn, “Statistics of Group Delays in Multimode Fiber With Strong Mode Coupling,” J. Lightwave Technol. 29(21), 3119–3128 (2011). [CrossRef]  

19. P. Sillard, D. Molin, M. Bigot-Astruc, K. de Jongh, and F. Achten, “Low-differential-mode-group-delay 9-LP-mode fiber,” in Optical Fiber Communication Conference 2016, OSA Technical Digest (online) (Optical Society of America, 2016), paper M2C.2.

20. S. Berdagué and P. Facq, “Mode division multiplexing in optical fibers,” Appl. Opt. 21(11), 1950–1955 (1982). [CrossRef]  

21. M. B. Shemirani, W. Mao, R. A. Panicker, and J. M. Kahn, “Principal modes in graded-index multimode fiber in presence of spatialand polarization-mode coupling,” J. Lightwave Technol. 27(10), 1248–1261 (2009). [CrossRef]  

22. K. I. Kitayama, S. Sikai, and N. Uchida, “Impulse response prediction based on experimental mode-coupling coefficient in a 10-km long graded-index fiber,” J,” Quantum Electron. 16(3), 356–362 (1980). [CrossRef]  

23. A.M. Tulino and S. Verdú, “Random matrix theory and wireless communications,” Commun. Inf. Theory, now Publishers Inc, 2004.

24. J. Zhou and Q. Hu, “Higher order statistics of the Mueller matrix in a fiber with an arbitrary length impacted by PMD,” Opt. Express 28(20), 30063–30074 (2020). [CrossRef]  

25. A. Mecozzi, C. Antonelli, and M. Shtaif, “Intensity impulse response of SDM links,” Opt. Express 23(5), 5738–5743 (2015). [CrossRef]  

26. C. Antonelli, A. Mecozzi, M. Shtaif, N. K. Fontaine, H. Chen, and R. Ryf, “Stokes-Space Analysis of Modal Dispersion of SDM Fibers With Mode-Dependent Loss: Theory and Experiments,” J. Lightwave Technol. 38(7), 1668–1677 (2020). [CrossRef]  

27. G. J. Foschini and C. D. Poole, “Statistical theory of polarization dispersion in single mode fibers,” J. Lightwave Technol. 9(11), 1439–1456 (1991). [CrossRef]  

28. Q. Lin and G. P. Agrawal, “Vector theory of stimulated Raman scattering and its application to fiber-based Raman amplifiers,” J. Opt. Soc. Am. B 20(8), 1616–1631 (2003). [CrossRef]  

29. Y. Li and A. Yariv, “Solutions to the dynamical equation of polarization-mode dispersion and polarization dependent losses,” J. Opt. Soc. Am. B 17(11), 1821–1827 (2000). [CrossRef]  

30. A. J. Laub, Matrix Analysis for Scientists and Engineers, (SIAM, 2004, Chap. 13).

31. Dan Klain, “The Matrix Exponential and Linear Systems of ODEs”, http://faculty.uml.edu/dklain/exponential.pdf

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (7)

Fig. 1.
Fig. 1. The evolution of the 2nd order transfer matrix element moments for the 2-mode system.
Fig. 2.
Fig. 2. The relative differences between the analytical and the numerical results for $\left\langle {{{|{{U_{11}}} |}^2}} \right\rangle$ , $\left\langle {{{|{{U_{12}}} |}^2}} \right\rangle$ , $\left\langle {{{|{{U_{21}}} |}^2}} \right\rangle$ , and $\left\langle {{{|{{U_{22}}} |}^2}} \right\rangle$ during the evolution.
Fig. 3.
Fig. 3. The evolution of the 4th order transfer matrix element moments for the 2-mode system.
Fig. 4.
Fig. 4. The evolution of the 4th order transfer matrix element moments for the 2-mode system.
Fig. 5.
Fig. 5. The evolution of the 2nd order transfer matrix element moments for the 6-mode system.
Fig. 6.
Fig. 6. The evolution of the 4th order transfer matrix element moments for the 6-mode system.
Fig. 7.
Fig. 7. The evolution of the 4th order transfer matrix element moments for the 6-mode system.

Equations (47)

Equations on this page are rendered with MathJax. Learn more.

| ψ ( z ) = U | ψ ( 0 ) ,
d U d z = j 1 2 N ( β 0 I + ( β Λ ) ) U ,
β = ( β 0 β ) Λ = ( I Λ ) ,
d U d z = j T U T = j 1 2 N ( β 0 I + ( β Λ ) ) = 1 2 N β Λ ,
β m ( z ) β n ( z ) = σ 2 4 N 2 δ m n δ ( z z ) β β = σ 2 δ ( z z ) ,
σ 2 = ω 2 γ 2 4 N 2 4 N 2 1 ,
t m n ( z ) t m n ( z ) = δ m m δ n n σ 2 8 N 3 δ ( z z ) ,
T d z T d z = T d z T H d z = σ 2 4 N 2 I d z ,
d U = j T d z U ( z + d z 2 ) ,
d U = j T d z U ( z ) 1 2 T d z T d z U ( z ) = j T d z U σ 2 8 N 2 I d z U ,
d U d U = T d z T d z ( U U ) ,
d ( U U ) = d U U + U d U + d U d U = ( j T d z U σ 2 8 N 2 I d z U ) U + U ( j T d z U σ 2 8 N 2 I d z U ) + T d z T d z ( U U ) = ( j T d z I I j T d z ) ( U U ) + ( σ 2 4 N 2 I d z + T d z T d z ) ( U U ) ,
d U U d z = S U U S = σ 2 4 N 2 I + T d z T d z d z ,
d ( ( U U ) ( U U ) ) = d ( U U ) ( U U ) + ( U U ) d ( U U ) + d ( U U ) d ( U U ) = ( ( j T d z I I j T d z + S d z ) ( U U ) ) ( U U ) + ( U U ) ( ( j T d z I I j T d z + S d z ) ( U U ) ) ( ( T d z I I T d z ) ( T d z I I T d z ) ) ( ( U U ) ( U U ) ) ,
d ( U U ) d ( U U ) = ( ( T d z I I T d z ) ( T d z I I T d z ) ) ( ( U U ) ( U U ) ) ,
d ( U U ) ( U U ) d z = F ( U U ) ( U U ) F = S I + I S ( T d z I I T d z ) ( T d z I I T d z ) d z ,
( T d z I I T d z ) ( T d z I I T d z ) d z = A + B C D A = ( T d z I ) ( T d z I ) d z B = ( I T d z ) ( I T d z ) d z C = ( T d z I ) ( I T d z ) d z D = ( I T d z ) ( T d z I ) d z ,
A ( ( m 1 ) ( 2 N ) 3 + ( n 1 ) ( 2 N ) 2 + ( l 1 ) 2 N + k , ( m 1 ) ( 2 N ) 3 + ( n 1 ) ( 2 N ) 2 + ( l 1 ) 2 N + k ) = { σ 2 m = l , l = m , n = n , k = k 0 o t h e r w i s e B ( ( m 1 ) ( 2 N ) 3 + ( n 1 ) ( 2 N ) 2 + ( l 1 ) 2 N + k , ( m 1 ) ( 2 N ) 3 + ( n 1 ) ( 2 N ) 2 + ( l 1 ) 2 N + k ) = { σ 2 m = m , l = l , n = k , k = n 0 o t h e r w i s e C ( ( m 1 ) ( 2 N ) 3 + ( n 1 ) ( 2 N ) 2 + ( l 1 ) 2 N + k , ( m 1 ) ( 2 N ) 3 + ( n 1 ) ( 2 N ) 2 + ( l 1 ) 2 N + k ) = { σ 2 n = n , l = l , m = k , m = k 0 o t h e r w i s e D ( ( m 1 ) ( 2 N ) 3 + ( n 1 ) ( 2 N ) 2 + ( l 1 ) 2 N + k , ( m 1 ) ( 2 N ) 3 + ( n 1 ) ( 2 N ) 2 + ( l 1 ) 2 N + k ) = { σ 2 m = m , k = k , n = l , n = l 0 o t h e r w i s e ,
U ( z ) U ( z ) = exp ( S z ) U ( 0 ) U ( 0 ) ( U ( z ) U ( z ) ) ( U ( z ) U ( z ) ) = exp ( F z ) ( U ( 0 ) U ( 0 ) ) ( U ( 0 ) U ( 0 ) ) ,
σ 2 8 N 3 ( 2 N , 2 N , 2 N , 0 ) ,
lim z U ( z ) U ( z ) = lim z exp ( S z ) = lim z ( l = 1 4 N 2 1 v l exp ( σ 2 4 N 2 z ) v l H + v exp ( 0 z ) v H ) = v v H ,
v ( k ) = { 1 2 N k = ( m 1 ) 2 N + m 0 o t h e r w i s e ,
lim z | U i j | 2 = 1 2 N ,
lim z ( U ( z ) U ( z ) ) ( U ( z ) U ( z ) ) = lim z exp ( F z ) = v 1 exp ( 0 z ) v 1 H + v 2 exp ( 0 z ) v 2 H = v 1 v 1 H + v 2 v 2 H ,
v 1 ( ( m 1 ) ( 2 N ) 3 + ( n 1 ) ( 2 N ) 2 + ( l 1 ) 2 N + k ) = { 1 2 ( 2 N 1 ) 2 N m = n , l = k , m l 1 2 ( 2 N 1 ) 2 N m = k , n = l , m n 0 o t h e r w i s e v 2 ( ( m 1 ) ( 2 N ) 3 + ( n 1 ) ( 2 N ) 2 + ( l 1 ) 2 N + k ) = { 1 ( 2 N + 1 ) N m = n = l = k 1 ( 2 N + 1 ) 4 N m = n , l = k , m l ; m = k , n = l , m n 0 o t h e r w i s e
lim z | U i j | 4 = 2 2 N ( 2 N + 1 ) lim z | U i k | 2 | U i l | 2 = 1 2 N ( 2 N + 1 ) ( k l ) lim z | U i k | 2 | U j l | 2 = 1 ( 2 N 1 ) ( 2 N + 1 ) ( i j , k l ) lim z U i k U j l U i l U j k = 1 ( 2 N 1 ) 2 N ( 2 N + 1 ) ( i j , k l ) ,
N ^ k = j γ l = 1 N m = 1 N n = 1 N 2 ( | U m l | 2 | U n k | 2 + U m l U n k U m k U n l ) ( 1 + d k l ) ( 1 + d m n ) f m m n n | A l | 2 A k ,
τ 2 = γ 2 z = σ 2 z ω 2 4 N 2 ( 4 N 2 1 ) ,
σ 1 = ( 1 0 0 1 ) , σ 2 = ( 0 1 1 0 ) , σ 3 = ( 0 j j 0 ) ,
N ( 0 1 1 0 ) ,
Λ 4 N 2 N + n = N n ( n + 1 ) ( 1 1 n n 0 0 ) , n = 1 N 1
T = 1 2 N β Λ = 1 2 N n = 0 4 N 2 1 β n Λ n ,
t m n ( z ) t m n ( z ) = δ m m δ n n σ 2 8 N 3 δ ( z z ) , m n , m n
t m n ( z ) t m m ( z ) = 0 , m n
1 2 N ( n = 0 N β n Λ n + n = 4 N 2 N + 1 4 N 2 1 β n Λ n ) ,
t ( 2 n + 1 ) ( 2 n + 1 ) ( z ) = 1 2 N ( β 0 + N β n + k = n + 1 N 1 N k ( k + 1 ) β 4 N 2 N + k + N n ( n + 1 ) n β 4 N 2 N + n ) t ( 2 n + 2 ) ( 2 n + 2 ) ( z ) = 1 2 N ( β 0 N β n + k = n + 1 N 1 N k ( k + 1 ) β 4 N 2 N + k + N n ( n + 1 ) n β 4 N 2 N + n ) ,
t ( 2 n + 1 ) ( 2 n + 1 ) ( z ) t ( 2 n + 1 ) ( 2 n + 1 ) ( z ) = t ( 2 n + 2 ) ( 2 n + 2 ) ( z ) t ( 2 n + 2 ) ( 2 n + 2 ) ( z ) = σ 2 4 N 2 1 4 N 2 ( 1 + N + n N ( n + 1 ) + k = n + 1 N 1 N k ( k + 1 ) ) δ ( z z ) = σ 2 4 N 2 1 4 N 2 2 N δ ( z z ) = σ 2 8 N 3 δ ( z z ) t ( 2 n + 1 ) ( 2 n + 1 ) ( z ) t ( 2 n + 2 ) ( 2 n + 2 ) ( z ) = σ 2 4 N 2 1 4 N 2 ( 1 N + n N ( n + 1 ) + k = n + 1 N 1 N k ( k + 1 ) ) δ ( z z ) = 0 ,
t m n ( z ) t m n ( z ) = δ m m δ n n σ 2 8 N 3 δ ( z z ) ,
T ( z ) T H ( z ) = σ 2 4 N 2 I δ ( z z ) ,
T Δ z T H Δ z 0 Δ z T ( z ) d z 0 Δ z T H ( z ) d z = 0 Δ z d z 0 Δ z T ( z ) T H ( z ) d z = σ 2 4 N 2 I 0 Δ z d z 0 Δ z δ ( z z ) d z = σ 2 4 N 2 I Δ z ,
T d z T d z = T d z T H d z = σ 2 4 N 2 I d z ,
d U = j TU ( z + d z 2 ) d z = j T d z U ( z ) + j T d z d U 2 + O ( d z ) ,
j T d z d U 2 = j 1 2 T d z ( j T d z U ( z ) + j T d z d U 2 + O ( d z ) ) = 1 2 T d z T d z U ( z ) + O ( d z ) ,
d U = j T d z U ( z ) 1 2 T d z T d z U ( z ) + O ( d z ) ,
d U = j T d z U ( z ) 1 2 T d z T d z U ( z ) ,
Λ 1 = 3 ( 1 0 0 0 0 0 0 1 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 ) , Λ 2 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) , Λ 3 = 3 ( 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 0 0 0 0 0 0 1 ) Λ 4 = 3 ( 0 1 0 0 0 0 1 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 0 ) , Λ 5 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) , Λ 6 = 3 ( 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 0 1 0 0 0 0 1 0 ) Λ 7 = 3 ( 0 j 0 0 0 0 j 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 0 ) , Λ 8 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) , Λ 9 = 3 ( 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 0 j 0 0 0 0 j 0 ) Λ 10 = 3 ( 0 0 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) , Λ 11 = 3 ( 0 0 j 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) , Λ 12 = 3 ( 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) Λ 13 = 3 ( 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) , Λ 14 = 3 ( 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 ) , Λ 15 = 3 ( 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 ) Λ 16 = 3 ( 0 0 0 0 0 1 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 0 0 0 0 0 ) , Λ 17 = 3 ( 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 ) , Λ 18 = 3 ( 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) ,
Λ 19 = 3 ( 0 0 0 0 0 0 0 0 j 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) , Λ 20 = 3 ( 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) , Λ 21 = 3 ( 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) Λ 22 = 3 ( 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 ) , Λ 23 = 3 ( 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 ) , Λ 24 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 ) Λ 25 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 ) , Λ 26 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 ) , Λ 27 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 ) Λ 26 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 ) , Λ 27 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 ) , Λ 28 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0 0 0 0 0 ) Λ 29 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 j 0 0 0 0 0 0 0 0 ) , Λ 30 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 1 0 0 ) , Λ 31 = 3 ( 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 j 0 0 0 0 0 0 0 0 0 j 0 0 ) Λ 32 = 3 ( 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 0 0 0 0 1 0 0 0 0 0 ) , Λ 33 = 3 ( 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 j 0 0 0 0 j 0 0 0 0 0 ) , Λ 34 = 3 2 ( 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 ) Λ 35 = 1 2 ( 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 2 0 0 0 0 0 0 2 ) , Λ 0 = ( 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 1 ) .
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.