## Abstract

We study theoretically, numerically and experimentally the nonlinear propagation of partially incoherent optical waves in single mode optical fibers. We revisit the traditional treatment of the wave turbulence theory to provide a statistical kinetic description of the integrable scalar NLS equation. In spite of the formal reversibility and of the integrability of the NLS equation, the weakly nonlinear dynamics reveals the existence of an irreversible evolution toward a statistically stationary state. The evolution of the power spectrum of the field is characterized by the rapid growth of spectral tails that exhibit damped oscillations, until the whole spectrum ultimately reaches a steady state. The kinetic approach allows us to derive an analytical expression of the damped oscillations, which is found in agreement with the numerical simulations of both the NLS and kinetic equations. We report the experimental observation of this peculiar relaxation process of the integrable NLS equation.

©2011 Optical Society of America

## 1. Introduction

As far as the influence of dissipative effects such as stimulated Raman scattering can be ignored, light propagation in single-mode optical fibers is well described by the one-dimensional (1D) scalar nonlinear Schrodinger (NLS) equation (see Eq. (1) below) [1]. This equation is integrable and has a class of special solutions called bright and dark solitons, which are sustained in the anomalous (focusing) and normal (defocusing) dispersion regimes respectively. During the past fifty years, the question of the interaction among solitons has been extensively studied by using the method of the inverse scattering transform [2, 3]. Recently the formation of schock-waves in the defocusing and nonlinear regime of the soliton dynamics has been studied [4, 5]. The evolution of a dense gas of uncorrelated NLS solitons has been also examined in Ref. [6]. In particular a general method to derive kinetic equations describing the evolution of the spectral distribution function of solitons has been proposed [6]. This method seems promising for future investigations of nonlinear propagation of incoherent optical waves in single mode optical fibers.

These questions are of practical importance for the fiber propagation of light beams delivered by highly multimode lasers or intense ASE light sources. The evolution of the optical power spectrum of incoherent waves has been analyzed in the framework of coherence theory for zero group velocity dispersion [7]. Taking into account second-order dispersion of the fiber, the statistical properties of optical wave systems ruled by NLS-like equations have been mainly examined by making use of the wave turbulence (WT) theory [8–15].

WT theory is based on a perturbation expansion procedure in which linear dispersive effects are supposed to dominate nonlinear effects. In this way a large separation of the linear and nonlinear length scales takes place, so that the statistics of the field can be assumed to be approximately Gaussian. This allows one to achieve the closure of the infinite hierarchy of moment’s equations, thus leading to a WT kinetic description of the evolution of the power spectrum of the field (i.e. the second-order moment of the gaussian field) [9, 10]. It is now well established that WT theory provides a detailed description of the process of thermalization and condensation of optical waves that occur in nonintegrable wave systems ruled by 2D or 3D NLS equations [11, 13]. WT theory has been also considered in purely 1D wave systems in various circumstances [16–19]. In particular, it has been shown to properly describe the asymmetric evolution of the spectrum of an incoherent light wave propagating near by the zero dispersion wavelength of an optical fiber [14]. In this specific case, the integrability of the 1D NLS equation is broken by the third-order dispersion term and the spectral asymmetry arises from a phenomenon of anomalous thermalization originating in degenerate resonances of the wave system [15].

Non-integrability of the model equation is usually considered as a prerequisite for the applicability of WT theory, because it implies a process of irreversible diffusion in phase space that is consistent with the formal irreversibility of the kinetic equation. On the other hand, the dynamics of integrable systems is essentially periodic in time, reflecting the underlying regular phase-space of nested-tori [11, 20]. In this respect, applying the WT theory to the integrable 1D scalar NLS equation, all collision terms in the kinetic equation vanish identically at any order [8]. Accordingly, the WT theory predicts that the spectrum of an incoherent light wave propagating in a single-mode optical fiber does not evolve during the propagation.

This conclusion is in contradiction with our experiments and numerical simulations reported in 2006 [21], in which evidence is provided that a significant evolution of the power spectrum of the optical field does occur during its propagation in single-mode optical fibers [21]. In particular, the width of the power spectrum has been shown to change in a non-monotonic way with the propagation distance [21]. It is important to underline that such spectral evolution cannot be simply ascribed to the fact that a real-world optical fiber experiment necessarily includes some small effects that break the integrability of the NLS equation, or that the split-step numerical scheme does not preserve the integrability of the equation [22]. The experimental and numerical results reported in [21] reveal that there exists an inadequate understanding of incoherent light propagation in single mode fibers. This issue was recently addressed by Soh *et al.*, who proposed to modify the traditional approach of the WT theory [23]. Using numerical simulations, they analyzed the behavior of the kinetic equation in the early stage of propagation of the optical field in the single mode fiber. An interesting agreement has been obtained between the numerical integration of their kinetic equation and the numerical simulations of the 1D NLS equation.

