Non-Manakovian Propagation in Optical Fiber

Abstract

Solving the nonlinear Schrödinger equation or similar calculus is essential for designing today’s long-haul optical communication systems. Associated numerical and analytical approaches were extensively studied over the past four decades; simplifications and adaptions for various applications and purposes have been introduced. Optical fibers installed in long-haul systems possess nonideal features such as birefringence, which some of these adaptions partially address to improve the simulation accuracy. But as the fiber birefringence frequently and randomly changes along a link, finding a mathematical solution is a more complex problem and beyond the task of dealing with the nonlinear nature of the corresponding equations. Recently, a novel propagation phenomenon related to the polarization evolution of a signal has been observed. In links with considerable length, i.e., bridging transatlantic or transpacific distances, the polarization state of a light wave is impacted by neighboring signals via the Kerr nonlinearity in glass. Established formalisms for describing polarization effects based on the nonlinear Schrödinger equation cannot fully capture this phenomenon. Here we discuss a first-order calculus for this problem. We start with high-level reviews of experimental observations to introduce the phenomenon and ways to model regular nonlinear propagation. Then we present a first-order calculus to describe the statistics behind the phenomenon by specifically discussing the interplay between fiber birefringence and fiber nonlinearities.

1. Introduction

Long-distance optical telecommunications using standard single-mode fiber (SSMF) is economically [1] one of the most important industrial applications, which requires an accurate solution to the nonlinear (NL) Schrödinger equations for developing competitive products. Modern fiber communications, the basis of all backbone networks, enable global long-reach and high-capacity data exchange such as the WWW and state-of-the-art systems provide 10s of Tb/s capacity per fiber over transpacific distances without electrical signal regeneration [2, 3].

By simulating the NL Schrödinger equation or similar approaches, the industry assesses potential commercial link design options for network operations [4, 5]. The chosen optical power levels balance the generated NL signal distortions with sufficiently high optical signal-to-noise ratios (OSNR) to guarantee certain bit error rates (BERs) on the receive side.

Throughout the half-century-long history of SSMF, various mathematical approaches for describing NL signal propagation have been proposed where all of them are somewhat related to the (coupled) NL Schrödinger equation(s) [6, 7]. These models vary in terms of complexity, their intended purposes, and user friendliness. For example, a more refined analysis of light propagation in SSMF for telecom application requires a polarization dependent field representation. Today’s coherent signaling technique, which exploits two orthogonal polarization states at same optical frequency to encode data information and the unavoidable birefringence of SSMF, requires a vector field representation of the light mode. But simultaneously considering both aspects, the polarization of light and its random changes along a regular link due to fiber birefringence lead to an extended set of equations that are for most applications impractical to solve exactly.

The high end of the “manageable” equation sets, the “Manakov-PMD equation” addresses NL propagation in a fiber with randomly varying birefringence [8]. For the current generation of communication systems, it provides quite accurate results for expected NL signal distortions. But recently, a fiber phenomenon, to which we refer to as nonlinear depolarization (NLDP) of light in fiber, has been observed and does not conform to the Manakov-polarization mode dispersion (PMD) equation. Here we discuss some relevant experimental aspects of this phenomenon and its impact on the mathematical description of NL propagation in realistic, i.e., birefringent fiber. The Manakov-PMD equation is in some ways an advanced form of the NL Schrödinger equation as it additionally considers PMD effects [9, 10], and on the other hand, a simplification as it uses averaged quantities over distance to describe the randomly changing fiber birefringence.

In this chapter, we describe an algorithm for solving the NL Schrödinger equation in vector form when the field variables are randomly and rapidly alternated by the stationary and linear properties (birefringence) of the channel. In other words, some coefficients of the coupled NL Schrödinger equation (CNLS) become distance dependent and describe the changing features of the glass medium along the propagation path.

We start with a brief overview on commonly used modeling for NL signal propagation (Manakovian propagation) in SSMF, report on a high-level view of experimental results that do not conform to such formalism (non-Manakovian propagation) and discuss an algorithm that yields some analytic quantities for a theoretical description of the later. Although not fully technically correct, we use the SSMF terminology in our chapter [11] even when we also mean other fiber types such as large effective area fiber (LEAF [12]) that are often installed in long-haul communication systems.

Our focus is to develop a mathematical formalism for NLDP that yields quantities that are observable in typical industrial test beds. Certainly, more sophisticated experimental setups can be built to characterize other features of NLDP. We report the experimental conditions to a degree that produces a qualitative understanding of the phenomenon and will reference further details in the literature.

2. Nonlinear propagation equations

Back in 1972 and far before commercial applications of optical communications became relevant, Sergey Manakov¹ suggested that a careful consideration of NL pulse propagation is required for accurate signal modeling in fiber. He proposed a set of coupled differential equations that to a large degree can capture the impact of the fiber’s Kerr nonlinearity on a signal’s evolution [13]. A slowly varying envelope constituted by the two orthogonal variables A_x and A_y represents the pulse in space and time domains:

$$\frac{\partial A_x}{\partial z}+{\beta }_1\frac{\partial A_x}{\partial t}+\frac{j{\beta }_2}{2}\frac{{\partial }^2{A}_x}{\partial t^2}+\frac{\alpha }{2}{A}_x=\mathit{j\gamma }\left({\left|A_x\right|}^2+{\left|A_y\right|}^2\right)A_x$$E1

$$\frac{\partial A_y}{\partial z}+{\beta }_1\frac{\partial A_y}{\partial t}+\frac{j{\beta }_2}{2}\frac{{\partial }^2{A}_y}{\partial t^2}+\frac{\alpha }{2}{A}_y=\mathit{j\gamma }\left({\left|A_x\right|}^2+{\left|A_y\right|}^2\right)A_y$$E2

where z and t are the propagation direction and time, respectively. The optical features of the fiber are characterized by $a,{\beta }_1,{\beta }_2,\text{and}\gamma$, which refer to the fiber attenuation, propagation constant, group velocity, and Kerr nonlinearity [14], respectively. By multiplying A_x and A_y with the optical carrier, which typically resides around 194 THz (1.5 μm wavelength), the optical field can be obtained.

In the field of telecom the above listed formulas are called the Manakov Equation (ME). Literature uses this name for both its scalar version and a two-dimensional version for polarization representation of a signal. The ME was heuristically found in the sense that the left side, known from linear transmission theory, has been extended by a source term on its right side describing the Kerr nonlinearity. It assumes that the nonlinearities are weak and proportional to the signal’s intensity.

The ME has been highly successfully applied to NL signal propagation and can explain phenomena such as optical solitons and nonlinear polarization rotation (NLPR [15, 16]) that were subjects of intensive research until about the end of the 1990 [17, 18, 19, 20]. But the ME in the form outlined above is restricted in twofold aspects compared with a more accurate propagation modeling that modern telecom applications require:

(a) It assumes a rotational symmetry of the fiber, i.e., without birefringence. But commercial SSMFs, even those with the lowest amounts of birefringence, demand for modeling that adapts the ME’s left side (the linear propagation features) to simulate typical evolutions of a signal’s state of polarization (SOP). We address this issue by substituting ${\beta }_1$ and ${\beta }_2$ with polarization-dependent constants ${\beta }_{1x},{\beta }_{1y},{\beta }_{2x}$ and ${\beta }_{2y}$ and

(b) other than the ME hypothesizes, the fiber Kerr NL of glass is significantly polarization dependent. The impact of an orthogonal polarization component needs to be weighted by a factor of 2/3.

Both adjustments in the ME lead to a simplified version of the so-called coupled NL Schrödinger equations (CNLSs):

$$\frac{\partial A_x}{\partial z}+{\beta }_{1x}\frac{\partial A_x}{\partial t}+\frac{j{\beta }_2}{2}\frac{{\partial }^2{A}_x}{\partial t^2}+\frac{\alpha }{2}{A}_x=\mathit{j\gamma }\left({\left|A_x\right|}^2+\frac{2}{3}{\left|A_y\right|}^2\right)A_x$$E3

$$\frac{\partial A_y}{\partial z}+{\beta }_{1y}\frac{\partial A_y}{\partial t}+\frac{j{\beta }_2}{2}\frac{{\partial }^2{A}_y}{\partial t^2}+\frac{\alpha }{2}{A}_y=\mathit{j\gamma }\left({\left|A_y\right|}^2+\frac{2}{3}{\left|A_x\right|}^2\right)A_y$$E4

The extended version differs by terms on the right side, which model fast oscillations of the nonlinearity (not shown). In our further considerations, those terms can be neglected; and the polarization dependence of ${\beta }_2$ is small and can be ignored. Strictly speaking, the CNLSs hold only for a single piece of fiber with linear birefringence, referred to in the literature as waveplate. Different velocities of light in both optical axes (principal states) of a waveplate lead in general to SOP changes for a cw tone passing through a single plate and complicate the calculus.

A model of a realistic fiber link would typically involve a concatenation of many waveplates, each with their optical axes randomly oriented. Solving the CNLS becomes for such a scenario a cumbersome process. On the other hand, it had been observed that NL differential equations of Manakovian type (Eqs. (1) and (2)) sufficiently capture experimental results such as soliton propagation. This led to the assumption that the ME can be derived from the CNLS when its polarization-dependent NL source term is averaged in a way that it simulates polarization scrambling caused by fiber birefringence [21, 22]. Several derivations have been published, which yield an ME with an 8/9 reduced Kerr nonlinearity as fitting parameter:

$$\frac{\partial A_x}{\partial z}+{\beta }_1\frac{\partial A_x}{\partial t}+\frac{j{\beta }_2}{2}\frac{{\partial }^2{A}_x}{\partial t^2}+\frac{\alpha }{2}{A}_x=\mathit{j\gamma }\frac{8}{9}\left({\left|A_x\right|}^2+{\left|A_y\right|}^2\right)A_x$$E5

$$\frac{\partial A_y}{\partial z}+{\beta }_1\frac{\partial A_y}{\partial t}+\frac{j{\beta }_2}{2}\frac{{\partial }^2{A}_y}{\partial t^2}+\frac{\alpha }{2}{A}_y=\mathit{j\gamma }\frac{8}{9}\left({\left|A_y\right|}^2+{\left|A_x\right|}^2\right)A_y$$E6

