Turbulence-induced rogue waves in Kerr resonators

Spontaneous emergence of self-organized patterns and their bifurcations towards a regime of complex dynamics in non-equilibrium dissipative systems is a paradigm of phase transition. Indeed, the behavior of these patterns in the highly nonlinear regime remains less explored, even in recent high-quality-factor resonators such as Kerr-nonlinear optical ones. Here, we investigate theoretically and experimentally the alteration of the resulting Kerr frequency combs from the weakly to the highly nonlinear regime, in the frameworks of spatiotemporal chaos, and dissipative phase transitions. We reveal the existence of a striking and easily accessible scenario of spatiotemporal chaos, free of cavity solitons, in a monostable operating regime, wherein a transition to amplitude turbulence via spatiotemporal intermittency is evidenced. Moreover, statistics of the light bursts in the resulting turbulent regime unveils the existence of rogue waves as extreme events characterized by long-tail statistics.


I. INTRODUCTION
The concept of order parameter description has played a key role in understanding dissipative structures and selforganized pattern in nonlinear systems.
Indeed, this description has been revealed to be a powerful theoretical tool for classifying and describing dissipative structures in weakly nonlinear regimes of nonequilibrium systems as the ones subject to a moderate external driving power. However, high external driving strengths lead dissipative systems to strongly nonlinear regimes where they exhibit extremely complicated dynamics such as spatiotemporal chaos and turbulence. In such a case, the degree of complexity of the highly nonlinear problem excludes any attempt to finding the appropriate order parameter that remains an elusive task. The development of a complete process for the study of the complex dynamics, including spatiotemporal chaos and turbulence, occurring in highly nonlinear regimes of dissipative systems is one of the most challenging open problems in nonlinear science [1][2][3]. Many recent advances in the understanding of complex dynamics have been driven by experimental and theoretical studies in modern optics, in fields as diverse as optical fiber cavities with quadratic or cubic nonlinear materials subject to externally injected radiations [4][5][6]. The field of optics, therefore, is ideally suited to the investigation of spatiotemporal chaos and turbulence in dissipative systems far from thermodynamic equilibrium (cavity and laser systems). Within a few decades, optical Kerr-nonlinear resonators have emerged as the paradigmatic setup for the study of externally driven nonlinear systems [7][8][9]. With length scales ranging from the meter to micrometer [10,11], their applications range from high-demand telecommunications [12] to precision spectroscopy and LIDAR (Light Detection And Ranging) systems [13] based on the coherent optical frequency * saliya.coulibaly@univ-lille.fr combs that can be delivered. Kerr resonators are also known for the property to continuously switch between a monostable (single-valued transmission curve) and bistable (S-shape transmission curve) regimes. Operating out of equilibrium, Kerr resonators can exhibit nontrivial outputs such as cavity solitons (localized coherent solution) [14][15][16] and the modulation instability [17,18] (process by which a homogeneous state breaks up into a periodic state). Unlike cavity solitons requiring the system to be bistable, the modulation instability (MI) was reported both in bistable and monostable regimes. There is currently a renewed interest for the MI in nearlyconservative physical systems after being tightly linked to the recent research activity and new developments in rogue wave phenomena [19]. Optical analogs of the hydrodynamic rogue waves are rare and short intense light pulses, characterized by long-tail statistics in the probability distribution of the intensity profile. They have been observed in numerous optical systems and may result from the emergence of coherent states or a collision of them. The Peregrine soliton, and Akhmediev breathers in conservative systems described by the nonlinear Schrödinger equation are the renowned examples [20]. Out of equilibrium systems, to which the ring resonators belong, have also been reported to exhibit rogue waves as extreme events. In this context, they result from deterministic processes leading to complex dynamical evolutions such as temporal chaos [21][22][23][24], spatiotemporal chaos and intermittency [25,26]. But from the general point of view, the identification of the necessary ingredients for the emergence of rogue waves and extreme events in dissipative systems remains a challenging open problem, including predictability [27], mechanisms of formation, and high sensitivity to noise sources [28] of these giant waves.
Returning to Kerr resonators, despite the tremendous interest on their dynamics, to the best of our knowledge, the only reference to rogue waves generation is a theoretical observation in the bistable regime with anomalous dispersion [29]. Moreover, more generally, it is MI of the monostable regime which has received much less attention. However, the MI has been demonstrated to seed an extremely interesting complex dynamics [30][31][32] when the pertinent control parameter is brought far from the instability onset. Examples include phase instability, spatiotemporal intermittency, and turbulence.
In this work, we investigate experimentally and theoretically the transition from periodic patterns with a triangular comb spectrum induced by MI in ring cavities towards more complex dynamics using a combination of three quantities: (i) The Lyapunov dimension density-The proof of the matters of the spatiotemporal chaos in an extended system is the existence of a continuous set of positive Lyapunov exponents. Besides, this set may increase with the number of degree of freedom [33]. Therefore, any quantity directly defined from a set of Lyapunov exponents yield to an extensive quantity. For example, all definition of the dimension of the attractor (the effective degree of freedom of the system) produces a quantity proportional to the volume of the system. The inverse of this proportionality coefficient is an intensive quantity with the same dimension of the volume of the system. Here, we compute the Lyapunov dimension density from the Kaplan-Yorke conjecture [33], which allows us to determine the characteristic size of the independent subsystems generated by the spatiotemporal chaos.
(ii) The two-point correlation length-A consequence of the MI is a diverging correlation length (i.e. a correlation that spawn overall a finite size system). Hence, any destabilization process of the MI may cause a meaningful reduction of the correlation lengths in the system [34][35][36][37][38]. In that case, short range correlation length given by the exponential decay rate of the correlation function is of particular interest to detect the phase transition-like processes when varying the control parameter of the system.
(iii) Intrinsic laminar length-Finite correlation range can be the result of a mixed state where spatiotemporal coherent (laminar) subsystems coexist with incoherent (or chaotic) ones. In this case, the probability distribution function of the laminar subsystem gaps can provide a significant signature about the nature of the complex dynamics. Namely, a power law distribution entails the spatiotemporal intermittency and an exponential distribution is a characteristic of the fully developed turbulence. However, it is more likely to find a range of the probability distribution of laminar regions that fits a decaying exponential law [2,30,31]. We will refer to this quantity as the intrinsic laminar length.
Note that, although the use of these quantities is widespread in hydrodynamics [2,30,31], only very few examples of their use are known in optics despite the many analogies between these areas. Besides, if the correlation length has been already associated with the Lyapunov dimension correlation or to the intrinsic exponential decay, the three quantities have not been used together yet. In the following, we numerically demonstrate that a monostable Kerr resonator can exhibit spatiotemporal chaotic regime. We have verify experimentally and numerically that the transition to this a chaotic behavior coincides with a diverging two-point correlation length in the system when varying the external driving strength. The correlation length allows us to identify a second transition in the dynamics. The numerical computation of the intrinsic laminar length helps to reveal this change as the transition from spatiotemporal intermittency to fully developed or amplitude turbulence. At this transition, we have identified the appearance of extreme events in the intensity profile. Finally, the temporal profile of the rogue waves is obtained through the analysis of the correlation function. Theoretical results and experimental observations are in a good agreement.
The paper is organized as follows. In section II the experimental setup, based on a coherently driven passive optical fibre ring cavity, is illustrated and the main ingredients for achieving the experiments are described. We also present both the governing equations of the experimental optical device and the reduced Lugiato-Lefever model. In section III spatiotemporal chaos resulting from modulation instability is analytically and numerically studied and the optimal range of parameters are determined. The results of our investigations on the complex dynamics appearing in the system, including spatiotemporal chaos, intermittency, and turbulence are contained in Sect. IV. Our concluding remarks are summarized in Section VI.