In this paper, we study theoretically, numerically and experimentally the nonlinear evolution of incoherent light waves whose propagation is governed by the integrable 1D scalar NLS equation. Starting from the treatment proposed by Soh *et al.* [23], we revisit the traditional WT approach of the 1D scalar NLS equation by considering that the fourth-order moment of the field is not necessarily a stationary quantity. A kinetic equation governing the evolution of the second-order moment of the field is obtained. Contrary to the conventional WT kinetic equation, the collision term does not vanish, but relaxes rapidly to zero. In spite of such a fast relaxation, the spectrum of the field may exhibit significant changes depending on the initial conditions. Considering Gaussian-shaped initial conditions, the evolution of the spectrum is characterized by a rapid growth of the spectral tails, their subsequent damped oscillations, and their ultimate relaxation toward the stationary state. The kinetic approach allows us to derive an analytical expression of evolution of the spectrum, whose damped oscillations have been found in agreement with the numerical simulations of both the NLS and kinetic equations. We report for the first time the experimental observation of the rapid growth of frequency components located in the tails of the spectrum.

## 2. Theoretical Treatment

The starting point of our analysis is the dimensionless integrable 1D scalar NLS equation:

*z*represents the distance of propagation in the optical fiber while

*t*measures the time in a reference frame moving with the incoherent wave

*ψ*(

*z,t*). The time variable

*t*has been normalized to

*τ*

_{0}= 1/Δ

*ω*, where Δ

*ω*is related to the width of the power spectrum of the incoherent light field launched inside the optical fiber. In this paper, we will particularly consider initial conditions corresponding to gaussian power spectra:

*n*(

_{ω}*z*= 0) =

*n*

_{0}

*exp*(−

*ω*

^{2}/Δ

*ω*

^{2}). For this particular shape, Δ

*ω*represents the half-width at 1/

*e*of the power spectrum. The space variable

*z*has been normalized with respect to the linear dispersion length

*L*= 2/(

_{D}*β*

_{2}Δ

*ω*

^{2}), where

*β*

_{2}is the second-order dispersion coefficient of the optical fiber. With these conventions, the amplitude of the field

*ψ*(

*z,t*) is expressed in units of $1/\sqrt{{L}_{D}\cdot \gamma}$ where

*γ*is the nonlinear Kerr coefficient of the fiber. The function

*ψ*and the variables

*z*and

*t*can be rephrased in physical units through the transformations

*z*→

*zL*,

_{D}*t*→

*t*

*τ*

_{0}=

*t*/Δ

*ω*and $\psi \hspace{0.17em}\to \hspace{0.17em}\psi /\sqrt{\gamma {L}_{D}}$.

*σ*=

*sign*(

*β*

_{2}) represents the sign of the second-order dispersion coefficient (

*σ*= +1 for the defocusing case (normal dispersion) and

*σ*= −1 for the focusing case (anomalous dispersion)). The linear dispersion relation simply reads

*k*(

*ω*) =

*σω*

^{2}. For incoherent waves considered in this paper, results that will be determined from WT theory are not sensitive to the sign of

*σ*. For the sake of clarity, let us mention that numerical simulations of Eq. (1) presented in this paper have been made only in the normal dispersion regime (

*σ*= +1) corresponding to the experiments presented in Sec. 3. The influence of the sign of

*σ*on mean spectra numerically calculated from integration of Eq. (1) is briefly discussed at the end of Sec. 2.

The nonlinear length *L _{NL}* is defined in the usual way as

*L*= 1/(

_{NL}*γP*

_{0}), where

*P*

_{0}represents the power carried by the incoherent wave. The NLS equation [Eq. (1)] conserves the power

*N*= ∫|

*ψ*(

*z,t*)|

^{2}

*dt*and the total energy

*H*=

*H*+

_{L}*H*that has a linear (kinetic) contribution

_{NL}*H*= ∫

_{L}*k*(

*ω*)

*$\tilde{\psi}$*(

*z*,

*ω*)

*dω*and a nonlinear contribution ${H}_{NL}\hspace{0.17em}=\hspace{0.17em}\frac{1}{2}\hspace{0.17em}\int {|\psi (z,\hspace{0.17em}t)|}^{4}dt$.

Considering a random wave characterized by a stationary statistics, the second-order moment of the field is

*$\tilde{\psi}$*(

*z*,

*ω*) is the Fourier transform of

*ψ*(

*z,t*) defined as $\tilde{\psi}(z,\hspace{0.17em}\omega )\hspace{0.17em}=\hspace{0.17em}\frac{1}{\sqrt{2\pi}}\hspace{0.17em}{\int}_{-\infty}^{+\infty}\psi (z,t){e}^{-i\omega t}dt$. Following the standard Hasselman derivation of kinetic equations for a random wave with stationary statistics, the sixth-order moment can be factored in products of second-order moments and we obtain one equation for the second-order moment and one equation for the fourth-order moment that read respectively [24]:

*𝒩*(

*z*) =

*n*

_{ω1}(

*z*)

*n*

_{ω3}(

*z*)

*n*

_{ω4}(

*z*) +

*n*

_{ω2}(

*z*)

*n*

_{ω3}(

*z*)

*n*

_{ω4}(

*z*) –

*n*

_{ω1}(

*z*)

*n*

_{ω2}(

*z*)

*n*

_{ω3}(

*z*) –

*n*

_{ω1}(

*z*)

*n*

_{ω2}(

*z*)

*n*

_{ω4}(

*z*) and Δ

*k*=

*k*(

*ω*

_{1}) +

*k*(

*ω*

_{2}) –

*k*(

*ω*

_{3}) –

*k*(

*ω*

_{4}).