Note, the fiber birefringence does not appear anymore in the form of polarization-dependent propagation constants, and the x, y-indexes should not be identified with the coordinates of the lab frame. Nevertheless, advanced effects such as NL polarization rotation of two cw tones can be quantitatively described by such equations.

The heuristically determined ME with and without reduced Kerr nonlinearity has become an integral part of soliton theory and is often used in today’s system modeling due to its relatively simple form and highly accurate results. Based on this success, the published attempts for deriving the ME from the more fundamental CNLSs have gained wide acceptance [21]. Here, however, we challenge their formalism as they rely on an elusively compelling mathematical argument and show that the solution spaces of the ME with reduced Kerr nonlinearity and of the section-wise solved CNLS for long birefringent fibers diverge from each other. This difference matters in case of NLDP.

The Manakov-PMD equation developed toward the end of the 1990 is cumbersome but allows to separate the principal effects by terms such as the Kerr nonlinearity, chromatic dispersion, nonlinear PMD, and linear PMD that contribute to the signal’s evolution. Its nonlinear PMD term can be ignored in applications for today’s regular telecom fibers, thus the remaining NL source term of the Manakov-PMD equation simplifies to the same form as known from the ME (Eqs. (5) and (6)). It differs from the latter by an additional term that accounts for linear PMD contributions. Nonlinear PMD refers to a situation where the fiber birefringence does not sufficiently scramble the SOP of a cw signal when it linearly propagates over distances that are comparable to the NL propagation length of the path. In our following considerations, we can (to good approximation) assume that a signal’s SOP gets strongly scrambled represented by a homogeneous coverage on the Poincare sphere while it propagates just a few kilometers. We can therefore neglect NL PMD.

Nevertheless, due to the same method of averaging that had been applied when deriving from the CNLS the ME with an 8/9 reduced Kerr nonlinearity, the Manakov-PMD equation does not represent NLDP fully correctly, regardless if it includes NL PMD or not.

3. Comparative polarimetry for detecting NLDP

The undersea communication industry continuously aims for a more precise modeling of signal propagation in fiber to enable longer unrepeated system spans and higher transport capacities, hence more cost-effective solutions. However, a noticeable discrepancy in channel capacity between measurable performance and predictions from the most advanced simulation tools remains, which attracts significant research in propagation modeling mostly focusing on pure optical effects [23] but lately also considering weak acousto-optic interactions [24, 25]. Here we discuss NL depolarization of light (NLDP), a recently observed non-Manakovian transmission phenomenon as a potential candidate to narrow down this discrepancy [26]. Unpolarized optical noise rapidly changes via the fiber’s Kerr nonlinearity the SOP from a fully polarized cw light by inducing antisymmetric phase noise in both of its orthogonal polarization states. These fluctuations become resolvable with a new generation of high-speed polarimeters [27] and do not average out over wide noise bandwidths but grow with propagation distance.

Prior to applying the CNLS to NLDP, we introduce the phenomenon by providing some experimental evidence that will prompt the discussion of quantities needed for developing our algorithm. Back in 2017, two measurements were performed on an undersea cable that connects Brazil and Florida to demonstrate NLDP (Figure 1a) on a real communication link [26, 28]. A cw light launched from Brazil and looped back in Florida propagates via the opposite fiber link to its origin. There are about 200 optical amplifiers in each direction to compensate for fiber attenuation but cause weak NL propagation in the first few tens of kilometer fiber length directly after each repeater.

Figure 1.

(a) Comparative polarimetry on a commercial undersea cable system. A cw tone and unpolarized ASE propagate together from Brazil to the USA and back. The probe’s SOP speed is measured after propagation and compared with a local reference (not transmitted through the cable, dashed path); (b) SOP speed histogram recorded over 24 h showing environmental impact on SOP stability; (c) spectra of transmitted and received signals; (d) matched optical spectra of the probe and reference, taken with an OSA, look similar, but the corresponding SOP speed histograms (e) differ significantly and prove the existence of NLDP (noise boost explained in text).

Undersea cables on the ocean floor undergo environmental impact such as motion due to water currents, seismic vibrations, and temperature fluctuations, etc. These environmental factors lead to small SOP motions on the receive side, trackable with a low-speed polarimeter, and represented as velocities in histograms. We define the amount of the time derivative from the corresponding normalized Stokes vector as “SOP speed.” Typical SOP speeds range on the order of a few rad/s (Figure 1b). In this experiment, a cw light (ECL, External Cavity Laser) also referred to as probe (or signal) is launched together with unpolarized Amplified Spontaneous Emission (ASE, loading) into the cable input. The ASE spectrum covers the whole repeater bandwidth (∼4.5 THz) except for a narrow central gap with ∼100 GHz width where the probe resides (Figure 1c). Without loading, the probe would pick up almost the entire repeater output power resulting in strong NL propagation, mainly self phase modulation (SPM) [29], which can produce unstable SOPs. With loading, the probe contributes less than 1% of the total repeater output power and NL effects such as Brillouin scattering [30] are avoided. The electrical detection bandwidth of the low-speed polarimeter (a few 10s of Hz) averages out all fast SOP motions; while increasing its bandwidth (a few 10s of MHz) produces much wider SOP speed histograms under identical experimental conditions (Figure 1d).

It would be incorrect to attribute the entire width of this histogram solely to NL interactions between the loading and the probe, as any noise during the detection corrupts the polarimeter, leading to artifacts in SOP speeds. In our experiment, signal-ASE beat noise is the main contributor and [31] biases the SOP detection. Every repeater adds small amounts of ASE that raise the noise floor within the spectral gap. This ASE then mixes with the probe on the four photodiodes of the polarimeter. Even a hypothetical constant SOP of the probe on the receive side will appear in a histogram with nonzero width due to the omnipresent signal-ASE beat noise. We refer to this noise-induced SOP speed (NISS) as an artifact since improving the optical signal-to-noise ratio (OSNR) at the polarimeter input would reduce the width of its histogram. However, probe power constraints, inevitable added repeater noise during transmission, and practical limitations on tighter filtering yield OSNR levels, which result in artificially broadened SOP speed histograms that partially obscure NLDP.

Instead, we demonstrate NLDP-induced fast SOP changes by means of a comparison technique. We contrast the SOP speed histograms from the transmitted probe with one from a reference signal that possesses equal power, equal ASE, and equal OSNR but bypasses the undersea cable (reference, btb). This reference is obtained by superimposing the transmitter signal with the noise output of the cable in Brazil and launched btb into the receiver (ECL switched to reference path, Figure 1a). Within the resolution capabilities (∼2 GHz) of our optical spectrum analyzer (OSA), both the transmitted probe and the reference spectra are identical (Figure 1d). However, the SOP speed histogram for the probe is significantly wider compared with the histogram for the reference (Figure 1e) as NL interactions fluctuated the probe’s SOP. Even when the optical noise floor in the btb experiment is subject to added ASE (∼5 dB, Figure 1d), the corresponding histogram does not expand to the same width of the probe’s plot (Figure 1e).

It is reasonable to assume a scaling of the NLDP magnitude with the probe’s transmission length and the loading power, as both determine the strength of the Kerr nonlinearity. Such parametric NLDP studies become feasible with a lab test bed based on a recirculating fiber loop (RFL, Figure 2a) [32]. In RFLs, a fast optical switch allows to launch a signal into a fiber link whose output is coupled back into its input. After a programmable number of round trips inside the loop, the signal is released via the switch and analyzed.

Figure 2.

Characterizing NLDP in a lab setup using a recirculating fiber loop: (a) block diagram of an RFL design and its fiber specifics. (b) Setup for experimental comparative polarimetry embedding an RFL. The RFL allows economically emulating differently long transmission paths, which impacts the magnitude of NLDP. The propagation time is determined by a central clock that synchronizes the gates of the RFL, the high-speed polarimeter, and the OSA. (c) Visualization of a typical NLDP dependence on transmission distances; the histograms widen for the probe and reference with increased transmission distances. Their difference σ²_NLDP shows the NLDP magnitude and grows nearly linearly with the transmission distance. (d) Repeater output power dependence of NLDP at 10 mm transmission. The power-dependent Kerr effect mediates the cross-phase modulation between loading and the probe. Lowering the repeater output power reduces the width of the probe’s histograms but deforms the spectral shape of the loading as the gain shape of the individual EDFA changes.

Our RFL comprises 11 transmission spans each with an output power-adjustable repeater followed by SSMF (specs see Figure 2a). The transmitter (Figure 2b) inserts the probe via a variable optical attenuator (VOA) in path I into the RFL (VOA_I open, VOA_II blocked). After leaving the RFL, a narrow bandwidth amplified filter cascade selects the probe from the wide band loading prior to its detection by the polarimeter. To acquire the reference for determining NISS (btb measurement), VOA_I blocks the access to path I for the ECL while the reference reaches the receiver via path II. In both scenarios, loading enters the RFL via path I. We fine-tune VOA_II to match the spectra from probe and reference within 0.1 dB across a bandwidth of approx. 50 GHz. A central clock synchronizes the gates from loop switch, kicker (loop synchronous polarization controller), OSA, and polarimeter. It periodically activates the latter’s data acquisition at 100 MS/s typically for a ∼3 ms long interval per transmission. A spectrally programmable ASE source delivers unpolarized loading (Degree of Polarization (DOP) <1%) with flat spectrum (constant power density) across the repeater gain bandwidth of about 4.5 THz except for a narrow gap (100 GHz wide, ∼50 dB ASE suppression) centered at 193.9 THz to contain the probe. Here we chose a flat launch spectrum for the loading and a flat transfer function for the RFL, adjusted via its GEF (Gain Equalizing Filter, Figure 2a), to simplify a mathematical treatment of NLDP [26]. In reality, every system component such as repeaters, fibers, etc., possesses small amounts of PDL [33] and PMD. Due to the large number and the random orientation of their optical axes in Stokes space, the results follow a “smooth distribution.” But in an RFL, a signal passes the same components several times, which causes filtering and different PDL and PMD statistics. To counteract this effect and emulate scenarios such as those in real systems, a triggered loop synchronous polarization controller [34] (kicker) randomizes PMD and PDL inside the loop. After the signal has traversed the device and before its return, the kicker is set to a different mode within a few nanoseconds to randomize the output SOP. The kicker does not directly impact the SOP speed as it is stationary when the probe is passing through, but it emulates more realistic propagation paths.