II. EXPERIMENTAL SETUP AND MODELIZATION
Our experimental setup based on a coherently driven passive optical fibre ring cavity is depicted in Fig. 1 . It is based on a resonant passive fibre ring cavity mainly made of 26.5 mlong segment of highly nonlinear optical fibre. The fibre combines a low group-velocity dispersion with a high nonlinear coefficient so as to enhance MI gain (β 2 = −0.89 ps 2 km −1 , γ = 10 W −1 km −1 at 1552.4 nm). It also exhibits a low thirdorder dispersion (β 3 0.01 ps 3 km −1 ) that can be neglected (the zero dispersion wavelength is located below 1500 nm far from the pump wavelength). A 90/10 input coupler is used to close the fibre loop cavity, whereas a 99/1 output coupler permits to extract and analyse the intracavity field. The two couplers are made of SMF28 fibre with a total length of 1.5 m that belongs to the cavity (i.e., the total cavity length is equal to 28 m). The SMF28 fibre exhibits the following usual parameters: β 2 = −21.7 ps 2 km −1 and γ = 1.2 W −1 km −1 . We split a continuous-wave (CW) laser at 1552.4 nm (linewidth <1 KHz) into two parts. The first part, the control wave, is injected into the cavity to stabilize and fix the linear detuning at the pump frequency ω 0 , with the help of a PID controller that finely tunes the laser wavelength. The second part, the pump wave, is intensity modulated to generate 2 ns square pulses at 7.36 MHz repetition rate (corresponding to the cavity's FSR). This stage simultaneously enables the increase of the pump peak powers and the circumvention of Brillouin backscattering within the cavity. The resulting quasi-CW pump is amplified by an erbium-doped fibre amplifier (EDFA), and launched into the cavity through the 90/10 coupler. To minimize their mutual interaction, the control and pump beams are counterpropagating. Input polarization states are controlled via polarization controllers to excite a neutral axis of the cavity fibre. Temporal and spectral characterizations of intracavity field are provided by an intensity autocorrelator (with a temporal resolution of 10 fs and a full time window limited to 80 ps) and a high-resolution (2.5 GHz) optical spectrum analyzer. The light field circulating in the effective fibre ring cavity (yellow box of Fig. 1) can be described by the following set of equations coupling the successive round-trip propagation described by the modified nonlinear Schrödinger equation (1a) to the synchronously coherent injection (1b) [39][40][41]: Here t R stands for the round-trip time which is the time taken by the pulse to propagate along the cavity with the group velocity, Φ 0 is the linear phase shift, θ (ρ) is the mirror transmission (reflection) coefficient, and L is the cavity length.
The complex envelope of the electric field inside the cavity at the m-th round trip is A m . Each of the coefficients β n is responsible for the n-th order dispersion, γ is the nonlinear coefficient, and α f is the attenuation along the fiber. The independent variable z refers to the longitudinal coordinate while T is the time in a reference frame moving with the group velocity of the light, and R(T ) is the nonlinear response including both instantaneous (Kerr effect) and delayed contributions (Raman effect). From transmission measurements with a continuous wave (CW) pump the cavity at 1552.4 nm (linewidth <1 kHz) the finesse F = 2π/θ eff is nearly 19. Hence, the total power losses of the cavity θ eff including fiber absorption and coupler losses are about 30% (θ eff = 0.3). From the fiber parameters, we obtain that our cavity is equivalent to a unique fiber ring cavity with β 2eff = −2 ps 2 km −1 and γ eff = 9.6 W −1 km −1 . Therefore, without loss of generality, the evolution of the electric field inside the cavity is well described by the Lugiato-Lefever equation (LL model) [41,42]: is the detuning with respect to the nearest cavity resonance k. The coefficient η = ±1 is the sign of the group velocity dispersion term and T n = |β 2 L| /(2α) and t ≡ mt R accounts for the continuous time introduced to account for the evolution of the intracavity field over the round trips.
In the bistable regime, the latter equation can exhibit soliton-like solutions (i.e., cavity solitons) [12,43,44]. Increasing the driving strength this soliton first undergoes an Andronov-Hopf instability yielding to self-pulsating localized state. For larger values of the pump intensity, this oscillatory localized state, in turn, becomes unstable. The evolution of the resulting complex state has been demonstrated to be of a spatiotemporal chaotic nature [45]. In this region of parameters, some of the emerging erratic spatiotemporal localized states have been statistically shown to satisfy the criteria of extreme events [29]. From a general point of view, the extensive literature in the highly nonlinear bistable regime contrasts with those of the monostable case that we are interested in. The reason is twofold, (i) bistability makes more accessible the experimental observation of solitons, (ii) the complex nonlinear dynamics is believed richer than in the monostable regime. Indeed, in the monostable regime, most studies are dominated by the weakly nonlinear analysis. In that case, only MI is present and the system remains stable and evolves in a regular way [42,[46][47][48][49][50]. The corresponding periodic pattern is characterized by a triangular comb spectrum [46][47][48][49][50] whose behavior in the highly nonlinear regime, however, remained unexplored.
To fill this lack, we have then performed a set of spectral and temporal measurements for a normalized linear detuning set to ∆ = 2(2πk − φ 0 )/θ eff = 0.55. Here k refers to the nearest resonance and θ eff /2 = α = π/F represents the power lost per roundtrip, with F the effective finesse of the cavity. For this value of the detuning, the setup is said to be in the monostable regime where no cavity soliton or complex dynamics related to bistability can take place, which prevents the formation of dissipative Kerr solitons. Only the Turing pattern or modulation instability (MI) can be observed at the onset of the cavity emission.  Figure 2 depicts the evolution of the intracavity spectrum while tuning the input peak power. The first sub-figure Fig. 2(a) represents the spectrum recorded for an input peak power P 0 = |A in | 2 = 0.16 W (slightly above the MI threshold, P th = 0.15 W). Two weak MI sidebands detuned by 407 GHz appears in the spectrum. By increasing the pump power, we observe their amplification and frequency detuning ( Fig. 2(b)). Also, the dynamics of cascade of MI sidebands appears progressively forming a triangular Kerr frequency comb (in a log-scale) composed by four harmonics of the MI frequency, equally spaced and detuned from the pump by 530 GHz (Fig. 2(c)). For higher pumping power ( Fig. 2(d)), intermediary spectral peaks centered between the MI bands start to growth.
While increasing further the pump power, a more complex and broader spectrum develops (Fig. 2(e)). For P 0 = 20 W, a broad continuous triangular spectrum is formed (Fig. 2(f)). We also observe the emergence of the broad Raman Stokes component detuned by −13 THz from the pump. We cannot further increase the power because of the formation of a Brillouin wave inside the cavity which destabilizes the locking loop.
Moreover, we performed temporal measurements of the intracavity field by means of a background-free secondharmonic autocorrelator. Figure 3 displays the distinct recorded autocorrelation traces corresponding to the spectra from Fig. 2. In the beginning, for a power just upon the MI threshold, we recorded a nearly flat trace corresponding to a quasi-continuous intensity field ( Fig. 3(a)). With the appearance and growth of MI sidebands, the trace becomes regularly modulated ( Fig. 3(b)). The contrast of this modulation is maximum when the triangular comb is reached (see Fig. 3(c)) indicating the formation of a short-pulse pattern on a finite background inside the cavity. The period of the pulse train is about 1.85 ps in agreement with the expected 530 GHz repetition rate driven by the MI frequency. The central peak corresponds to the autocorrelation of a single pulse, whereas the adjacent peaks are the cross-correlations between neighboring pulses. With the growth of the additional spectral peaks between the MI harmonics, the contrast seen in the trace starts to be degraded (Fig. 3(d)). Since the cross-correlations between neighboring pulses is highly sensitive to neighboring pulse differences, we can infer the degradation of the formed pulses train because of these new sidebands. For high pump powers, the modulation completely disappears and only a central peak of coherence remains (Fig. 3(e-f)). This reveals the strong intensity incoherence of the intracavity field (over a single round-trip).
We have also performed numerical simulations based on Eqs. (1). Figure 4 shows the different simulated intracavity spectra obtained after 10,000 cavity round trips (i.e., a stationary state circulating in the cavity) as a function of the input peak power. The numerical results are found to be in a good qualitative agreement with experiments. At first, MI sidebands and their harmonics appear and get amplified while increasing the pump power ( Fig. 4(a-c)). Then, inter- mediary sidebands start to grow and the spectrum broadens progressively until a continuous triangular shape is reached (Fig. 5(f)). Temporal analysis of the simulated intracavity field was also carried out (see Fig. 5). In particular, the autocorrelation traces (red lines, right axis) are found to be similar to the experimental ones. Fig. 5 also shows the corresponding field evolution (black lines, left axis). Details from the field evolution allow confirming the experimental results. At first, the appearance of the MI sidebands induces the modulation of the intracavity field until the formation of a train of short pulses (of about 360 fs) on a finite background (Fig. 5(a-c)).
Next, the growth of intermediary bands progressively degrades the regular pulse train, first, by doubling the period (Fig. 5(d)) and then by strongly decreasing the intensity coherence, i.e., an irregular wave forms with the appearance of some rare pulses with high peak powers compared to the average. To characterize this behavior from a statistical point of view, we computed the probability density distribution of the intensity peaks of the circulating intracavity field (using data from the last 1,000 roundtrips). Fig. 6 displays the probability density function versus the pulse peak power. We remark that for an input power of ∼ 1 W (or below), we mostly observe a distribution centered on a unique value of peak power that evolves towards a bimodal distribution. This behavior is directly related to the initial regular pulse train on a finite background whose period becomes doubled due to the growth of intermediary spectral bands (Fig. 4(d)). After that, the temporal pattern is wholly degraded (Fig. 5(d-e)), the corresponding statistics broadens and becomes progressively more tailed (i.e., towards a right-skewed distribution, see Fig. 6(c)). For further increased power, the statistics develops a long tail indicating the presence of rare events with extremely high peak powers (Fig. 6(d)). These events can reach a peak power ten times higher than the average intracavity power. Note that the calculated rogue wave intensity threshold, defined as I RW = 2I SI (I SI being the significant intensity equal to the mean of the upper third of events in the distribution), is I RW = 5.1, 9.2, 17.6, and 34.9, respectively for subplots from Fig. 6. As a consequence, this criterion confirms that extreme events then appear for input power powers ∼ 2 W. To the best of our knowledge, extreme events have never been reported in this monostable regime of a passive Kerr-type cavity. For a better understanding of the mechanism behind this prediction, the rise of the observed spatiotemporal complexity needs to be deeply characterized. However, it is an elusive task to apply many of tools of dynamical systems on Eqs. (1). Hence, for the sake of simplicity and without loss of generality, in what follows, we will use the reduced equation (2). For the accuracy of our prediction based on this model (2), it is worthy to determine the effective parameter after the experimental measurements.
Close to the emission threshold, from Eq. (2), the dimensionless frequency shift of the MI is given by Ω c = √ 2I s − ∆/ (2π) when η = −1. Ω 2 c is then a linear function of the intracavity field with the y-intercept corresponding to −∆. By setting P in the mean intracavity intensity and ω c the physical frequency shift we can infer that ω 2 c = aP in + b with a and b to be determined. P in can be also given in ratio of its value at the MI threshold P th in , such that I s = P in /P th in . Hence, it follows: Consequently, On the other hand, in the weakly nonlinear regime, the amplitude of the MI varies linearly with the distance at the threshold. From the experimental data, we have extracted the ω c and the corresponding magnitude. The result is given in Fig. 7. From Fig.7b the threshold is obtained as P th in = 1.985/3.615 = 0.55 W. By considering the effective length L = 28 m and group velocity dispersion parameter β 2eff = −2 ps 2 km −1 we have calculated the effective parameters corresponding to the configuration of the setup. The most relevant are shown in Table I Interpolation of experimental data extracted from the spectrum. The values of these figures have been extracted from the experimental measured spectra.

III. FROM MODULATION INSTABILITY TO SPATIOTEMPORAL CHAOS
After the validation of expected parameters, we now follow the evolution of the modulation instability from the weak to the unexplored large driving strength regime. Figure 8a shows the evolution of the recorded intracavity power spectrum when the power of the driving beam (i.e., the pump) increases. We first recall the two significant and striking features that have drawn our attention to the evolution of these spectra. The first one is a transition from a purely triangular comb (in a log. scale) to a quasi-periodic-like spectrum formed by the appearance of secondary instability bands in between the primary comb teeth. The second change arises while a continuous spectrum is observed in place of the comblike profile. Later, this continuous spectrum only broadens as the pump power increases. A typical behavior is shown in Fig. 8a obtained from the numerical simulations based on the Lugiato-Lefever equation (2). The bottom panel (Fig. 8b) shows the corresponding density plot of the temporal traces for successive round trips. The key feature of this evolution is the apparent complexity of the dynamics when the spectrum goes beyond the triangular-shape comb (I s 2.0). The analysis of the structural changes in the optical spectrum is performed through the computation of the Lyapunov spectrum (see Appendix A for details) of this reduced model for an increasing pump power. The results are displayed in Fig. 9a. From this figure, we can observe a value of the pump power (I s = 2.0) below which the continuous Lyapunov spectra are composed of only negative exponents. Above this value, the spectra show a region with positive exponents that increases with pump power suggesting the emergence of a chaotic dynamics in the spatiotemporal complexity observed for I s > 2.0 in Fig. 8b. An indicative quantity of the number of degree of freedom needed to describe the chaotic dynamics is the dimension of the strange attractor. An estimator of this dimension, based on the Lyapunov spectrum, has been conjectured by Kaplan and Yorke [51]. The dependance of this Kaplan- Yorke dimension D KY on the control parameter is shown by Fig. 9c, which confirms the increasing of the complexity with the pump power. Besides this growth of the Kaplan-Yorke dimension with the pump power, we must also consider the conjecture that the dimension increases with the size of the system-expansive nature of the spatiotemporal chaos [33,[52][53][54][55]. To focus only on the effect of the pump power, it is necessary to use an intensive measurement of the complexity. This can be done through the Lyapunov dimension density (see Appendix A and B). Indeed, D KY may change linearly with the volume of the system [33,34]. That is, for a 1D system, D KY = ξ −1 δ ∆T where ∆T is the extension of the system and ξ δ represents the Lyapunov dimension density of the system for a fixed value of the pump power. Therefore, it is an intensive quantity that provides an estimation of the extension of the independent subsystems generated by the chaotic dynamics. We have performed a systematic calculation of the Lyapunov dimension density ξ δ and the results are given in Fig. 9b. It appears that the the characteristic range of chaotic fluctuations decreases as the driving power increases.
Considering that the spatiotemporal chaos produces a subsystem of the order of ξ δ , one might wonder how far these subsystems interact. This can be addressed by characterizing the average temporal disorder provided by the equal-time two point correlation length ξ 2 obtained from the exponential decay of the following integral [34,35,38]: The brackets · stand for the average process (see Appendix C for practical computation). For our system, the result is shown in Fig. 9d. We clearly observe a non-monotonic evolution of the correlation time in a subsequent range of pump power. This behavior is well confirmed by corresponding ex-perimental data, suggesting that short-range fluctuations measured by ξ δ may be decoupled from the long-range order measured by ξ 2 [38]. In general, for unidimensional (1D) systems, ξ δ and ξ 2 are proportional, specifically when the later one is computed from the magnitude of the considered field. For now, only one study has reported different evolutions of ξ δ and ξ 2 . However, in that case ξ 2 has been computed from the phase [37]. Hence, the decoupling between ξ δ and ξ 2 measured here from the same dynamical component-the magnitude of the intracavity field-is an unexpected observation. From Figs. 9 b) and c) we can see that the correlation length ξ 2 can be larger up to ten times the Lyapunov density dimension. Therefore, a significant number of the subsystems may be dynamically correlated. Two central questions arise from this observation. Are the spatiotemporal chaotic subsystems contiguous or not and if not what is the characteristic time in between? An appropriate method to answer these questions is to address the dynamics of the different points of the system in terms of laminar (regular) and turbulent (irregular) periods according to a cut-off threshold [56]. If the chaotic subsystems are not contiguous, the dynamics may lead to a fluctuating mixture of turbulent and laminar periods [5,31], the so-called spatiotemporal intermittency. If so, the probability distribution of the laminar periods should be a mixture of the power law (long-range correlation) and the exponential law (short range) [31]. In this mixture distribution, the exponent of the power law (µ) is expected to be insensitive to the input parameter or the value of the cut-off threshold. By contrast, the decay rate (m) of the exponential law depends on both parameters. Since the dependence on the cut-off threshold is only a decaying exponential, the y-intercept (m 0 ) contains the dependence on the control parameter S through I s . We will refer to m 0 as the intrinsic exponential decay.

IV. TRANSITION TO AMPLITUDE TURBULENCE VIA SPATIOTEMPORAL INTERMITTENCY
Before going further, it may be helpful to describe with details how we have detected the laminar and turbulent domain in the temporal traces. This characterization was done following the process explained in Ref. [31]. As reminded by authors in the latter reference, the distinction between laminar and turbulent regions solely in terms of the amplitude does not have a rigorous justification. However, as we will see later, the domain boundaries obtained using the amplitude also correspond to a local change in the periodicity of the patterns. Notice that a rigorous distinction condition should be based on local frequency fluctuations over a roundtrip.
According to the aforementioned process separating our temporal profile in a laminar and turbulent region is done as follows: 1. Detecting hills and valley in the temporal tracer ( Fig. 10a) 2. Computing the amplitude of the local peak to peak amplitude (Fig. 10b) 3. Fixing a threshold below which a local peak is said laminar or turbulent otherwise ( Fig. 10c and 10d).
After the last step of this process, the statistical study of the laminar duration is done. According to the theory, spatiotemporal intermittency is characterized by a power-law decay of the number of laminar time duration, while an exponential decay governs a fully amplitude turbulent regime. However, in case of Ising-type transition the probability density distribution of the laminar time should follow a mixture function as suggested by authors in [31]. Figure 11 shows the loglog plots of the probability density function of laminar region when increasing the pump power. As expected, the slope of the best-linear fit exhibits almost no dependence neither on the detection threshold nor the control parameter. Considering the mean value µ = 4.8 ± 0.4, we were able always to find a portion of the statistics that satisfies the law PDF(x) ∝ x −µ , characteristic of a spatiotemporal intermittent evolution. Elsewhere, there is also some portions of the statistics that may decay exponentially: PDF(x) ∝ e −mx (see Fig. 12a). However, this decay rate depends on c, and the control parameter I s such that m (I s , c) = m 0 (I s ) e −c/c 0 as shown in Fig. 12b. Finally, as it can be seen from Fig. 9e m 0 (I s ) itself presents a power law increase on the control parameter. Indeed, different size of the system show the same linear evolution of m 2 0 upon I s . This figure also allows to identify two critical pointsvalues of the control parameter I s from which m 0 scales as a power law. At the first transition which occurs at I s 2, the distribution of laminar periods follows a power law. Therefore, we conclude that I s 2 , that we set to I s = I ST I s , corresponds to the onset of the spatiotemporal intermittency. Since I s = I ST I s 2 also coincides with the onset of the spatiotemporal chaos, we can conclude that the spatiotemporal intermittency is the route to this chaos [31,56]. The second critical phenomenon which occurs at I s = I T s 3 reveals the onset to the fully developed amplitude turbulence [57], given that the distribution of laminar periods follows an exponential law [56]. For I ST I s < I s < I T s we were able to fit the probability distribution as P(τ) = (Aτ −µ + B) exp(−mτ) (see Fig. 13ad), typical of a phase transition process [31]. In this range of parameters, from the numerical data, we can observe the presence of intensity bursts in the intracavity field temporal profile. Figure 14a shows the typical profile of 100 largest bursts when increasing the pump power. It reveals that the largest peaks are always surrounded by smaller ones located at a position corresponding roughly to ξ 2 . This evidences that 2ξ 2 quantifies the region needed by the bursts to rise and disappear. We have also plotted over these profiles the  experimental two-points correlation function |C (∆τ)| 2 . We emphasize that the correlation function gives an accurate estimation of the temporal profiles of the intensity bursts, specifically inside the turbulent regime. We have then considered the statistical analysis of these intensity peaks in connection with extreme event studies [58]. Figure 9c shows the evolution of the ratio of bursts that can be considered as extreme events using the criterion defined in hydrodynamics [58] and Fig. 14b gives some of the corresponding probability distributions. From these two results, we observe that the emergence of extreme events appears to coincide with the onset of the turbulent regime at I s = I T s . Furthermore, in the vicinity of a local maximum of the correlation length and just before the emergence of the extreme events, we notice the bimodal shape of the probability density function (PDF). This suggests the existence of two main sub-populations in the intensity profile like in the metal-insulator phase transition [59,60].

V. CONCLUSION
In summary, we have shown that the dynamics of the wellknown periodic Turing patterns associated with the triangularshape frequency comb in a Kerr-nonlinear ring cavity can be subject to transition to spatiotemporal chaos. Our findings show that combining the study of different order parameters instead of trying to get the better one might be the sound approach to describe the spatiotemporal complexity. Once established that the whole system can be split into independent subsystems through the computation of the Lyapunov dimension length, the equal-time correlation length has provided the range over which their fluctuations are dynamically connected even not contiguous. Indeed, chaotic subsystems (turbulent) are shown to be separated by coherent subsystems (laminar). By analogy with fluids, the resulting mixture undergoes a transition to turbulence via spatiotemporal intermittency. Our findings also show that non-excitable passive resonators can generate spatiotemporal pulses whose probability distribution exhibits a long tail (i.e., extreme events). Our results highlight how experiments in optics with controllable complexity using technology-driven components can be used to develop understanding of fundamental nonlinear dynamics in dissipative systems. We anticipate applications in establishing links between different branches of nonlinear science for which the occurrence of critical phenomena and extreme events is universal. where d is the dimension of the system and n the number of Lyapunov exponents to be computed. After a time increment dt, the matrix L evolves to L (t 0 + dt) =ÛL (t 0 ) wherê U = e J * dt . Using the modified Gram-Schmidt QR decomposition on L (t 0 + dt), the diagonal elements of R account for the Lyapunov exponentsλ i (i = 1, . . . , n) at time t 0 + dt, that is Repeating this procedure several time, after a large number of iterations N, the Lyapunov exponents can be approximated by ln (R ii (t 0 + kdt)) . (A3)