The usual way to proceed from this point in WT theory consists in integrating Eq. (5). The exact solution of Eq. (5) reads:

*k*= 0). For Δ

*k*≠ 0, the contribution of the fast oscillating function

*e*

^{iΔkz}is considered as being non essential for propagation distances larger than 1/Δ

*k*. With these assumptions, one obtains the following standard WT kinetic equation [10, 12, 24]

*z*= 0 has been supposed gaussian so that the fourth-order moment ${J}_{1,2}^{3,4}\hspace{0.17em}(z=0)$ can be factorized into a product of second-order moments (〈

*$\tilde{\psi}$*(0,

*ω*

_{1})

*$\tilde{\psi}$*(0,

*ω*

_{2})

*$\tilde{\psi}$*

^{*}(0,

*ω*

_{3})

*$\tilde{\psi}$*

^{*}(0,

*ω*

_{4})〉 =

*n*

_{ω1}(0)

*n*

_{ω2}(0)[

*δ*(

*ω*

_{1}–

*ω*

_{3}).

*δ*(

*ω*

_{2}–

*ω*

_{4}) +

*δ*(

*ω*

_{1}–

*ω*

_{4}).

*δ*(

*ω*

_{2}–

*ω*

_{3})]). The fourth-order moment is therefore a real quantity having a vanishing contribution when Eq. (6) is substituted into Eq. (4).

As underlined by Zakharov in Ref. [8], the collision term found in the right-hand side of Eq. (7) is identically equal to zero because the term *𝒩* (*z*) vanishes as a result of the integration over the Dirac delta functions. In opposition with two- or three-dimensional geometry, the phase-matching conditions found in 1D scalar NLS equation only permit trivial interactions among frequency components (i.e. *ω*
_{3} = *ω*
_{1,2} and *ω*
_{4} = *ω*
_{2,1}). This results in a collision term that is identically equal to zero. Therefore the only possible way to describe an evolution of *n*
_{ω1} in a one-dimensional geometry from WT theory consists in taking into account non-phase matched interactions among frequency components (i.e. Δ*k* ≠ 0).

In Ref. [23], Soh *et al.* have proposed a theoretical treatment of nonlinear propagation of an incoherent light wave inside an optical fiber. Although full details of their calculation are not explicitly given in Ref. [23], the way through which Soh *et al.* establish their kinetic equation (Eq. (5) of Ref. [23]) can be easily understood from our previous considerations. Equation (6) is simply substituted into Eq. (4) yielding the following kinetic equation with a collision term including a cosine oscillating function (see also Eq. (5) of Ref. [23]):

*et al.*in Ref. [23] therefore gives a kinetic equation including a collision term that does not vanish in a straightforward way, contrary to the collision term usually obtained from a standard treatment of WT theory [10, 12, 24].

Before undertaking our own analytical treatment of WT equations, let us examine some preliminary results from numerical simulations of Eq. (1). Our numerical simulations of Eq. (1) are made by using periodic boundary conditions and a standard pseudo-spectral step-adaptative numerical scheme. In all the simulations of the NLS equation [Eq. (1)] presented in this paper, the random phases *ϕ*(*ω*) of the Fourier components of the field taken as initial condition (*ψ*(*z* = 0, *ω*) = |*ψ*(*z* = 0, *ω*)|*e*
^{iϕ(ω)}) are uniformly distributed between 0 and 2*π* [21]. Moreover the power spectra are computed from an average procedure made over an ensemble of 100 numerical integrations of Eq. (1).

Generally speaking numerical simulations of Eq. (1) show that the power spectrum of the incoherent wave launched at *z* = 0 transiently changes from nonlinear propagation. It reaches a statistical steady state after a transient stage whose duration and amplitude critically depend on the shape of the power spectrum taken as initial condition. The most significant evolutions are observed when the tails of the power spectrum decay more sharply than an exponential function (i.e. for power spectra such as *n _{ω}* (

*z*= 0) =

*n*

_{0}

*exp*(−|

*ω*/Δ

*ω*|

*) with*

^{α}*α*> 1). This is illustrated in Fig. 1 that shows power spectra (black lines) computed from numerical simulations of Eq. (1) with the following gaussian power spectrum (dashed black lines) taken as initial condition:

*ω*at 1/

*e*of the power spectrum is unity and with our normalization, the parameter

*n*

_{0}becomes a measure of the ratio between the linear dispersive length

*L*and the nonlinear length

_{D}*L*: ${n}_{0}\hspace{0.17em}=\hspace{0.17em}2\sqrt{\pi}{L}_{D}/({L}_{NL})$. As shown in Fig. 1, numerical simulations of Eq. (1) evidence an irreversible evolution of the power spectrum of the incoherent light field. It reaches a statistical steady state after a short propagation distance that is typically of the order of 0.1 in our dimensionless variables (i. e.

_{NL}*z*∼ 0.1

*L*∼ 0.2/(

_{D}*β*

_{2}Δ

*ω*

^{2}) in physical units). The central part of the spectrum does not change much from nonlinear propagation whereas the the wings of the spectrum are significantly modified. Their shape does not remain gaussian but it becomes exponential. Let us emphasize that the power carried by frequency components lying in the wings of the spectrum increases by several orders of magnitude. Note that this phenomenon has already been illustrated in Ref. [23].