By programming the loop’s timing gate, SOP speed histograms are determined (Figure 2c) for propagation lengths of approx. 1023 km, 3069 km, 5115 km, 10,230 km, and 20,460 km further referred to as 1 Mm, 3 Mm, …, 20 Mm transmissions and correspond to 1, 3, 5, 10, and 20 loop circulations, respectively. For longer transmissions, the corresponding btb measurements show wider NISS histograms as the received OSNRs decay. The probe’s SOP speed histograms broaden with propagation length due to NISS and NLDP. We hypothesize both as statistically independent processes and visualize NLDP by subtracting the reference’s NISS from the probe’s SOP speed variance (Figure 2c, inset). This quantity monotonically increases, indicating a growing NLDP with transmission distance. A linear fit reasonably resembles the measured NLDP variances as a function of the propagation length. Each recording consists of about 100 MSamples and was repeated multiple times to verify stable measurement readings. We found experimentally for the short-term reproducibility of all measured SOP speed variances a relative error smaller than 1%.

Without NL interaction with the loading, the probe would linearly propagate. Altering the nominal probe power of ∼−5.2 dBm at the repeater output by ±3 dB has shown an insignificant dependence of σ²_NLDP on SPM-based effects at 10 Mm transmission distance. Hence, Brillouin scattering or modulation instability (known to be the weakest NL fiber process) can be ruled out as origin for NLDP.

The repeater output power inside our RFL can be to some degree controlled without strongly tilting and distorting the spectral shape of the loading. At maximum repeater output power (20.9 dBm into transmission fiber) and at 1 and 2 dB down, we record SOP speed histograms over 10 Mm transmission distance following the same methodology. As three different phenomena need to be onsidered, NLDP, receive OSNR degradation, and spectral distortions of the loading, the data analysis becomes more complex. Increasing widths of the reference’s NISS histograms with decreasing repeater power (Figure 2d) can be explained by receive OSNR reduction. To a good approximation, the ASE that an optical amplifier adds to a boosted signal is proportional to the repeater’s gain. Since our repeaters operate in constant gain mode (as the span loss stays constant the repeater gain must do the same), altering the output power does not affect their ASE contributions. Consequently, lower repeater output power means smaller probe power and degraded OSNR after transmission (at constant ECL power), which enhances NISS. Remarkably, after transmission, and despite improved receive OSNR, the corresponding SOP speed histograms trend oppositely and widen with stronger repeater power. This broadening stems from NLDP, which even exceeds a theoretically expected narrowing of the histograms when only the enhanced receive OSNRs are considered as in case of the reference. Other than in the NLPD versus transmission distance—study, the received spectra diverge from an ideal boxcar shape and deform at different test conditions (Figure 2d, inset). While around 194.4 THz, the spectral density remains about constant, it drops disproportionately toward the spectral edges with reduced repeater power, making a quantitative analysis more challenging. But the noise floor surrounding the probe (zoom-in Figure 2d) remains relatively constant, which supports our previous OSNR argument.

As discussed in Section 4.7, NLDP is formed on long-range nonlinear optical interactions that become observable in the spectral domain. Distortions mainly generated within the first section of a transmission fiber connected to a repeater output substantially interfere among each other when propagating along a multi span link. Qualitatively spoken, the spectral features of the probe’s Stokes vector depend on the transmission distance—more precisely, the further the signal propagates, the leaner the Stokes vector spectrum becomes. In our study, we define as the probe’s spectrum the Fourier transform of the autocorrelation of its normalized, three-dimensional Stokes vector.

In our experimental verification, we utilize the same RFL as in the aforementioned setup. But instead of using a polarimeter, a polarization scrambling interferometer connected to a photodetector followed by an RF spectrum analyzer is deployed (Figure 3a). It can be shown that averaging many RF spectra yields the Stokes vector spectrum [35]. For emulated transmission distances between 1 and 20 Mm spectral width (FWHM), between 16.1 and 3.3 MHz, respectively, were observed. This demonstrates the spectral contraction over propagation length (Figure 3b).

Figure 3.

(a) Recirculating fiber loop (RFL) test bed with a polarization scrambling interferometer (PSI). A probe (ECL), embedded by unpolarized flat ASE, propagates through the RFL and the PSI records the power density spectrum of its stokes vector. (b) Stokes vector spectra for 1–20 mm transmissions and corresponding spectra for references showing flat noise floors. (c) Typical launch spectrum and receive spectra at the RFL output.

An analogy to Fabry-Perot etalons [36] whose filter width narrows with the number of interfering rays inside its cavity maybe helpful for understanding the spectral shaping of the Stokes vector spectrum with transmission distance. While a higher reflectivity of the etalon’s mirrors produces more interfering rays, more link spans lead to more coherently superimposing distortions in the receiver plane that cause spectral compression.

4. Applying the coupled NL Schrödinger equations to NLDP

Techniques for solving the CNLS for a short piece of fiber with constant linear birefringence have been extensively discussed in the literature [37]. SSMF can be envisioned as a concatenation of many fiber pieces with linear birefringence. The length of each piece and the orientation of its optical axes follow known statistics. Our goal is to calculate closed-form solutions for quantities that describe the underlying statistics of the reported NLDP effects based on the aforementioned fiber models. For the sake of simplicity, we take advantage of certain experimental conditions that justify a significantly reduced formalism. The conceptional simplifications relate to NL propagation, the repeater functionality, and the fiber model.

We assume NL signal propagation but consider its impact as relatively weak, which effectively addresses the operational range of today’s telecom systems. This allows us to model NL distortions as first-order perturbations that propagate linearly through the path once they have been generated.

Every repeater adds small amounts of ASE to the loading. While this extra noise does not significantly change the spectral shape of the loading, it slightly impacts its temporal correlation features, which will be ignored.

We will partially diverge from the picture of discrete and concatenated waveplates that form a fiber. For SSMF, it is appropriate to imagine smooth transitions between the single waveplates. However, we will consider a discrete waveplate model to discuss local SOP rotations and apply the theory of PMD statistics to cope with long-distance SOP correlations. These assumptions should not impact our main conclusions.

4.1 Sorting the Kerr nonlinearity in even and odd operators

For simplicity, we consider a weak cw field $a_{x\left(y\right)}$ (probe) residing in a narrow spectral gap of a surrounding, fully unpolarized and co-propagating strong ASE field $A_{x\left(y\right)}$ (loading) with a boxcar-shaped spectrum (Figure 4a). As stated, significantly less than 1% of the total signal power stems from the probe. The known CNLS [38] for propagation in z direction within a waveplate, given here in complete form, can approximate this scenario well:

Figure 4.

(a) Spectral grid for boxcar-shaped ASE spectrum with a probe residing in its center gap. (b) Single and multiple span link design for NLDP study.

$$\frac{\partial A_x^{\mathrm{\Sigma }}}{\partial z}+{\beta }_{1x}\frac{\partial A_x^{\mathrm{\Sigma }}}{\partial t}+\frac{j{\beta }_2}{2}\frac{{\partial }^2{A}_x^{\mathrm{\Sigma }}}{\partial t^2}+\frac{\alpha }{2}{A}_x^{\mathrm{\Sigma }}=\mathit{j\gamma }\left({\left|A_x^{\mathrm{\Sigma }}\right|}^2+\frac{2}{3}{\left|A_y^{\mathrm{\Sigma }}\right|}^2\right)A_x^{\mathrm{\Sigma }}+j\mathrm{\gamma }\frac{1}{3}{A}_x^{\mathrm{\Sigma }\ast }{A}_y^{{\mathrm{\Sigma }}^2}{e}^{-j2\mathrm{\Delta \beta }z}$$E7

$$\frac{\partial A_y^{\mathrm{\Sigma }}}{\partial z}+{\beta }_{1y}\frac{\partial A_y^{\mathrm{\Sigma }}}{\partial t}+\frac{j{\beta }_2}{2}\frac{{\partial }^2{A}_y^{\mathrm{\Sigma }}}{\partial t^2}+\frac{\alpha }{2}{A}_y^{\mathrm{\Sigma }}=\mathit{j\gamma }\left({\left|A_y^{\mathrm{\Sigma }}\right|}^2+\frac{2}{3}{\left|A_x^{\mathrm{\Sigma }}\right|}^2\right)A_y^{\mathrm{\Sigma }}+j\mathrm{\gamma }\frac{1}{3}{A}_y^{\mathrm{\Sigma }\ast }{A}_x^{{\mathrm{\Sigma }}^2}{e}^{+j2\mathrm{\Delta \beta }z}$$E8

$$\text{with}{A}_{x\left(y\right)}^{\mathrm{\Sigma }}=A_{x\left(y\right)}+a_{x\left(y\right)};\left|a_{x\left(y\right)}\right|\ll \left|A_{x\left(y\right)}\right|,$$E9

where the wavelength-independent α, β_1x(y)⁻¹, β₂, and γ denote the attenuation, polarization-dependent group velocities, dispersion coefficient, and the Kerr nonlinearity of a waveplate, respectively. As previously mentioned, the high modal birefringence Δβ of regular SMF induces fast oscillating of the terms at the far-right side, which leads to ineffective NL interference. This has been discussed in great detail for the derivation of the Manakov-PMD equation [8]; however, it is negligible in our analysis and thus left out in the following.

We assume that the loading modulates the probe, but the probe has no impact on the loading. In the limit of a negligibly small probe power, such interaction can be justified and simplifies the Kerr nonlinearity. Two sets of equation pairs follow, for the two cases $a_{x\left(y\right)}=0$ and $1\gg \left|a_{x\left(y\right)}\right|\ne 0$. Subtracting both sets, neglecting second-order terms, and considering coupling conditions yield an equation set that describes the motion of the probe induced by the loading

$$\frac{\partial a_x}{\partial z}+{\mathrm{\beta }}_{1x}\frac{\partial a_x}{\partial t}+\frac{{j\mathrm{\beta }}_{2x}}{2}\frac{{\partial }^2{a}_x}{\partial t^2}+\frac{\mathrm{\alpha }}{2}{a}_x=j\mathrm{\gamma }\left(2{\left|A_x\right|}^2+\frac{2}{3}{\left|A_y\right|}^2\right)a_x$$E10

$$\frac{\partial a_y}{\partial z}+{\mathrm{\beta }}_{1y}\frac{\partial a_y}{\partial t}+\frac{{j\mathrm{\beta }}_{2y}}{2}\frac{{\partial }^2{a}_y}{\partial t^2}+\frac{\mathrm{\alpha }}{2}{a}_y=j\mathrm{\gamma }\left(2{\left|A_y\right|}^2+\frac{2}{3}{\left|A_x\right|}^2\right)a_y.$$E11

The remaining stochastic perturbation on the right side is decomposed into a symmetric and an antisymmetric term (Eqs. (12) and (13)) with respect to the loading’s field

$$\frac{\partial a_x}{\partial z}+{\mathrm{\beta }}_{1x}\frac{\partial a_x}{\partial t}+\frac{{j\mathrm{\beta }}_{2x}}{2}\frac{{\partial }^2{a}_x}{\partial t^2}+\frac{\mathrm{\alpha }}{2}{a}_x=j\mathrm{\gamma }\left[\frac{4}{3}\left({\left|A_x\right|}^2+{\left|A_y\right|}^2\right)a_x+\frac{2}{3}\left({\left|A_x\right|}^2-{\left|A_y\right|}^2\right)a_x\right]$$E12

$$\frac{\partial a_y}{\partial z}+{\mathrm{\beta }}_{1y}\frac{\partial a_y}{\partial t}+\frac{{j\mathrm{\beta }}_{2y}}{2}\frac{{\partial }^2{a}_y}{\partial t^2}+\frac{\mathrm{\alpha }}{2}{a}_y=j\mathrm{\gamma }\left[\frac{4}{3}\left({\left|A_y\right|}^2+{\left|A_x\right|}^2\right)a_y+\frac{2}{3}\left({\left|A_y\right|}^2-{\left|A_x\right|}^2\right)a_y\right].$$E13