In addition to the simulations of NLS equation [Eq. (1)], we have made numerical simulations of Eqs. (4) and (5). Figure 1 shows that Eqs. (4) and (5) describe the changes occurring in the wings of the power spectrum in a quantitative way over a great number of decades whatever the value of *n*
_{0} (see spectra plotted in red lines in Fig. 1). As previously mentioned, *n*
_{0} is proportional to the ratio between linear and nonlinear lengths. Changes in *n*
_{0} are therefore associated with changes in the ratio *ɛ* between the *H _{NL}* and

*H*(

_{L}*ɛ*=

*H*/

_{NL}*H*). Figure 1(a) shows that numerical simulations of Eqs. (4) and (5) are in quantitative agreement with simulations of Eq. (1) over more than 20 decades in the linear (kinetic) regime in which

_{L}*ɛ*=

*H*/

_{NL}*H*≃ 0.05. As shown in Fig. 1(b), this quantitative agreement is preserved over ∼ 15 decades even in the nonlinear regime in which

_{L}*ɛ*=

*H*/

_{NL}*H*≃ 0.5. Whatever the interaction regime explored in our numerical simulations (i.e. from

_{L}*ɛ*∼ 0.01 to

*ɛ*∼ 1), the central part of the spectrum that carries the essential of the power of the incoherent wave is not significantly modified from nonlinear propagation.

In the theoretical analysis made in the present paper, we use the previous result and we suppose that the power spectrum *n _{ω}* (

*z*= 0) of the incoherent field only slightly changes with propagation distance

*z*. With this assumption, the term

*𝒩*(

*z*′) in the integral of Eq. (6) is approximated by its value

*𝒩*(

*z*= 0) at

*z*= 0 and it is extracted from the integral. Substituting Eq. (6) into Eq. (4), we obtain the following kinetic equation

*z*≫ 1/Δ

*k*, the function $\frac{\text{sin}(\mathrm{\Delta}kz)}{\mathrm{\Delta}k}$ in Eq. (10) tends to the Dirac delta function

*δ*(Δ

*k*) =

*δ*(

*k*(

*ω*

_{1}) +

*k*(

*ω*

_{2}) –

*k*(

*ω*

_{3}) –

*k*(

*ω*

_{4})) and the collision term found in the kinetic equation [Eq. (7)] vanishes. However as long as

*z*∼ 1/Δ

*k*, non-phase-matched interactions among spectral components cannot be neglected. The collision integral found in Eq. (10) is not identically equal to zero and changes in the power spectrum

*n*

_{ω1}can be observed despite the integrability of the wave system.

The theoretical analysis can be pursued by simplifying the kinetic equation [Eq. (10)] from the integration over *ω*
_{2}. This gives:

*k*=

*k*(

*ω*

_{1}) +

*k*(

*ω*

_{3}+

*ω*

_{4}–

*ω*

_{1}) –

*k*(

*ω*

_{3}) –

*k*(

*ω*

_{4}) = 2

*σ*(

*ω*

_{1}–

*ω*

_{3})(

*ω*

_{1}–

*ω*

_{4}). The term

*ℳ*(

*z*= 0) is given by $\mathcal{M}(z=0)\hspace{0.17em}=\hspace{0.17em}{n}_{{\omega}_{1}}^{0}{n}_{{\omega}_{3}}^{0}{n}_{{\omega}_{4}}^{0}\hspace{0.17em}+\hspace{0.17em}{n}_{{\omega}_{3}}^{0}{n}_{{\omega}_{4}}^{0}{n}_{{\omega}_{3}+{\omega}_{4}-{\omega}_{1}}^{0}\hspace{0.17em}-\hspace{0.17em}{n}_{{\omega}_{1}}^{0}{n}_{{\omega}_{3}}^{0}{n}_{{\omega}_{3}+{\omega}_{4}-{\omega}_{1}}^{0}\hspace{0.17em}-\hspace{0.17em}{n}_{{\omega}_{1}}^{0}{n}_{{\omega}_{4}}^{0}{n}_{{\omega}_{3}+{\omega}_{4}-{\omega}_{1}}^{0}$ where ${n}_{{\omega}_{i}}^{0}\hspace{0.17em}=\hspace{0.17em}{n}_{{\omega}_{i}}(z=0)$ (

*i*= 1, 3, 4). Further simplification of the kinetic equation can be achieved by considering that two parts of the collision term found in Eq. (11) balance each other. We thus obtain:

Power spectra plotted in blue lines in Fig. 1 are computed from numerical simulations of Eq. (12). Figure 1(a) shows that our approximation is very good in the linear regime (*n*
_{0} = 0.1, *ɛ* = *H _{NL}*/

*H*= 0.05). In the slightly nonlinear regime (Fig. 1(b) computed for

_{L}*n*

_{0}= 1,

*ɛ*=

*H*/

_{NL}*H*= 0.5, the approximation is less effective but a quantitative agreement with the simulation of Eq. (1) is nevertheless preserved over ∼ 8 decades. The first integral found in the right-hand side of Eq. (12) mainly contributes to change occurring in the wings of the spectrum whereas the second integral essentially contributes to changes occurring around the center of the spectrum (i. e. around

_{L}*ω*

_{1}∼ 0). The aim of this paper is to study the deep quantitative changes found in the tails of the spectrum (see Fig. 1). We now restrict our theoretical analysis to the only study of this specific point (

*ω*

_{1}>> Δ

*ω*). The problem thus reduces to

*x*=

*ω*

_{3}–

*ω*

_{1}/3 and

*y*=

*ω*

_{4}–

*ω*

_{1}/3. Using the variables

*x*and

*y*, the function ${n}_{{\omega}_{3}}^{0}\hspace{0.17em}{n}_{{\omega}_{4}}^{0}\hspace{0.17em}{n}_{{\omega}_{3}+{\omega}_{4}-{\omega}_{1}}^{0}$ is simply equal to ${n}_{0}^{3}\hspace{0.17em}\text{exp}(2({x}^{2}\hspace{0.17em}+\hspace{0.17em}{y}^{2}\hspace{0.17em}+xy))\hspace{0.17em}\text{exp}(-{\omega}_{1}^{2}/3)$, which corresponds to a peaked function localized around

*x*=

*y*= 0. With a gaussian power spectrum taken as initial condition, the growth rate of the power carried by a spectral component at a frequency

*ω*

_{1}falling in the tails of the spectrum is determined by the interaction among spectral components at frequencies ±

*ω*

_{1}/3 falling in the center of the gaussian spectrum. This result is not intuitive : the evolution of the component at a frequency

*ω*

_{1}is not driven by degenerate four-wave mixing among the pairs of frequencies (0, 0) and (−

*ω*

_{1}, +

*ω*

_{1}). Indeed, the dominant contribution corresponds to degenerate four-wave mixing among the pairs of frequencies (+

*ω*

_{1}/3, +

*ω*

_{1}/3) and (−

*ω*

_{1}/3, +

*ω*

_{1}). Note that this result is valid for any initial spectrum with an hypergaussian shape (i. e.

*n*(

_{ω}*z*= 0) =

*n*

_{0}

*exp*(−

*ω*

^{2}

*) where*

^{p}*p*> 0 is an integer)

Using an approximation commonly made with integral of oscillating functions, we consider that the term Δ*k* found in the denominator of Eq. (13) can be taken constant
$(\mathrm{\Delta}k\hspace{0.17em}=\hspace{0.17em}8\sigma {\omega}_{1}^{2}/9)$ and we extract it from the integral. We thus obtain

*k*in the sine function of Eq. (14) by $\mathrm{\Delta}k\hspace{0.17em}\simeq \hspace{0.17em}\sigma \frac{4{\omega}_{1}}{3}({\omega}_{4}\hspace{0.17em}-\hspace{0.17em}{\omega}_{3})$. We finally obtain

*ω*>> Δ

*ω*) with propagation distance

*z*. For the sake of clarity, we rephrase its expression in physical units:

*β*=

*β*

_{2}/2.

Equation (15) shows that the power *n*
_{ω1} of a spectral component taken in the wings of the gaussian spectrum grows with the propagation distance *z*. Moreover Eq. (15) shows that *n*
_{ω1} will reach a steady value after damped oscillations. This phenomenon is illustrated in Fig. 2 that shows the decaying oscillations of the power of two different spectral components. As illustrated in Fig. 2, there is a good quantitative agreement between the results obtained from the numerical integrations of Eq. (1) and of Eq. (15). Note that simulations of Eq. (1) have been realized over an ensemble of 100 realizations and that the curves plotted with black lines in Fig. 2 represent an averaged result. Note also that the decaying oscillations plotted in Fig. 2 are those leading to the steady-state spectrum plotted in Fig. 1(a).

As shown by Eq. (15), the spatial period
$\mathrm{\Lambda}\hspace{0.17em}\simeq \hspace{0.17em}\frac{9\pi}{4{\omega}_{1}^{2}}$ of the decaying oscillations is inversely proportional to
${\omega}_{1}^{2}$. Rephrasing Λ in physical units, the spatial period of the oscillations is:
${\mathrm{\Lambda}}^{phys}\hspace{0.17em}\simeq \hspace{0.17em}\frac{9\pi}{4}\hspace{0.17em}{\left(\frac{\mathrm{\Delta}\omega}{{\omega}_{1}}\right)}^{2}$
*L _{D}*. Let us recall that these oscillations have their physical origin in transient and non-phase-matched interactions among frequency components of the incoherent light field (see Eq. (10)). The period of the oscillations is determined by the fact that the spectral component at

*ω*

_{1}dominantly interact with the frequency components ${\omega}_{3}\hspace{0.17em}\simeq \hspace{0.17em}{\omega}_{4}\hspace{0.17em}\simeq \frac{{\omega}_{1}}{3}$ and ${\omega}_{2}\hspace{0.17em}\simeq \hspace{0.17em}-\frac{{\omega}_{1}}{3}$. In these conditions, the dominant spatial frequency Δ

*k*is around $\mathrm{\Delta}k\hspace{0.17em}=\hspace{0.17em}8\sigma {\omega}_{1}^{2}/9$. The damping of the oscillations is given by the gaussian function in Eqs. (15) and (16). It shows that the damping length measuring the distance needed for the wave system to reach its steady state scales as (

*βω*

_{1}Δ

*ω*)

^{−1}.

We now carefully compare our results with those established by Soh *et al* in Ref. [23]. First of all, Eqs. (8) and (9) of Ref. [23] only provide an upper bound to the change of the power of one frequency component. Rephrasing in a simplifying way the result obtained by Soh *et al.* in Ref. [23], we have indeed:

*B*in a parameter depending on

*γ*and

*P*

_{0}in a way that does not need to be explicited for our discussion. According to Soh

*et al.*, $\overline{\Omega}$

^{2}is a constant quantity defined as the “averaged root-mean squared full optical bandwidth”. With the averaging procedure proposed by Soh

*et al*, the upper-bound |

*Bsinc*($\overline{\Omega}$

^{2}

*z*)| is an oscillating and decaying function of

*z*that is independent of

*ω*

_{1}. In our detailed analysis, we have shown that spectral components at frequencies

*ω*

_{1}falling in the wings of the spectrum relax to the statistical steady state from damped oscillations whose properties (period and damping rate) fundamentally depend on the frequency

*ω*

_{1}.

Let us conclude Sec. 2 by a brief discussion about the influence of the sign of *σ* on the mean spectra found from numerical integration of Eq. (1) (see Fig. 1). The results determined above from WT theory are inherently insensitive to the sign of *σ* in linear (kinetic) regime (*ɛ* << 1). In the linear regime, coherent structures such as solitons do not emerge from the propagation of incoherent waves. Moreover the frequencies at which modulational instability can appear in anomalous dispersion regime are localized in the vinicity of the center of the spectrum. In other words, modulational instability does not influence the growth of the tails of the spectrum in the linear regime. We have checked this behavior from numerical integration of Eq. (1) and mean spectra of Fig. 1(a) do not depend of the sign of *σ* in this linear regime (*ɛ* = 0.05).

In the slightly nonlinear interaction regime in which *ɛ* = 0.5 (see Fig. 1(b)), numerical simulations of Eq. (1) contrarily show that the mean spectrum *n _{ω}* (

*z*= 1) depends on the sign of

*σ*. Numerical simulations made in the anomalous dispersion regime (

*σ*= −1) show that the mean spectrum

*n*(

_{ω}*z*= 1) reached at

*z*= 1 is broader than the one plotted in black line in Fig. 1(b) for the normal dispersion regime (

*σ*= +1). This feature can be interpreted from the fact that soliton-like structures and modulational instability begin to play a role in the anomalous dispersion regime at

*ɛ*= 0.5. We have shown that our analytical results obtained from kinetic theory are still in good agreement with numerical integration of NLS equation [Eq. (1)] at

*ɛ*= 0.5 (see Fig. 1(b)). In other words, the range of values of

*ɛ*=

*H*/

_{NL}*H*over which WT theory provides results that are quantitatively correct is wider in the normal dispersion regime than in the anomalous dispersion regime.

_{L}## 3. Experiments

In this section we present experimental results showing the transient changes in the spectrum that have been theoretically predicted and studied in Sec. 2. Our experimental setup is schematically shown in Fig. 3. The light source delivering a partially-coherent field is a homemade Nd:YVO4 continuous-wave laser. The laser is linearly-polarized and its central emission wavelength of ∼ 1064 nm is far from the zero-dispersion wavelength *λ*
_{0} of the optical fiber (*λ*
_{0} ∼ 1400 nm). The laser optical power spectrum has a full-width at half maximum (FWHM) of ∼ 0.15 nm (i.e. ∼ 36 GHz). Observing the time evolution of the laser power with a photodiode having a bandwidth much greater than the free spectral range Δ*ν _{FSR}* of the laser cavity (Δ

*ν*≃ 150 MHz), we did not notice any significant changes in the laser power on time scales of the order of 1/Δ

_{FSR}*ν*. We conclude that the laser output is made of approximately 200 longitudinal modes with phases that can be considered as being uncorrelated.

_{FSR}Let us explain why we have chosen to use a Nd:YVO4 laser rather than a more highly multimode laser like a Ytterbium-doped fiber laser [21]. As discussed in Sec. 2, the duration and the amplitude of the transitory evolution of the light power spectrum significantly change with the shape of the power spectrum of the incoherent light wave launched inside the optical fiber. Our Ytterbium-doped fiber lasers have power spectra that do not decay in a sharp way. On the other hand, our homemade Nd:YVO4 laser has a power spectrum that decay sharply enough to observe transient changes discussed in Sec. 2.

The power of the laser light is controlled by adjusting a half-wave plate (HWP1) placed between the Nd:YVO4 laser and a Faraday isolator. The laser light is launched inside a 1.5 km-long single-mode polarization-maintaining fiber (PMF). The polarization direction of the input light is adjusted along one of the birefringence axes of the PMF by using another half-wave plate (HWP2). In our experiments, the light wave remains linearly polarized all along the PMF and the extinction ratio between the two axes is greater than 20 dB. The nonlinear Kerr coefficient of the PMF is *γ* = 6 W^{−1}km^{−1} and its second-order dispersion coefficient is *β*
_{2} = 20 ps^{2} km^{−1} at 1064 nm.