While both describe weak and independently treatable NL interactions in a first-order perturbation calculus, the polarization-dependent sign of the antisymmetric perturbation (right sides of Eqs. (12) and (13) root-causes opposite phase noises in both principal axes that manifests experimentally as NLDP. Identical phase changes in both polarizations as produced by the symmetric perturbation do not alter the probe’s SOP but lead to NL phase noise. In the following, we discuss solutions for the pair

$$\frac{\partial a_x}{\partial z}+{\mathrm{\beta }}_1\frac{\partial a_x}{\partial t}+\frac{{j\mathrm{\beta }}_2}{2}\frac{{\partial }^2{a}_x}{\partial t^2}+\frac{\mathrm{\alpha }}{2}{a}_x=j\mathrm{\gamma }\frac{2}{3}\left({\left|A_x\right|}^2-{\left|A_y\right|}^2\right)a_x$$E14

$$\frac{\partial a_y}{\partial z}+{\mathrm{\beta }}_1\frac{\partial a_y}{\partial t}+\frac{{j\mathrm{\beta }}_2}{2}\frac{{\partial }^2{a}_y}{\partial t^2}+\frac{\mathrm{\alpha }}{2}{a}_y=j\mathrm{\gamma }\frac{2}{3}\left({\left|A_y\right|}^2-{\left|A_x\right|}^2\right)a_y.$$E15

The impact of birefringence has been disregarded in Eqs. (14) and (15) by replacing β_1x(y) with a polarization-independent group velocity β₁⁻¹ but will be revisited at a later stage. That is, Eqs. (14) and (15) captures nonlinear SOP changes along a waveplate but does not include the much larger SOP changes caused by linear birefringence. The linear SOP changes are static and do not contribute to the measured SOP speeds. In practical applications, the right side represents a very fast fluctuating term as the noise components $A_x$and $A_y$ stem from stochastically independent processes with bandwidths in the range of several THz. Since we assume unpolarized loading, the quantity $\left({\left|A_{x\left(y\right)}\right|}^2-{\left|A_{y\left(x\right)}\right|}^2\right)$ is on average zero. For our algorithm, we represent it on an evenly spaced frequency grid, which yields experimentally observable quantities to describe some statistical features of the probe $a_{x\left(y\right)}$.

4.2 Asymmetric phase noise in first-order approximation

The right side of Eqs. (14) and (15) weakly perturbs the probe by adding a first-order correction term a_1x(y)(z,t) to it. We write the perturbation by means of the undistorted fields of the probe a_x(y)(z) and the loading A_x(y)(z,t). The latter’s components with amplitudes A_x(y) form a comb Eqs. (16) and (17) on an evenly spaced grid with an infinitesimally small angular frequency pitch ω (Figure 4a). For a single span system, the probe a_x(y)(z,t) in first-order development and the ASE field read

$$A_{0x}\left(z,t\right)=\sum _{m=-N}^N{A}_x^m{e}^{j\left(k_{m}z-\mathit{m\omega t}\right)}{e}^{-\frac{\mathrm{\alpha }}{2}z}$$E16

$$A_{0y}\left(z,t\right)=\sum _{m=-N}^N{A}_y^m{e}^{j\left(k_{m}z-\mathit{m\omega t}\right)}{e}^{-\frac{\mathrm{\alpha }}{2}z}$$E17

$$a_x\left(z,t\right)\approx \left(a_{0x}+a_{1x}\left(z,t\right)\right)e^{-\frac{\mathrm{\alpha }}{2}z}=\left(a_{0x}+{\sum }_l{a}_{1x}^l\left(z\right)e^{-\mathit{\text{jlωt}}}\right)e^{-\frac{\mathrm{\alpha }}{2}z}$$E18

$$a_y\left(z,t\right)\approx \left(a_{0y}+a_{1y}\left(z,t\right)\right)e^{-\frac{\mathrm{\alpha }}{2}z}=\left(a_{0y}+{\sum }_l{a}_{1y}^l\left(z\right)e^{-\mathit{\text{jlωt}}}\right)e^{-\frac{\mathrm{\alpha }}{2}z}$$E19

$$\text{with}{A}_{x\left(y\right)}^m=0\text{for}\left|m\right|<N_u,\left|m\right|>N,z\ge 0\text{and}{A}_{0x}\left(z,t\right)=a_x\left(z,t\right)=0\text{for}z<0$$E20

where $k_m={\beta }_1\mathit{m\omega }+\frac{{\beta }_2}{2}{\left(\mathit{m\omega }\right)}^2$ stands for the propagation constant of a component at ‘$\mathit{m\omega }$’. We synthesize the Kerr nonlinearity in Eqs. (14) and (15) as a sum to address the impact of low-frequency beat noises among its terms. Due to phase matching conditions, this noise alone can efficiently interact with the probe and is used to redefine the perturbation term in Eqs. (14) and (15) as

$${\left|A_x\right|}^2-{\left|A_y\right|}^2\stackrel{\scriptscriptstyle\mathrm{def}}{=}\sum _{m,l}{\mathcal{A}}_x^{l,m}{e}^{-\mathit{\alpha z}}{e}^{j\left(\left({\beta }_1\mathit{l\omega }+{\beta }_2\mathit{lm}{\omega }^2\right)z-\mathit{l\omega t}\right)}$$E21

$${\left|A_y\right|}^2-{\left|A_x\right|}^2\stackrel{\scriptscriptstyle\mathrm{def}}{=}\sum _{m,l}{\mathcal{A}}_y^{l,m}{e}^{-\mathit{\alpha z}}{e}^{j\left(\left({\beta }_1\mathit{l\omega }+{\beta }_2\mathit{lm}{\omega }^2\right)z-\mathit{l\omega t}\right)}$$E22

$$\mathit{\text{with}}{\mathcal{A}}_{x\left(y\right)}^{l,m}=A_{x\left(y\right)}^{m+\frac{l}{2}}{A_{x\left(y\right)}^{m-\frac{l}{2}}}^{\ast }-A_{y\left(x\right)}^{m+\frac{l}{2}}{A_{y\left(x\right)}^{m-\frac{l}{2}}}^{\ast },\left|l\right|\ll N_u.$$E23

Here $\mathit{m\omega }$ stands for the angular frequency spacing between the two beating ASE tones at $\omega \left(m\pm \frac{l}{2}\right)$ and the probe, which is typically in the THz range. Due to coupling inefficiency, any beating among the two noise tones can be neglected when its frequency $\mathit{l\omega }$ resides beyond a few tens of MHz. Eqs. (21) and (22) hold separately for every $l$ and $m$, and its solution for a single span system can be written by means of a Green’s function

$$a_{1x\left(y\right)}^{l,m}\left(z,t\right)=j\frac{2}{3}\gamma \underset{0}{\overset{z}{\int }}{e}^{j{k}_l\left(z-z_i\right)-\alpha z_i}{\mathcal{A}}_{x\left(y\right)}^{l,m}{e}^{j\left(\left({\beta }_1\mathit{l\omega }+{\beta }_2\mathit{lm}{\omega }^2\right)z_i-\mathit{l\omega t}\right)}{a}_{0x\left(y\right)}d{z}_i$$E24

According to this Ansatz, the NL distortions are generated as a kind of “wave packet” in fiber sections of incremental length $d{z}_i$, propagate thereafter linearly through the span, and coherently superimpose in the receiver plane. We portion out the integral of Eq. (24) into a sum of infinitesimal short waveplates $d{z}_i$ (Eq. (25)) to analyze NL interactions in the presence of Polarization Mode Dispersion (PMD), which originates from fiber birefringence. At the span output at $L_0$ holds

$${\delta }_{1x\left(y\right)}^{l,m}\left(L_0,t\right)=\frac{a_{1x\left(y\right)}^{l,m}\left(L_0,t\right)}{a_{0x\left(y\right)}}=j\frac{2}{3}\gamma e^{j\left({\beta }_1\mathit{l\omega }{L}_0-\mathit{l\omega t}\right)}\sum _i{\mathcal{A}}_{x\left(y\right)}^{m,l}{e}^{-\alpha z_i}{e}^{j{\beta }_2\mathit{lm}{\omega }^2{z}_i}{e}^{j\frac{1}{2}{\beta }_2{l}^2{\omega }^2\left(L_0-z_i\right)}d{z}_i.$$E25

Alternating the sign of$l$ conjugates its right side, except for its last and typically negligible small exponent ($a_{1x\left(y\right)}^{l,m}\left(L_0,t\right)\approx a_{1x\left(y\right)}^{-l,m}{\left(L_0,t\right)}^{\ast })$. Hence, pairing contributions at $\pm l$ results into a correction with a 90^○ phase offset relative to the undistorted probe. Therefore, all pairs of NL distortions stemming from single waveplates generate pure phase oscillations in the receiver plane $\epsilon \left(L_0,t\right)=\left({\delta }_{1x\left(y\right)}^{l,m}\left(L_0,t\right)+{\delta }_{1x\left(y\right)}^{-l,m}\left(L_0,t\right)\right)$ at an angular frequency $\mathit{l\omega }$. As ${\mathcal{A}}_{x\left(y\right)}^{m,l}=-{\mathcal{A}}_{y\left(x\right)}^{m,l}$ holds, the oscillations in both orthogonal polarizations are 180^○ out of phase and cause SOP fluctuations. We define for later purposes a temporal autocorrelation as

$${\phi }_{i,k}^{l,m}\left(\tau \right)=⟨{\delta }_{1x,i}^{l,m}\left(L_0,t+\tau \right){\delta }_{1x,k}^{l,m}{\left(L_0,t\right)}^{\ast }⟩+⟨{\delta }_{1y,i}^{l,m}\left(L_0,t+\tau \right){\delta }_{1y,k}^{l,m}{\left(L_0,t\right)}^{\ast }⟩,$$E26

where the indices $i,k$ denote different waveplates and $⟨⟩$ denotes the averaging over time and fields, which involves reestablishing the birefringent fiber features in our model as detailed below.

4.3 Correlations of asymmetric phase noise in birefringent fiber

For deriving correlations between phase noises generated at different propagation distances of the probe, we focus on a single span system and extend the results to a multiple span link. Modeling a birefringent fiber as a concatenation of waveplates uses Eq. (25) to determine the correlation among incremental distortions stemming from two different plates. Our x(y)-coordinate system is congruent with the fast (slow) axis of a waveplate, i.e., it rotates and follows the plates’ orientations along the propagation path. Our model incorporates birefringence, originating from, e.g., axis-specific group velocities (β_1x ≠ β_1y), by means of a Jones matrix that transforms the input SOPs from the probe and the noise components $A_{x\left(y\right)}^m$ when traversing a waveplate. A Jones matrix ${\overline{R}}_i$ of a waveplate shall be given by a unitary matrix

$${\overline{R}}_i=\left[\begin{array}{cc}{R}_{11}^i& R_{12}^i\\ -R_{12}^{i\ast }& R_{11}^{i\ast }\end{array}\right]\text{with}{\left|R_{11}^i\right|}^2+{\left|R_{12}^i\right|}^2=1.$$E27

With$a_{0x\left(y\right)}≔\overrightarrow{a_0}=\left(\begin{array}{c}{a}_{0x}\\ a_{0y}\end{array}\right)$ and $\overline{I}=\left[\begin{array}{cc}1& 0\\ 0& -1\end{array}\right]$, an NL distortion generated in waveplate $k$ with length $d{z}_k$ concisely reads at the span output

$$\overrightarrow{a_{1,k}^{m,l}}\left(L_0,t\right)=\left(\begin{array}{c}{a}_{1x,k}^{m,l}\left(L_0,t\right)\\ a_{1y,k}^{m,l}\left(L_0,t\right)\end{array}\right)=j\frac{2}{3}{\mathrm{\gamma }e}^{j\left({\beta }_1\mathit{l\omega }{L}_0-\mathit{l\omega t}\right)}{e}^{j{\beta }_2\mathit{lm}{\omega }^2{z}_k}d{z}_k\times \left(\prod _{k=i}^{M_{L_0}}{\overline{R}}_i\right){\mathcal{A}}_{x,k}^{m,l}\overline{I}{e}^{-\alpha z_k}\left(\prod _{i=1}^{k-1}{\overline{R}}_i\right)\overrightarrow{a_0},$$E28

where ${\mathcal{A}}_{x,k}^{m,l}$ und $M_{L_0}$stand for the undepleted noise components inside of waveplate $k$constituted by Eq. (23), and the total number of waveplates, respectively. The correlation between contributions from two consecutive waveplates follows as

$$⟨\overrightarrow{a_{1,k+1}^{m,l}}\left(L_0\right)|\overrightarrow{a_{1,k}^{m,l}}\left(L_0\right)⟩=\frac{4}{9}{\mathrm{\gamma }}^2{e}^{j{\beta }_2\mathit{lm}{\omega }^2\left(z_k-z_{k+1}\right)}{e}^{-\alpha \left(z_k+z_{k+1}\right)}d{z}_{k}d{z}_{k+1}\times ⟨{\mathcal{A}}_{x,k+1}^{m,l}\overline{I}{\overline{R}}_k\overrightarrow{a_{0,k}}|{\mathcal{A}}_{x,k}^{m,l}{\overline{R}}_k\overline{I}\overrightarrow{a_{0,k}}⟩,$$E29

where $\overrightarrow{a_{0,k}}$ and $⟨|⟩$ represent the undistorted and undepleted probe field within waveplate $k$ and the field-averaged scalar product, respectively. We restore PMD in our fiber model by using matrices ${\overline{R}}_i,{{\overline{R}}_i}^{\prime }$ obeying Eq. (27) and transforming wavelength-dependent the probe and ASE fields:

$$⟨{\mathcal{A}}_{x,k+1}^{m,l}\overline{I}{\overline{R}}_k\overrightarrow{a_{0,k}}|{\overline{R}}_k{\mathcal{A}}_{x,k}^{m,l}\overline{I}\overrightarrow{a_{0,k}}⟩=⟨\left({\left|R_{11}^k\right|}^2-{\left|R_{12}^k\right|}^2\right)\left({\left|R_{11}^{k\prime }\right|}^2-{\left|R_{12}^{k\prime }\right|}^2\right)⟩{\mathcal{H}}^m{\left|\overrightarrow{a_{0,k}}\right|}^2$$E30

$$\text{with}{\mathcal{H}}^m=⟨{\left|A_x^{m"}\right|}^2{\left|A_x^{m"}\right|}^2+{\left|A_y^m\right|}^2{\left|A_y^{m"}\right|}^2⟩,$$E31

where double primes indicate statistically independent noises. Since $m\gg l$ holds, we can evaluate Eq. (30) just at $m$ but must treat the original two noise components at frequencies $\left(m\pm \raisebox{1ex}{$l$}\!\left/ \!\raisebox{-1ex}{$2$}\right.\right)\omega$ as uncorrelated. The theory of PMD statistics [39] specifies the distance-dependent decorrelation from two SOPs of two co-propagating cw tones at different wavelengths, which provides a correlation between the matrix elements

$$⟨\left({\left|R_{11}^k\right|}^2-{\left|R_{12}^k\right|}^2\right)\left({\left|R_{11}^{k\prime }\right|}^2-{\left|R_{12}^{k\prime }\right|}^2\right)⟩=\frac{1}{3}{e}^{-\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2∆L},$$E32

where $m\omega$, ${\tau }_{\rho }$, and $∆L=\left|z_{k+1}-z_k\right|$ are the angular frequency spacing between the probe and the two noise components, the mean fiber DGD per √length (Differential Group Delay), and the propagation distance, respectively. PMD effects across frequency intervals of size ∼$l\omega$ such as the spectral width of the received and distorted probe or spacings between two efficiently beating noise components are negligibly small. Hence, the sum of all incremental phase oscillations correlates as defined by Eq. (26) at the span output like

$${\phi }^{l,m}\left(\tau \right)=\frac{4}{27}{\gamma }^2{e}^{-\mathit{\text{jlωτ}}}{\mathcal{H}}^m\sum _{i,k}{e}^{j{\beta }_2\mathit{lm}{\omega }^2\left(z_i-z_k\right)-\alpha \left(z_i+z_k\right)}{e}^{-\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2\mid z_i-z_k\mid }d{z}_{i}d{z}_k.$$E33

The above outlined calculus assumes two consecutive birefringent waveplates $k,k+1$. But it holds for any pair of further spaced waveplates, indexed $i,k$ with $i\ne k+1$, as well. Since a matrix product $\overline{R}={\overline{R}}_i\dots {\overline{R}}_{k+2}{\overline{R}}_{k+1}$ of intermediately located waveplates can be expressed by a single unitary matrix that fulfills Eq. (27), the conclusion from Eq. (30) will equally hold and leads to Eq. (33).

4.4 Fiber PMD constitutes NLDP

System PMD imposes cutoff conditions via the Gaussian for the number of interacting waveplates addressed by the double sum of Eq. (33). Without this constrain, the sum tends to zero as its complex exponential function causes averaging for sufficiently small $\alpha$. For a single span system with $L_0\gg {\alpha }^{-1}$ (typically tens of kilometers) and relatively short waveplates $(\ll {\alpha }^{-1}$, typically tens of meters), we will replace the double sum with an integral.

Experimentally observed NLDP-caused SOP features such as scattering angles and speed as their time derivatives are detected after O/E conversion of the optical fields and conveniently reported in Stokes space (Section 3). To derive such quantities, we will confine the optical autocorrelation density by introducing electrical low-pass filtering, which represents the detection process, and then convert the result into Stokes space. In Jones space, the density of the optical autocorrelation (Eq. (33)) at $\mathit{l\omega }$ reads for sufficiently long propagation distances $\left(L_0\to \infty \right)$

$$\begin{array}{rc}{\phi }^{l,m}\left(\tau \right)& =\frac{4}{27}{\gamma }^2{\mathcal{H}}^m\frac{e^{-\mathit{\text{jlωτ}}}}{2\alpha }{\int }_{-\infty }^{+\infty }{e}^{-\left(\alpha +\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2\right)\left|z\right|}{e}^{j{\beta }_2\mathit{lm}{\omega }^{2}z}\mathit{dz}\\ & =\frac{4}{27}{\gamma }^2{\mathcal{H}}^m\frac{e^{-\mathit{\text{jlωτ}}}}{\alpha }\frac{\alpha +\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2}{{\left(\alpha +\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2\right)}^2+{\left({\beta }_2\mathit{lm}{\omega }^2\right)}^2}\end{array}$$E34

and its integration over $m,l^{\prime }$ yields the autocorrelation for the total phase noise.

So far, we have considered a single span system, i.e., one optical amplifier boosts the probe and loading prior of their launch into an as infinite long assumed transmission fiber. Today’s undersea communication cables can consist of a few hundred spans to bridge transpacific distances of up to about 15 Mm length. In contrast to terrestrial systems, they are strictly modularly designed, which eases our modeling. To keep the calculation effort at an introductory level, we further assume a large enough fiber attenuation$\alpha ,$ thus all NL propagation fades away far before the span end. Additionally, our amplification is assumed to be a noiseless process. For a cable with N_S spans, the autocorrelation then reads

$${\phi }_{N_s}^{l,m}\left(\tau \right)=\frac{4}{27}{\gamma }^2{e}^{-\mathit{\text{jlωτ}}}{\mathcal{H}}^m\times \sum _{p=0,q=0}^{N_s,N_s}\sum _{i,k}{e}^{j{\beta }_2\mathit{lm}{\omega }^2\left(z_i-z_k+L_0\left(p-q\right)\right)-\alpha \left(z_i+z_k\right)}{e}^{-\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2\left|z_i-z_k+L_0\left(p-q\right)\right|}d{z}_{i}d{z}_k$$E35

with $z_i+L_{0}p$ and $z_k+L_{0}q$ as the positions of two interacting waveplates.

For approximating the interleaved summations, we take advantage of the assumed large fiber attenuation, which forms cutoff conditions for the $z_i$ and $z_k$. When one of the two coordinates or both are large enough but still significantly smaller than $L_0$, the corresponding summands do not substantially contribute to the overall sum. Thus, a replacement of the amount in the exponents by

$$\mid z_i-z_k+L_0\left(p-q\right)\left|\to |z_i-z_k\right|+L_0\left|p-q\right|$$E36

can be justified and yields

$$\begin{array}{l}{\phi }_{N_s}^{l,m}\left(\tau \right)\cong \frac{4}{27}{\gamma }^2{e}^{-\mathit{\text{jlωτ}}}{\mathcal{H}}^m\times \sum _{i,k}{e}^{j{\beta }_2\mathit{lm}{\omega }^2\left(z_i-z_k\right)-\alpha \left(z_i+z_k\right)}{e}^{-\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2\left|z_i-z_k\right|}d{z}_{i}d{z}_k\\ \times \sum _{p=0,q=0}^{N_s,N_s}{e}^{j{\beta }_2\mathit{lm}{\omega }^2{L}_0\left(p-q\right)}{e}^{-\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2{L}_0\left|p-q\right|}\end{array}$$E37

$$\begin{array}{l}=\frac{4}{27}{\gamma }^2\frac{e^{-\mathit{\text{ilωτ}}}}{\alpha }{\mathcal{H}}^m\times \frac{\alpha +\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2}{{\left(\alpha +\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2\right)}^2+{\left({\beta }_2\mathit{lm}{\omega }^2\right)}^2}\\ \times N_s\sum _{n=-N_s}^{N_s}{e}^{j{\beta }_2\mathit{lm}{\omega }^2{L}_{0}n}{e}^{-\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2{L}_0\left|n\right|}\mathrm{\Lambda }\left(\frac{n}{N_s}\right),\end{array}$$E38

where $\mathrm{\Lambda }\left(.\right)$ denotes the triangular function. For the degree of accuracy, we follow in our modeling, it is sufficient to approximate it with an exponential function of even area and substitute both sums by integrals to obtain at least trends that show how NLDP qualitatively depends on system parameters.

$$\begin{array}{l}{\phi }_{N_s}^{l,m}\left(\tau \right)\approx \frac{4}{27}{\gamma }^2\frac{e^{-\mathit{\text{ilωτ}}}}{\alpha }{\mathcal{H}}^m\frac{\alpha +\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2}{{\left(\alpha +\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2\right)}^2+{\left({\beta }_2\mathit{lm}{\omega }^2\right)}^2}\\ \times N_s\sum _{n=-\infty }^{\infty }{e}^{j{\beta }_2\mathit{lm}{\omega }^2{L}_{0}n}{e}^{-\left(\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2{L}_0+\frac{2}{N_s}\right)\left|n\right|}.\end{array}$$E39

Further we approximate the sum with a Fourier integral and find for the density of the optical correlation

$${\phi }_{N_s}^{l,m}\left(\tau \right)\approx \frac{8}{27}{\gamma }^2{\mathcal{H}}^m{N}_s\frac{e^{-\mathit{\text{jlωτ}}}}{\alpha L_0}\frac{1}{\alpha +\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2}\frac{2}{\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2+\frac{2}{N_s{L}_0}}\frac{1}{1+{\left(\frac{{\beta }_2\mathit{lm}{\omega }^2}{\alpha +\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2}\right)}^2}\frac{1}{1+{\left(\frac{{\beta }_2\mathit{lm}{\omega }^2}{\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2+\frac{2}{N_s{L}_0}}\right)}^2}.$$E40

4.5 Detecting optical fields in the electrical domain

For the sake of simplicity, we consider two extreme cases for the electrical detection bandwidth ${\mathrm{\Omega }}_e$ of our polarimeter. In one scenario its bandwidth tends to infinity and in the other, it strongly filters the photocurrents. In both cases${\tau }_S^{-1}\gg {\mathrm{\Omega }}_e$ should hold for its sampling rate. In typical lab experiments, the electrical detection bandwidth is the only parameter that can be practically tuned over a larger range without the need for readjusting other model assumptions. For example, reducing the optical bandwidth of a system (∼4.5 THz) at launch or the launch power into a span usually distorts the assumed boxcar-shaped spectrum of the loading.

We model the low-pass characteristic of the polarimeter with a Lorentzian curve, thus the autocorrelations in the optical and electrical domain interrelate as

$${\phi }_{N_s}^{\mathit{\text{elec}}}\left(\tau \right)={\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }\frac{1}{1+{\left(\frac{\mathit{l\omega }}{{\mathrm{\Omega }}_e}\right)}^2}{\phi }_{N_s}^{l,m}\left(\tau \right)\mathit{\text{dωldωm}},$$E41

Provided sufficiently short sampling periods ${\tau }_S$, compared with the inverse of the autocorrelation’s spectral width², the exponential in Eq. (40)$\left(e^{-\mathit{\text{jlωτ}}}\approx 1-\mathit{\text{jlωτ}}-\frac{1}{2}{\left(\mathit{l\omega \tau }\right)}^2\right)$ can be approximated in $\tau$. Its linear term vanishes after the integration (Eq. (41)) due to symmetry aspects. The second order in ${\tau }_S$ of ${\phi }_{\left({\tau }_S\right)}^{\text{elec}}$ determines the variance of the stochastic SOP speed as shown below and reads

$${\phi }_{N_s}^{\text{elec}}\left(0\right)-{\phi }_{N_s}^{\text{elec}}\left({\tau }_S\right)={\int }_{-\infty }^{\infty }{\int }_{-\infty }^{\infty }\frac{1}{1+{\left(\frac{\mathit{l\omega }}{{\mathrm{\Omega }}_e}\right)}^2}\left({\phi }_{N_s}^{l,m}\left(0\right)-{\phi }_{N_s}^{l,m}\left({\tau }_s\right)\right)\mathit{\text{dωldωm}},$$E42

which tends for the case of wide bandwidth detection (${\mathrm{\Omega }}_e\to \infty$) and under the assumptions of typical system parameters as explained in Section 3 toward

$${\phi }_{N_s}^{\text{elec}}\left(0\right)-{\phi }_{N_s}^{\text{elec}}\left({\tau }_S\right)\approx \frac{\pi }{81}{\gamma }^2\frac{N_s}{\alpha L_0}{\beta }_2^{-3}{\tau }_S^2{\tau }_{\rho }^2{P}_{\mathrm{rep}}^2\frac{1}{{\left({\mathrm{\Omega }}_{\max }-{\mathrm{\Omega }}_{\min }\right)}^2}ln\left(\frac{{\mathrm{\Omega }}_{\max }}{{\mathrm{\Omega }}_{\min }}\right),$$E43

where $P_{\mathrm{rep}}$ represents the total launch power of the loading at the fiber input.

For the case of detection at a reduced electrical bandwidth (${\mathrm{\Omega }}_e\ll \frac{1}{3}{\mathrm{\Omega }}_{\min }\frac{{\tau }_{\rho }^2}{{\beta }_2}$), we find

$${\phi }_{N_s}^{\text{elec}}\left(0\right)-{\phi }_{N_s}^{\text{elec}}\left({\tau }_S\right)\approx \frac{\pi }{27}{\gamma }^2\frac{N_s}{\alpha L_0}\frac{{\Omega }_e^2}{{\mathrm{\Omega }}_{\max }^2}\frac{{\tau }_S^2}{{\beta }_2\alpha }{P}_{\mathrm{rep}}^{2}ln\left(\sqrt{\frac{3\alpha }{2}}\frac{1}{{\tau }_{\rho }{\mathrm{\Omega }}_{\min }}\right).$$E44

As Eq. (44) rapidly declines for small ${\mathrm{\Omega }}_e^2$, NLDP stays hidden in the study of environmentally driven SOP fluctuations in undersea cables (Figure 1b), which had been performed at a detection bandwidth below 100 Hz.

4.6 Transforming phase noise from Jones space to stokes space

So far, NLDP has been characterized in Jones space for which the CNLSs hold. But for convenience in general, quantities such as the SOP speed and the SOP scattering angle [40] are usually discussed and experimentally obtained in Stokes space. To transform the incremental field distortions, determined by means of the Jones calculus and the CNLS, into Stokes space, we represent the assumed normalized probe’s Jones vector $\overrightarrow{a_{\mathrm{norm}}}\left(N_s{L}_0,t\right)={\left[a_x\left(N_s{L}_0,t\right),a_y\left(N_s{L}_0,t\right)\right]}^T$ at the fiber output by

$$\overrightarrow{a_{\mathrm{norm}}}\left(N_s{L}_0,t\right)=\left[\begin{array}{c}{a}_x\left(N_s{L}_0,t\right)\\ a_y\left(N_s{L}_0,t\right)\end{array}\right]=\left[\begin{array}{c}cos\vartheta cos\theta +jsin\vartheta sin\theta \\ sin\vartheta cos\theta -jcos\vartheta sin\theta \end{array}\right]$$E45

with $\vartheta ={\vartheta }_0+\mathit{\delta \vartheta }\left(t\right),\left|2{\vartheta }_0\right|<\pi$and $\theta ={\theta }_0+\mathit{\delta \theta }\left(t\right),\left|4{\theta }_0\right|<\pi$where ${\vartheta }_0,{\theta }_0$ obey known distributions [41] to uniformly cover the Poincare sphere. For the SOP density $\Psi$ holds

$$\Psi \left(2{\vartheta }_0,2{\theta }_0\right)d2{\vartheta }_{0}d2{\theta }_0=\frac{1}{4\pi }\text{rect}\left(\frac{2{\vartheta }_0}{2\pi }\right)\text{rect}\left(\frac{2{\theta }_0}{\pi }\right)cos\left(2{\theta }_0\right)d2{\vartheta }_{0}d2{\theta }_0.$$E46

NLDP causes small temporal fluctuations $\mathit{\delta \vartheta }\left(t\right),\mathit{\delta \theta }\left(t\right)$ whose autocorrelations must equal $⟨\mathit{\delta \vartheta }\left(t+{\tau }_S\right)\mathit{\delta \vartheta }\left(t\right)⟩=⟨\mathit{\delta \theta }\left(t+{\tau }_S\right)\mathit{\delta \theta }\left(t\right)⟩={\phi }_J\left({\tau }_S\right)$ due to averaging caused by birefringence fluctuations along the long propagation path. The corresponding three-dimensional Stokes vector with unity length$\left(S_0=1\right)$ to Eq. (45) reads

$$\overrightarrow{S}\left(N_s{L}_0,t\right)={\left[\begin{array}{c}{\left|a_x\right|}^2-{\left|a_y\right|}^2\\ +2\mathfrak{R}ℯ\left\{a_x{a_y}^{\ast }\right\}\\ -2\mathfrak{Im}\left\{a_x{a_y}^{\ast }\right\}\end{array}\right]}_{\left(N_s{L}_0,t\right)}={\left[\begin{array}{c}cos\left(2\vartheta \right)cos\left(2\theta \right)\\ sin\left(2\vartheta \right)cos\left(2\theta \right)\\ -sin\left(2\theta \right)\end{array}\right]}_{\left(N_s{L}_0,t\right)}$$E47

and can be analyzed with respect to the impact of $\mathit{\delta \vartheta }\left(t\right)$ and $\mathit{\delta \theta }\left(t\right).$ One finds

$$⟨{\left|\overrightarrow{S}\left(t+{\tau }_S\right)-\overrightarrow{S}\left(t\right)\right|}^2⟩=\frac{40}{3}\left({\phi }_J\left(0\right)-{\phi }_J\left({\tau }_S\right)\right).$$E48

In Jones space we get

$$⟨{\left|\overrightarrow{a_{\mathrm{norm}}}\left(t+{\tau }_S\right)-\overrightarrow{a_{\mathrm{norm}}}\left(t\right)\right|}^2⟩=4\left({\phi }_J\left(0\right)-{\phi }_J\left({\tau }_S\right)\right).$$E49

Thus, combing both results add up to

$$⟨{\left|\overrightarrow{S}\left(t+{\tau }_S\right)-\overrightarrow{S}\left(t\right)\right|}^2⟩=\frac{10}{3}⟨{\left|\overrightarrow{a_{\mathrm{norm}}}\left(t+{\tau }_S\right)-\overrightarrow{a_{\mathrm{norm}}}\left(t\right)\right|}^2⟩$$E50

Identifying $⟨{\left|\overrightarrow{a_{\mathrm{norm}}}\left(t+{\tau }_S\right)-\overrightarrow{a_{\mathrm{norm}}}\left(t\right)\right|}^2⟩$=$2\left({\phi }_{N_s}^{\text{elec}}\left(0\right)-{\phi }_{N_s}^{\text{elec}}\left({\tau }_S\right)\right)$ yields the variance of the NLDP-induced SOP speed in Stokes space

$$\frac{⟨{\left|\overrightarrow{S}\left(t+{\tau }_S\right)-\overrightarrow{S}\left(t\right)\right|}^2⟩}{{\tau }_S^2}\approx \frac{2\pi }{81}{\gamma }^2\frac{N_s}{\alpha L_0}{\beta }_2^{-3}{\tau }_{\rho }^2{P}_{\mathrm{rep}}^2\frac{1}{{\left({\mathrm{\Omega }}_{\max }-{\mathrm{\Omega }}_{\min }\right)}^2}ln\left(\frac{{\mathrm{\Omega }}_{\max }}{{\mathrm{\Omega }}_{\min }}\right),$$E51

which depends on the PMD of the link and differentiates NLDP from NL polarization rotation (NLPR) [15]. NLPR is a phenomenon between two cw tones where one impacts the SOP of the other. However, its fundamental equations do not include fiber features such as PMD and chromatic dispersion.

4.7 NLDP-induced stokes vector spectrum of the probe

In analogy to the relationship between temporal and spectral features of an electrical signal in time and frequency domain, we define the spectrum of the probe’s Stokes vector by the Fourier transform of its temporal autocorrelation

$$⟨\overrightarrow{S}\left(t+{\tau }_S\right)\bullet \overrightarrow{S}\left(t\right)⟩=1-\frac{20}{3}{\phi }_J\left(0\right)+\frac{20}{3}{\phi }_J\left({\tau }_S\right)\approx 1+\frac{20}{3}{\phi }_J\left({\tau }_S\right)\circ -∎\Phi \left(f\right).$$E52

For our following derivations, we revisit Eq. (52). It does not include several approximations we have meanwhile applied in our calculus and can therefore be considered as a more accurate starting point. But we simplify the integration over the loading’s optical spectrum by neglecting its central gap and assume infinite integral bounds, i.e., ${\mathcal{H}}^m=\overline{{\mathcal{H}}^m}=\mathit{\text{const}}$ holds

$$\begin{array}{l}{\phi }_{N_s}^{l,m}\left(\tau \right)=\frac{4}{27}{\gamma }^2\frac{e^{-\mathit{\text{ilωτ}}}}{{\alpha }^2}\overline{{\mathcal{H}}^m}\frac{1}{1+\frac{1}{3\alpha }{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2}\frac{1}{1+{\left(\frac{{\beta }_2\mathit{lm}{\omega }^2}{\alpha +\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2}\right)}^2}\\ \times N_s\sum _{n=-N_s}^{N_s}{e}^{j{\beta }_2\mathit{lm}{\omega }^2{L}_{0}n}{e}^{-\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2{L}_0\left|n\right|}\mathrm{\Lambda }\left(\frac{n}{N_s}\right).\end{array}$$E53

For typical system parameters, Eq. (53) is dominated by the term $e^{-\frac{1}{3}{\left(\mathit{\omega m}\right)}^2{\tau }_{\rho }^2{L}_0\left|n\right|}$, which sets even for low fiber PMD (∼50 fs/√km) and small |n| > 0 a strong cutoff criterion for the integration. In this case the first two fractions on the right side will be ignored and the integration yields

$${\phi }_{N_s}^l\left(\tau \right)\sim e^{-\mathit{\text{jlωτ}}}\sum _{n=-N_s,n\ne 0}^{N_s}\frac{1}{\sqrt{\left|n\right|}}{e}^{-\frac{3{\beta }_2^2}{4{\tau }_{\rho }^2}{\left(\mathit{\omega l}\right)}^2{L}_0\left|n\right|}\mathrm{\Lambda }\left(\frac{n}{N_s}\right).$$E54

A narrowing of the probe’s Stokes vector spectrum with an increasing number of system spans can be qualitatively understood by defining its spectral width as

$$\frac{1}{{\phi }_{N_s}^0\left(0\right)}\int {\phi }_{N_s}^l\left(0\right)\mathit{d\omega l}\approx \frac{ln\left(N_s\right)-1}{\sqrt{N_s}-1}.$$E55

This quantity decreases for large enough and increasing $N_s$ . NL distortions generated along the transmission link interfere, which leads to spectral shaping of the Stokes vector spectrum. Such kind of spectral compression has been experimentally observed (Section 3). The more spans a link consists of, the more nonlinear distortions in form of wave packets superpose in the receiver plane and steepen the pedestal of the averaged RF spectra (Figure 3b) [35].

5. Non-Manakovian transmission

We continue to assume that the CNLSs provide a sufficiently accurate basis for signal modeling in a single birefringent waveplate and will compare in the following their nonlinearity with corresponding terms of often-used simplifications such as the Manakov equation [13], the Manakov-PMD equation [8], and GN theory [23]. It turns out that the nonlinearity used in our NLDP calculus differs significantly from simplified forms often exploited in today’s system simulations. For this purpose, we examine a known transformation of the CNLS into the ME for regular birefringent fiber, which leads to an inconsistency in the foundation of nonlinear propagation theory. This commonly applied assumption biases microscopic NL polarization effects. As for non-Manakovian effects such as the NLDP, we like to refer to when propagating signals cannot be sufficiently accurately described by applying equations of Manakovian type or the Manakov-PMD equation, which possess averaged and symmetric Kerr nonlinearities.

5.1 Shortfalls of an averaged Kerr nonlinearity

Here we review key steps for introducing an averaged Kerr nonlinearity into the CNLS and point out the associated shortfalls. The two main concepts [22, 41] to derive an averaged Kerr nonlinearity can be found in the literature. We follow a well-documented scheme in [41], which considers a piece of isotropic fiber (with zero birefringence) and neglect for now linear propagation effects such as attenuation and all forms of dispersion. Furthermore, like all published derivations we assume for $A_{x\left(y\right)}^{\mathrm{\Sigma }}$ monochromatic fields at $\omega =0$, which some papers paraphrase as fully polarized light³. In this case Eqs. (7) and (8) reduces to polarization-dependent differentials that equal small NL distortions

$$\frac{\partial A_x^{\mathrm{\Sigma }}}{\partial z}=j\mathrm{\gamma }\left({\left|A_x^{\mathrm{\Sigma }}\right|}^2+\frac{2}{3}{\left|A_y^{\mathrm{\Sigma }}\right|}^2\right)A_x^{\mathrm{\Sigma }}+j\mathrm{\gamma }\frac{1}{3}{A}_x^{{\mathrm{\Sigma }}^{\ast }}{A}_y^{{\mathrm{\Sigma }}^2}$$E56

$$\frac{\partial A_y^{\mathrm{\Sigma }}}{\partial z}=j\mathrm{\gamma }\left({\left|A_y^{\mathrm{\Sigma }}\right|}^2+\frac{2}{3}{\left|A_x^{\mathrm{\Sigma }}\right|}^2\right)A_y^{\mathrm{\Sigma }}+j\mathrm{\gamma }\frac{1}{3}{A}_y^{{\mathrm{\Sigma }}^{\ast }}{A}_x^{{\mathrm{\Sigma }}^2}.$$E57

For normalized variables

$${\left|A_x^{\mathrm{\Sigma }}\right|}^2+{\left|A_y^{\mathrm{\Sigma }}\right|}^2=1$$E58

Eqs. (56) and (57) can concisely be written as

$$\frac{\partial A_x^{\mathrm{\Sigma }}}{\partial z}=j\mathrm{\gamma }\left(2A_x^{\mathrm{\Sigma }}+A_x^{\mathrm{\Sigma }\ast }\left({A_x^{\mathrm{\Sigma }}}^2+{A_y^{\mathrm{\Sigma }}}^2\right)\right)$$E59

$$\frac{\partial A_y^{\mathrm{\Sigma }}}{\partial z}=j\mathrm{\gamma }\left(2A_y^{\mathrm{\Sigma }}+A_y^{\mathrm{\Sigma }\ast }\left({A_x^{\mathrm{\Sigma }}}^2+{A_y^{\mathrm{\Sigma }}}^2\right)\right).$$E60

Similar as in Section 4.6, we evaluate the magnitude of the polarization state-dependent NL distortion by expressing the field with a Jones vector

$$\left[\begin{array}{c}{A}_x^{\mathrm{\Sigma }}\\ A_y^{\mathrm{\Sigma }}\end{array}\right]=\left[\begin{array}{c}\text{cosϑ}cos\epsilon +\mathrm{j}\text{sinϑ}sin\epsilon \\ \text{sinϑ}cos\epsilon -\mathrm{j}\text{cosϑ}sin\epsilon \end{array}\right]{\mathrm{e}}^{\mathrm{i\psi }}.$$E61

But here, the time-independent variables $\vartheta ,ϵ,\psi$are functions of the distance $z$. As empirically observed, the SOP of a cw light propagating through randomly and rapidly changing fiber birefringence is uniformly distributed on the Poincare sphere with corresponding densities for the angles $\vartheta ,\epsilon$ as assumed for Eq. (46).

Inserting Eq. (61) into the left sides of Eqs. (59) and (60) leads to a differential term for the phase

$${\left[\begin{array}{cc}{A}_x^{\mathrm{\Sigma }}& A_y^{\mathrm{\Sigma }}\end{array}\right]}^{†}{\partial }_z\psi$$E62

that can be evaluated as scalar when multiplying both sides of the equation by ${\left[\begin{array}{cc}{A}_x^{\mathrm{\Sigma }}& A_y^{\mathrm{\Sigma }}\end{array}\right]}^{\ast }$. Averaging the differentials across the Poincare sphere by means of Eq. (46) must equate to zero for terms including ${\partial }_z\vartheta$and ${\partial }_z\epsilon$ due to symmetry considerations. Thus, on the left only a term in ${\partial }_z\psi$ remains. On the right, $\left({\left|A_x^{\mathrm{\Sigma }}\right|}^2+{\left|A_y^{\mathrm{\Sigma }}\right|}^2\right)$ is factored out, which is no longer assumed to be limited to unity, but represents time- and distance-dependent power levels of a signal. When adding the time-dependent differentials on the left sides that describe linear propagation effects and including the absorption term, then the ME as stated in form of Eqs. (5) and (6) with an averaged Kerr NL follows. As this derivation assumes an isotropic fiber, ${\mathrm{\beta }}_1,{\mathrm{\beta }}_2$ must be defined as x-, y-independent. Some closed-form solutions of Eqs. (5) and (6) are known (solitons). Its nonlinearity is symmetric in the field components. But strictly speaking, the above derivation holds only for cw light and so Eqs. (5) and (6). The reintroduction of Eq. (58) as time- and distance-dependent quantity without expanding term (62) by additional corresponding differentials violates the derivation’s assumptions and complicates it. Literature justifies the above derivation with lacking phase matching of mixing terms that can be neglected. However, if the field contains components that are closely separated, a low-frequency beat tone can occur. When sufficient phase matching is present (as in the case of NLDP), measurable modulation of the probe can appear. Additionally, the derivation does not mirror correlations among differentials in ${\partial }_z\vartheta$ or ${\partial }_z\epsilon$ taken at two different positions separated by PMD. But as such, they are essential for understanding the spectral properties of NLDP.

When choosing Eqs. (5) and (6) as starting point for the splitting of the overall NL operator as explicated by Eqs. (12) and (13), the symmetric and asymmetric NL perturbation terms change their relative weights. Hence, both phenomena NLDP and NL phase noise cannot be represented at the same time using the same effective average Kerr nonlinearity.

5.2 Manakovian simulators in telecom

A common technique in modern optical communications is to polarization multiplex two orthogonal channels at same wavelength, which maximizes the spectral efficiency of a system. When representing in a simplified picture transmitted data symbols per channel by optical wave packets, their instantaneous common receive SOP is equally blurred by NLDP in both azimuth and polar angles on the Poincare sphere. Especially, fast SOP motions in azimuthal direction (assuming the individual channels possess polarizations aligned with the x-y coordinates) impair coherent cross talk at high receive OSNR for advanced modulation formats and reduce established limits for the channel capacity. Current research on capacity limits of fiber channels has not explicitly considered NLDP. These theories apply Shannon’s theorem while computing NL signal distortions by means of the Manakov equation [42, 43, 44, 45, 46, 47, 48, 49]. They need to be revisited when more accurate estimates are required.

To ease computations, most of the industrial link simulations resort to a type of ME when estimating NL transmission penalties. For signals with low-density constellation (e.g., PM-QPSK), the small NLDP-induced SOP scattering has little performance impact and a Manakovian simulation can be a good approximation. As reference, for a link with moderate Kerr nonlinearity and transpacific transmission distance, we experimentally found an apex angle ∼11° for NLDP-induced scattering on the Poincare sphere [40].

The Manakov-PMD equation as defined in [8] includes the 8/9 factor for an averaged Kerr nonlinearity. For regular birefringent fiber, its NL PMD term becomes irrelevant and only its linear PMD term distinguishes it from an ordinary ME. However, the averaging of its Kerr nonlinearity that yields its 8/9 weighting questions its claim of universal acceptance as a governing model.

The GN model, another widely discussed approach for simulating optical communication channels, is a technique to analytically solve the ME for a weakly NL WDM system with D+ propagation [23, 50]. Its main advantage resides in the derivation of variances for NL field distortions that can be linearly added to noise powers (repeater ASE) in the SNR formula to determine a BER. Note, these NL field distortions should not be confused with polarization state noise generated by NLDP. The GN model does not include PMD effects, which lead to a linear SOP decorrelation between the probe and the beating ASE components [39]; but this decorrelation is essential for the foundation of the probe’s NLDP-induced Stokes vector spectrum, which is interferometrically formed by long-range NL interactions (Section 3 [35]).

6. Conclusions

Precise modeling of NL signal propagation in optical fibers is critical for maximizing the data capacity of long-haul communication systems. It balances signal powers and received OSNR to mitigate nonlinearities. Over the past five decades, simplified techniques have been developed to efficiently compute NL propagation in fiber. They adapt models for ideal or piece-wise linear birefringent fiber to simulate propagation paths with randomly varying birefringence.

Recently, a novel transmission phenomenon in fiber to which we refer to as NL DePolarization (NLDP) has been introduced. Unpolarized ASE depolarizes a co-propagating probe in long-haul communication systems and lab test beds due to the fiber Kerr nonlinearity. This phenomenon has proven elusive to simpler propagation modeling.

We have described NLDP by means of propagation-dependent SOP speed histograms. And under some simplifying assumptions, our outlined analytical model yields a closed-form solution for NLDP-induced SOP speed in single and multiple span systems. Although small compared with other polarization effects, this phenomenon leads to a qualitatively different microscopic understanding of nonlinear light propagation in fiber. An antisymmetric perturbation operator in the CNLS generates phase noises that produce the SOP fluctuations. A major aspect of our model forms the PMD dependence of NLDP, which fundamentally differentiates it from other NL polarization phenomena such as NL polarization rotation. NLDP is based on long-range NL interactions where contributions from Kerr nonlinearities interfere over long transmission distances. Counterintuitively, the NL-generated Stokes vector spectrum of a signal’s polarization narrows with increasing propagation length.

Our derivations show that in the case of NLDP (non-Manakovian propagation), the solution spaces of the CNLS and the Manakov equation do not converge as suggested by earlier work. Under consideration of NLDP, reassessing fiber channel capacity simulations that are utilizing Manakovian-type equations can be beneficial for scientific purposes and could show small performance offsets.