In our experiment, we want to investigate the linear (kinetic) regime in which only the wings of the spectrum are modified from nonlinear propagation inside the PMF (*ɛ* = *H _{NL}*/

*H*<< 1). Up to now, only the nonlinear interaction regime has been explored in experiments (

_{L}*ɛ*=

*H*/

_{NL}*H*>> 1) [21]. In the nonlinear interaction regime, the power spectrum of the incoherent wave undergoes deep changes affecting its wings but also its central part. In other words, the FWHM of the power spectrum significantly changes in nonlinear interaction regime (see Fig. 2 of Ref. [21]). For incoherent light waves with wavelength around 1

_{L}*μ*m and with spectra spreading over ∼ 1 nm, the nonlinear interaction regime in a standard singlemode fiber typically refers to optical powers around ∼ 1 W [21]. The linear regime that we want to explore corresponds to optical powers typically lower than ∼ 1 mW. As illustrated in Fig. 1(a), the changes occurring in the wings of the power spectrum in the linear interaction regime involve spectral components carrying an optical power that is typically 10

^{10}times (equivalently 100 dB) lower than the power of spectral components lying in the center of the spectrum. The possibility to explore the linear interaction regime therefore critically depends on the sensitivity and the dynamic range of the optical spectrum analyzer (OSA) used in the measurement.

In our experiment, the output end of PMF is directly connected to the OSA (Ando AQ6317B). With the setup arranged in this way and using the highest degree of sensibility of the OSA, we have observed changes that can be undoubtedly associated with the effects discussed in Sec. 2 for an optical power of ∼ 12 mW at the output of the PMF. At this level of optical power, our experiment corresponds to a weakly nonlinear regime in which *ɛ* = *H _{NL}*/

*H*≃ 1. To explore the linear interaction regime, the optical power should be decreased by at least one order of magnitude but the detection of the slight changes occurring in the wings of the spectrum then becomes impossible with our OSA.

_{L}Figure 4 shows that the power of spectral components falling in the wings of the spectrum has increased with propagation distance whereas the central part of the power spectrum has not changed. Let us emphasize that the experiment only permits to explore the first stage of transient regime corresponding to a monotonic growth of the power carried by spectral components falling in the tails of the spectrum. In particular, we have not been able to observe damped oscillations of the power of these frequency components by using our setup (see Fig. 2).

Figure 4(a) shows that behaviors observed in the experiment are also found in numerical simulations of Eq. (1). In these simulations, the power spectrum of the field taken as initial condition is approximated by an hyper-gaussian profile
${n}_{\omega}^{0}\hspace{0.17em}=\hspace{0.17em}{n}_{\omega}(z=0)\hspace{0.17em}=\hspace{0.17em}{n}_{0}exp\left(-{\left(\frac{\omega}{A}\right)}^{4}\right)$. The value of the parameter *A* is adjusted to get the best fit between the hyper-gaussian profile and the experimental profile. As shown in Fig. 4(a), the power spectrum found in numerical simulations of Eq. (1) is in good quantitative agreement with the power spectrum recorded in experiments. Despite the fact that the experiment is made in a slightly nonlinear regime (*ɛ* = *H _{NL}*/

*H*≃ 1), a quantitative agreement over approximately four decades is obtained between numerical results computed from the integration of simplified kinetic Eq. (12) (blue line in Fig. 4(b)) and numerical results computed from integration of Eq. (1) (red lines in Fig. 4(b)). We have also made additional numerical simulations in order to study the influence of the exact shape of the laser power spectrum on the shape of the spectrum found at the output end of the PMF. As shown in Fig. 4(a), the power spectrum of the NdYVO4 laser is well fitted by an hyper-gaussian profile over approximately two decades but its tails deviates from this profile. Numerical simulations show that these deviations do not significantly influence the spectrum found at the output of the PMF. This can be understood by recalling that the evolution of the tails of the spectrum is only determined by the interaction among spectral components lying in the center of the spectrum (see Sec. 2)

_{L}## 4. Conclusion

In this paper, we have studied the nonlinear evolution of incoherent light waves whose propagation is governed by the integrable 1D scalar NLS Eq. (1). We have revisited an approximation usually made in WT theory for the study of problems such as optical wave thermalization or condensation. We have considered that the fourth-order moment of the field is not necessarily a stationary quantity. With the additional assumption that the power spectrum of the field does not change much with propagation distance, we have derived a kinetic equation governing the evolution of the second-order moment of the field. Contrary to the conventional WT kinetic equation, the collision term does not vanish, but relaxes rapidly to zero. In spite of such a fast relaxation, the spectrum of the field may exhibit significant changes depending on the initial conditions. Considering Gaussian-shaped initial conditions, we have shown that the evolution of the spectrum is characterized by a rapid growth of the spectral tails, their subsequent damped oscillations, and their ultimate relaxation toward the steady state. Using WT theory we have derived an analytical expression of evolution of the spectrum, whose damped oscillations have been found in agreement with the numerical simulations of both the NLS and kinetic equations. In particular, we have established a simple analytical expression for the spatial period of the damped oscillations. Experiments made with a cw laser emitting multiple longitudinal modes have evidenced the rapid growth of frequency components located in the tails of the spectrum.

Various theoretical and experimental developments can be made from the present work. the theoretical analysis presented here has been restricted to incoherent optical waves having power spectra that are initially gaussian. Further works are needed to understand the exact role taken by the shape of the power spectrum in the transient evolution leading to the steady state. From the experimental point of view, additional work is also needed to design and make an experiment that would permit to explore the linear kinetic regime. In this respect, the experimental observation of the damped oscillations found in this regime represents a striking challenge.

## Acknowledgment

This research was supported by the ANR-COSTUME 08-SYSC-004-03.

## References and links

**1. **G. P. Agrawal, *Nonlinear Fiber Optics* (Academic, 2001).

**2. **M. J. Ablowitz and H. Segur, *Solitons and the Inverse Scattering Transform* (SIAM, Society for Industrial and Applied Mathematics, 1981). [CrossRef]

**3. **W. Zhao and E. Bourkoff, “Interactions between dark solitons,” Opt. Lett. **14**, 1371–1373 (1989). [CrossRef] [PubMed]

**4. **C. Conti, A. Fratalocchi, M. Peccianti, G. Ruocco, and S. Trillo, “Observation of a gradient catastrophe generating solitons,” Phys. Rev. Lett. **102**, 083902 (2009). [CrossRef] [PubMed]

**5. **A. Fratalocchi, C. Conti, G. Ruocco, and S. Trillo, “Free-energy transition in a gas of noninteracting nonlinear wave particles,” Phys. Rev. Lett. **101**, 044101 (2008). [CrossRef] [PubMed]

**6. **G. A. El and A. M. Kamchatnov, “Kinetic equation for a dense soliton gas,” Phys. Rev. Lett. **95**, 204101 (2005). [CrossRef] [PubMed]

**7. **J. T. Manassah, “Self-phase modulation of incoherent light,” Opt. Lett. **15**, 329–331 (1990). [CrossRef] [PubMed]

**8. **V. E. Zakharov, “Turbulence in integrable systems,” Stud. Appl. Math. **122**, 219–234 (2009). [CrossRef]

**9. **A. C. Newell, S. Nazarenko, and L. Biven, “Wave turbulence and intermittency,” Physica D **152–153**, 520–550 (2001). [CrossRef]

**10. **S. Dyachenko, A. C. Newell, A. Pushkarev, and V. E. Zakharov “Optical turbulence: weak turbulence, condensates and collapsing filaments in the nonlinear Schrodinger equation,” Physica D **57**, 96–160 (1992). [CrossRef]

**11. **C. Connaughton, C. Josserand, A. Picozzi, Y. Pomeau, and S. Rica, “Condensation of classical nonlinear waves,” Phys. Rev. Lett. **95**, 263901 (2005). [CrossRef]

**12. **A. Picozzi, “Towards a nonequilibrium thermodynamic description of incoherent nonlinear optics,” Opt. Express **15**, 9063–9083 (2007). [CrossRef] [PubMed]

**13. **P. Aschieri, J. Garnier, C. Michel, V. Doya, and A. Picozzi, “Condensation and thermalization of classical optical waves in a waveguide,” Phys. Rev. A **83**, 033838 (2011). [CrossRef]

**14. **C. Michel, P. Suret, S. Randoux, H. R. Jauslin, and A. Picozzi, “Influence of third-order dispersion on the propagation of incoherent light in optical fibers,” Opt. Lett. **35**, 2367–2369 (2010). [CrossRef] [PubMed]

**15. **P. Suret, S. Randoux, H. R. Jauslin, and A. Picozzi, “Anomalous thermalization of nonlinear wave systems,” Phys. Rev. Lett. **104**, 054101 (2010). [CrossRef] [PubMed]

**16. **V. E. Zakharov, F. Dias, and A. Pushkarev, “One-dimensional wave turbulence,” Phys. Rep. **398**, 1 (2004). [CrossRef]

**17. **U. Bortolozzo, J. Laurie, S. Nazarenko, and S. Residori, “Optical wave turbulence and the condensation of light,” J. Opt. Soc. Am. B **26**, 2280 (2009). [CrossRef]

**18. **B. Barviau, B. Kibler, A. Kudlinski, A. Mussot, G. Millot, and A. Picozzi, “Experimental signature of optical wave thermalization through supercontinuum generation in photonic crystal fiber,” Opt. Express **17**, 7392 (2009). [CrossRef] [PubMed]

**19. **B. Barviau, B. Kibler, and A. Picozzi, “Wave turbulence description of supercontinuum generation: influence of self-steepening and higher-order dispersion,” Phys. Rev. A **79**, 063840 (2009). [CrossRef]

**20. **M. Tabor, *Chaos and Integrability in Nonlinear Dynamics* (Wiley, 1989).

**21. **B. Barviau, S. Randoux, and P. Suret, “Spectral broadening of a multimode continuous-wave optical field propagating in the normal dispersion regime of a fiber,” Opt. Lett. **31**, 1696–1698 (2006). [CrossRef] [PubMed]

**22. **B. M. Herbst and M. J. Ablowitz, “Numerically induced chaos in the nonlinear Schrodinger equation,” Phys. Rev. Lett. **62**, 2065–2068 (1989). [CrossRef] [PubMed]

**23. **D. B. S. Soh, J. P. Koplow, S. W. Moore, K. L. Schroder, and W. L. Hsu “The effect of dispersion on spectral broadening of incoherent continuous-wave light in optical fibers,” Opt. Express **18**, 22393–22405 (2010). [CrossRef] [PubMed]

**24. **V. E. Zakharov, V. S. L’vov, and G. Falkovich, *Kolmogorov Spectra of Turbulence I* (Springer, 1992).