PHYSICAL REVIEW A 108, 013520 (2023)
Walk-off-induced dissipative breathers and dissipative breather gas in microresonators
A. Villois ,1,2 D. N. Puzyrev,3,4 D. V. Skryabin ,3,4 and M. Onorato 1,5
1Dipartimento di Fisica, Universitá degli Studi di Torino, 10125 Torino, Italy
2School of Mathematics, University of East Anglia, Norwich Research Park, Norwich NR4 7TJ, United Kingdom
3		Department of Physics, University of Bath, Bath BA2 7AY, United Kingdom
4		Centre for Photonics and Photonic Materials, University of Bath, Bath BA2 7AY, United Kingdom
5		INFN, Sezione di Torino, Via Pietro Giuria 1, 10126 Torino, Italy
(Received 13 March 2023; revised 20 June 2023; accepted 12 July 2023; published 28 July 2023)
Dissipative solitons in optical microcavities have attracted significant attention in recent years due to their direct association with the generation of optical frequency combs. Here, we address the problem of dissipative soliton breathers in a microresonator with second-order nonlinearity, operating at the exact phase matching for efficient second-harmonic generation. We elucidate the vital role played by the group-velocity difference between the first- and second-harmonic pulses for the breather existence. We report the dissipative-breather-gas phenomenon, in which multiple breathers propagate randomly in the resonator and collide nearly elastically. Finally, when the breather gas reaches an out-of-equilibrium statistical stationarity, we show how the velocity locking between first and second harmonics is still preserved, naming such phenomena turbulence locking.
DOI: 10.1103/PhysRevA.108.013520

I.		INTRODUCTION
Ultrahigh Q-factor optical cavities, such as whisperinggallery-mode (WGM) microresonators, offer a compact and power-efficient platform for generating various types of spectrally broadband wave forms, including dissipative Kerr solitons [1,2], Turing patterns [3], soliton molecules [4], and soliton crystals [5]. The spectra of these wave forms are associated with optical frequency combs. Recent developments have shown that a similar range of effects happens in microresonators with the second-order, χ (2), nonlinearity [6–10]. Coherent and equidistant microresonator combs can be used for various applications such as precision spectroscopy, optical clocks, and the search for exoplanets [11–13].
One of the prerequisites for efficient operation of a solitonbased microresonator device is the stability of dissipative solitons. Hence, accurate knowledge of the parameters responsible for instabilities is a fundamental problem relevant to technological applications. On the other hand, understanding soliton instabilities enables the experimental realization of breather combs [14–16] and also opens up the possibility of using microresonators as a platform to investigate the interplay between dissipative solitons, multimode chaos, and turbulence [17–19]. The study of dissipative breathers in χ (2) microresonators has so far been restricted to the regime of optical parametric down-conversion [20,21], where dissipative solitons become Hopf unstable, leading to the breather formation. This work will investigate whether this scenario takes place as well in a χ(2) microresonator set for second-harmonic generation (SHG).
Published by the American Physical Society under the terms of the Creative Commons Attribution 4.0 International license. Further distribution of this work must maintain attribution to the author(s) and the published article’s title, journal citation, and DOI.
Assuming that the second harmonic is far from the phasematching conditions, the nonlinear interactions resemble four-wave mixing (cascading regime), and one can expect the existence of dissipative breathers similar to those observed in Kerr microresonators. However, it remains an open question whether dissipative breathers can also exist for exact phase matching (SHG breathers). This question is not limited to microresonators; it also extends to bulk crystals. Recent studies demonstrated breather solutions only in the cascading regime and brought up close analogies with the Akhmediev and Kuznetsov-Ma breathers [22–25]. In this work, we show that dissipative solitons at the phase-matching point (SHG solitons) are linearly stable, provided that the group-velocity walk-off is zero and, therefore, SHG breathers do not exist. Inducing a nonzero walk-off, which is practically unavoidable in experiments, opens up a window of instability not only for the continuous-wave (cw) solution [26] but also for the solitons and, consequently, triggers SHG breathers. In addition, our findings give evidence that the walk-off-induced instability can also lead to the generation of multiple randomly moving and quasielastically interacting dissipative breathers. We refer to such a state as a SHG dissipative breather gas, and we compare it with its SHG dissipative-soliton-gas counterpart. Furthermore, we find that when the light inside the resonator reaches a statistically stationary turbulent state, the locking between the first- and second-harmonic components of the chaotically moving, disappearing, and emerging pulses remains intact. We refer to this phenomenon as the turbulencelocking regime.
II.		MODEL
In this work, we consider a WGM LiNbO3 ring microresonator like the one recently used to experimentally generate a second-harmonic frequency comb [7]. By pumping ordinary polarized light with a cw laser at 1065 nm, it is possible to
achieve natural phase matching and generate extraordinary polarized second-harmonic light. Dispersion relations around the first-harmonic cavity mode mp and second-harmonic cavity mode ms = 2mp are defined as
ωμ = ωmp + μD1p + 2 μ D2p, (1)
ωμ = ωms + μD1s + 2 μ D2s .	(2)
Here, ωmp and ωms denote the cold-cavity resonant frequencies, while μ = 0, ±1, ±2, ... represents the mode number offset with respect to mp,s. The group-velocity dispersion coefficients D2p/2π = -100 kHz and D2s /2π = -200 kHz are both normal, while the repetition rate is D1p/2π = 21 GHz. The walk-off parameter U = (D1p - D1s) is of the order of 1 GHz at the phase-matching point, 2ωmp = ωms . Throughout this work we treat the walk-off as a free parameter. The equations governing the evolution of the first- and second- harmonic envelopes ψp,s are given by
i∂tψp= δp-iD1p∂θ-12D2p∂θ2 ψp - i κ2p ψp+h - γpψsψp∗,
1κ
i∂t ψs =	δs	- iD1s∂θ	- D2s∂2	ψs-i sψs
t s	s 1s θ 2 2s θ	s 2	s

γs ψp2;

(3)
(4)

a formal derivation can be found in [27]. Here, δ = ω - ^ , p	mp represents the detuning of the cw pump laser frequency ^ from the cold-cavity resonance frequency ωp . Similarly, δs = ωms - 2^ denotes the detuning of the second harmonic. From the detuning definitions it is possible to rewrite the second- harmonic detuning as δs = 2δp - ε, where ε = 2ωmp - ωms is the frequency-mismatching parameter. The linewidths are given by κp/2π = 1 MHz and κs /2π = 4 MHz, setting resonator finesse F = D1p/κp ∼ 104. The nonlinear coefficients are defined as γp,s/2π = 300 MHz W-1/2. The power of the cw laser W is related to the pump parameter h = iκp/2 η/πFW, where η is the coupling coefficient, η = 0.5.
III.		QUADRATIC DISSIPATIVE SOLITONS
AND WALK-OFF-INDUCED INSTABILITY
One of the main effects that the walk-off has on dissipative quadratic solitons is the fact that localized structures cannot travel with linear group velocity, and their velocity is selected by the dissipative effects. The value of the velocity at which the first- and second-harmonic solitons are locked together is such that the linear momentum M is conserved in time, i.e.,
M =-(κ M +κM) =0,	(5) dt	p p s s ,
where
M =Mp+Ms,	(6)
π
(ψp∗∂θ ψp - c.c.)dθ, π π
(ψs∗∂θψs - c.c.)dθ.
π

(7)
(8)

The functional form of (7) and (8) is such that Mp + Ms is conserved in time in the Hamiltonian limit, κ = 0 (see [28]). p,s
Note that the momentum is not, in general, a conserved quantity, given the presence of dissipation; however, the system can support traveling-wave solutions moving with the common velocity V if the right-hand side (rhs) of Eq. (5) vanishes.
A natural question that arises is how to estimate the value of this velocity. In order to answer this, it is instructive to consider the case when κp,s = 0. At such a Hamiltonian limit, soliton solutions are not limited to a single value of velocity locking, and V can span a continuous range of values [28]. Since an analytical solution for second-harmonic χ2 cavity solitons does not exist to our knowledge, it is useful to start identifying the range of possible velocities at which a solitary traveling wave can propagate in a microresonator. To do so, one can rely on the so-called band-gap analysis already introduced in [29].
The idea behind the band-gap analysis is as follows: a localized or solitary wave moving with velocity V in an optical cavity can exist on top of a homogeneous cw solution ψp0,s only if the condition
ω
V ^= μ (9) μ
is true for each μ in the cavity, where ωμ represents the dispersion relation for weakly nonlinear waves propagating on top of the cw solution. Such a dispersion relation can be found by performing a standard Bogoliubov–de Gennes analysis. In order to emphasize the role of the walk-off parameter, we work in a frame of reference moving with angular velocity D1p. We can now substitute the ansatz
ψp,s = ψp0,s + up,seiμθ+λt +v∗p,se-iμθ+λ∗t, (10) with λ ∈ C, in Eqs. (3) and (4) and solve the eigenvalue problem
MA = λA (11)
arising after linearizing with respect to the low-amplitude
perturbations |u | and |v |, where A = (u , p,s	p,s	p	v p , us , vs )T and 0⎤	
		0 -γp ψs	∗ γp^ψp0^		
M=⎢⎢	Lγpp^-ψs0i^κ∗p	-L∗ - iκ pp	0	γpψ0		
	-2gsψ0	0	Ls - iκs	0	⎥⎥,	
	0	∗ 2gs^ψp0^	0	∗
-Ls - iκs		
(12)

with κp,s = 0 and the following definitions: Lp ≡ (δp + 1/2D2 p μ2) and Ls ≡ (δs + 1/2d2s μ2 + Uk). Figure 1 shows the four phase-speed branches ωμ/μ for low-amplitude waves numerically computed in the frame of reference moving with angular velocity D1p . Note that ωμ corresponds to Im(λ). Dashed lines in the plot highlight a band gap: a range of velocities such that low-amplitude waves cannot propagate and, hence, it is possible to prevent energy transferring from solitonic solutions to low-amplitude waves.
We can now focus on the role of the walk-off parameter in the band-gap analysis. One of the main effects of the walk-off is presented in Fig. 1(b), which shows how nonzero walk-off in the system modifies the dispersion of linear waves, leading to the narrowing of the band gap. As we increase the walk-off

FIG. 1. Numerically evaluated phase speed for low-amplitude waves on the top lower branch cw solution obtained for the following parameters: δp/2π =-26.64MHz, δs = 2δp, and W = 260 mW. In (a) the walk-off parameter U is zero, while in (b) U/2π =-5 MHz.

parameter beyond a critical value Uc , the band gap completely closes, thereby impeding the existence of soliton solutions.
The presence of walk-off also has a significant impact on the behavior of SHG solitons, specifically with respect to the velocity at which the first- and second-harmonic solitons become locked together. As already mentioned, by reintroducing losses in the system, first- and second-harmonic solitons will travel at a speed such that the rhs of Eq. (5) vanishes. This value will be, in general, selected by the values of dissipation parameters κp,s, as discussed in the Appendix, and it will lie within the range estimated using band-gap analysis. The specific value of velocity locking can be evaluated numerically, with a velocity-selective Newton-Raphson method (see red circles in Fig. 9 in the Appendix). Alternatively, it can be analytically shown that for small values of walk-off, the locking velocity V depends linearly on U :
V=αU. (13)
A derivation of the α coefficient can be found in the Appendix [see Eq. (A6)], where a comparison between the velocity-locking values obtained numerically and analytically is shown in Fig. 10(b).
IV.		QUADRATIC DISSIPATIVE BREATHERS
After presenting the influence of the walk-off on the velocity and the existence of SHG soliton solutions, we now examine its impact on their stability and how it contributes to the transition from solitons to breathers. It is well known that breather solutions can exist in the integrable 1D Nonlinear Schrödinger equation (NLS) in the form of time- [23,24] and space- [25] periodic solutions. In the case of quadratic nonlinearity, similar solutions have been found only in the cascading regime limit, where four-wave mixing is dominant [22]. In our work, such a limit would correspond to the case with δs ^ δp. For Kerr optical cavities, where driving and dissipations are added to the standard NLS equation (also known as the Lugiato-Lefever equation [30]), breathers can exist when solitons undergo a Hopf instability [20]. Numerical and experimental results concerning dissipative breathers can be found in [16,31]. The main problem concerning finding SHG breathers in the phase-matching regime is related to the stability of SHG solitons. In order to study the linear stability

FIG. 2. Existence region for soliton solutions for different values of detuning and cw-laser input power (a) in the absence of walk-off and (b) in the presence of walk-of. The existence region is studied up to -30 MHz, although it can extend beyond that value. The stability and instability regions are denoted in green and orange, respectively. The bistability line denotes where the homogeneous cw solution becomes multivalued.

of soliton solutions, one can modify Eq. (12) by substituting the soliton profile instead of the cw background solution; considering low-amplitude perturbations to be generic functions of θ, up,s(θ), and vp,s (θ); and redefining Lp ≡ (δp - 1/2D2p∂θ2 - iV∂θ) and Ls ≡ [δs - 1/2D2s∂θ2 - i(V -U)∂θ]. Different from dissipative solitons in the Kerr case, SHG solitons appear to be stable for all values of detuning and/or driving power [see the green region in Fig. 2(a)], as long as the walk-off parameter is neglected. Note that moving away from the ideal scenario of perfect group-velocity matching (zero walk-off) makes it harder to find SHG solitons due to the closure of the band gap. Only recently has the existence of quadratic dissipative solitons (DSs) with large walk-off been proven, but they are limited to the cascading regime [10]. However, it is possible to alter the stability of SHG solitons by considering a value of the walk-off parameter large enough to trigger a Hopf instability [see the orange region in Fig. 2(b)] but less than the critical value Uc responsible for the closure of the existence gap. By selecting detuning or cw input power within the orange region in Fig. 2(b), one can investigate the effects of such instability in the dynamical evolution of the SHG soliton. Figure 3(a) shows how the instability is responsible for the soliton-breather transition. The plot shows the dynamics of the first-harmonic component of the SHG soliton in a comoving frame of reference. This choice was made to emphasize how, after a few thousand round-trip times, the soliton starts radiating dispersive waves and stabilizes at a new velocity. The system reaches a metastable state characterized by random velocity shifts, as highlighted by the white circle in Fig. 3(a). Metastability has been tested numerically for up to 4 × 106 round-trip times. Given the finite size of the system, the emitted dispersive waves cannot escape, and they keep traveling within the cavity. We stress that this is not a numerical artifact, but an intrinsic property of microresonators. This process will trigger a continuous emission and reabsorption of the waves causing the breathing of the soliton [32]. The coherence of such a state can be observed in the space-time Fourier analysis of the dynamics. Specifically, from Fig. 3(b) it is possible to see how the dissipative breather is formed by the

FIG. 3. (a) Dynamical evolution of the unstable soliton solution. The detuning parameters and cw-laser power are δp/2π =-26.64MHz, δs = 2δp, W = 870 mW, and V =-1.62MHz. (b) Space-time Fourier transform of breather dynamics within the time interval [2.5, 2.6] × 106 round-trip times. The weakly nonlinear dispersion relation is superimposed by the white dashed lines. (c) rf spectrum for the first harmonic. The inset in (c) shows the first- and second-harmonic spectra (black and orange lines, respectively) of the dissipative breathers after 2.6 × 106 round-trip times.

superposition of several coherent structures, all moving with the same velocity. The crossing of the breather signal with the weakly nonlinear dispersion [see the white dashed lines in Fig. 3(b)] causes the resonant emission of dispersive waves with characteristic wave number μr ∼ 22 [see the white circle in Fig. 3(b)]. Note that the weakly nonlinear dispersion is the same as that plotted in Fig. 1(b) but in a frame of reference comoving with the soliton. Spectra evaluated after 2.6 × 106 round-trip times are shown in the inset of Fig. 3(c) for both first and second harmonics. In the first-harmonic spectrum, in addition to the main peak associated with the cw pumping, a smaller secondary peak is also present at μ = 22. In order to evaluate the oscillation period of the SHG breather, one can evaluate the dominant frequency in the system through the so-called rf spectrum, evaluated as Ft [^ψ |ψ ^θ], where Ft represents the Fourier transform with respect to time and ^ψ |ψ ^θ = |ψ|2dθ. Such a plot is shown in Fig. 3(c), where a peak is observed at 0.01 GHz. This value is remarkably close to half of the gap size, as pointed out in a recent study on dissipative Kerr breathers [33].
V.		DISSIPATIVE BREATHER GAS
The effects of the walk-off-induced instability on soliton dynamics can be further studied by examining the transition from a SHG soliton to a turbulent state This transition can be achieved by increasing the growth rate Re(λ) of the linearly unstable eigenstates responsible for the soliton instability, for instance, by maintaining the same detuning as shown in Fig. 3 while increasing the input cw power W . Figure 4(a) displays the time evolution of the first-harmonic component of an unstable SHG soliton inside the resonator. Due to the walk- off-induced instability, the soliton first transforms into an SHG breather, then into a state characterized by an increasing number of breathers, before finally settling into a turbulent state similar to the Kerr spatiotemporal chaos observed in [5]. This state is referred to as SHG dissipative breather gas, a term adopted from the integrable-turbulence community [34]. Unlike integrable systems in which solitons or breathers can interact only elastically with each other, our system can reach an out-of-equilibrium stationarity characterized by a multitude of pulses continuously emitting dispersive waves. The loss of coherence typical of a turbulent state is evident in both Figs. 4(b) and 4(c), where there is no longer a clear signal associated with a coherent structure or a well-defined peak in the rf spectrum.
A question that arises naturally concerns the difference between the SHG dissipative breather and soliton gas. As already pointed out, different from dissipative solitons in the Kerr case, SHG solitons are not affected by Hopf-type instability; as long as the walk-off parameter is negligible, this result is independent of the detuning value and the initial power. Therefore, we can exploit this feature to generate a gas of SHG dissipative solitons. To engineer such a state, we start by considering an initial condition characterized by 15 solitons in both harmonics, which we will refer to as a soliton crystal. Figures 5(a) and 5(b) display a portion of the crystal obtained numerically using the Newton-Raphson method. These plots allow for a visual comparison between the crystal (see the black line) and the single-soliton solution (see the orange line) for the first-harmonic component. Figure 5(c) demonstrates that the velocity locking of the crystal remains consistent with that of the single-SHG-soliton solution for different values of the walk-off parameter. Additionally, Fig. 5(c) confirms the linear trend as predicted in Eq. (13). A substantial difference between the crystal and the single-SHG-soliton solution concerns their linear stability. In Fig. 5(d) we plot the dimensionless value of the total momentum for the system

FIG. 4. (a) Dynamical evolution of the unstable soliton solution. The detuning parameters and cw-laser power are δp/2π = -26.64 MHz, δs = 2δp, W = 910 mW, and V =-1.63MHz. (b) Space-time Fourier transform within the time interval [3, 5] × 104 round-trip times. (c) rf spectrum evaluated as in Fig. 3(c).

with respect to different values of walk-off. Despite having the same velocity, the crystal exhibits higher momentum due to its higher power intensity. Dashed and solid lines in Fig. 5(d) denote unstable and stable solutions, respectively. More specifically, the crystal is unstable for small walk-off values, different from the single soliton (see the black and orange lines, respectively). By increasing the walk-off value, the crystal stabilizes and then destabilizes again, experiencing the same type of walk-off-induced instability experienced by the single SHG soliton which led to the formation ofa SHG breather, as discussed in the previous section. The different

FIG. 5. Space profile of the stationary DS solution and soliton crystal (orange and black lines, respectively) in the case of (a) zero walk-off and (b) U/2π =-5MHz. (c) Plot of velocity locking for both a single DS and soliton crystal versus walk-off. (d) Total linear momentum in dimensionless units versus walk-off. Dashed lines highlight unstable solutions. Fourier spectrum of the soliton crystal (black lines) superimposed on the Fourier spectrum of unstable crystal eigenstates (green lines) for (e) zero walk-off and (f) higher walk-off, |U |/2π = 5 MHz. For each panel detuning and cw input power are the same as those used in Fig. 4.

natures of the two instabilities affecting the crystal can be understood by considering the spectrum associated with the most linearly unstable eigenstates of the crystal [see Figs. 5(e) and 5(f)]. In the case of a small walk-off, the instability leads to the formation of a secondary comb with a different periodicity, while for large walk-off values, the instability leads to the creation of sidebands around the primary comb. As a final remark, the instability for small values of walk-off appears in the form of a single positive real eigenvalue λ, while for large walk-off it appears in the form of a pair of complex*conjugate eigenvalues with positive real parts [35].

FIG. 6. (a) Collection of spectra from the dynamical evolution of stationary crystal states for different values of the walk-off parameter. Each spectrum is measured after 5 × 105 round-trip times. Time evolutions of the soliton crystal for (b) |U |/2π = 1 MHz and (c) |U |/2π = 5 MHz. For each panel detuning and cw input power are the same as those used in Fig. 4.

We will now study the long-term effects of these two different types of instability on the temporal dynamics of the crystal. Figure 6(a) shows the first-harmonic component of a collection of spectra obtained from the dynamical evolution of stationary crystal states for different values of the walk-off parameter. Each spectrum is measured after 5 × 105 roundtrip times. As predicted from the linear-stability analysis [see dashed lines in Fig. 5(d)], the crystals exhibit instabilities in two separate walk-off regions, denoted by I and II in Fig. 6(a). The boundaries of these regions are highlighted by white dashed lines. Examples of the dynamical evolution of the first-harmonic component of the crystal for each unstable region are presented in Figs. 6(b) and 6(c). The instability present in region I causes the breaking of the crystal into its individual constituents, which are, however, stable according to the linear stability analysis presented in Fig. 2(a). For this reason, in Fig. 6(a) itis still possible to identify the trajectories of the individual solitons and their interaction. This type of dynamics is, in contrast, not present in the dynamics associated with region II, where the crystal seems to transition into a state similar to the one already reported in Fig. 4(a). For this reason, we are led to identify the type of turbulence characterizing region I as a SHG dissipative soliton gas, while the one characterizing region II is a SHG dissipative breather gas. To compare these two different types of dissipative gases in a more quantitative way, it is possible to study the probability distribution function (PDF) for the field amplitudes. To do this, we repeat the simulations presented in Figs. 6(b) and 6(c) while considering the same initial conditions but different random noises.
To construct a PDF, we follow a process that involves generating a histogram. In this histogram, each bin value, denoted as vi , is determined using the formula vi = ciN/wi , where ci represents the number of elements within the bin, N denotes the total number of data points in the histogram, and wi represents the bin width, which is set at 1 mW. From the PDF associated with region I, plotted in Fig. 7(a), one can see how the first harmonic has a peak around 0.05 W, while the second harmonic has a peak around 0.03 W (see the dark gray and light gray histograms, respectively). Such peaks correspond to the amplitude of the single-SHG-soliton solution [see Fig. 5(a)], justifying the assumption of a SHG dissipative soliton gas. Note the presence of a peak for low

FIG. 7. Probability density function or the wave amplitude for many realizations of (a) the dissipative soliton gas when walk-off is zero and (b) the dissipative breather gas when walk-off is U /2π = -5 MHz. Dark gray and light gray histograms are associated with the first and second harmonics, respectively. (c) Ensemble average of time evolution of first-harmonic peak counts higher than 0.04 W.

FIG. 8. Plot of the ratio between ensemble average of the linear momentum for first and second harmonic over different times. tch stand for the time when the system fully transits into a turbulent state. tRT denotes the round trip time defined as the inverse of the free spectral range of the resonator.

power too, which can be associated with the homogeneous cw background on which the solitons lie. In contrast, the PDF for region II, plotted in Fig. 7(b), appears to be flat for both harmonics. The absence of well-defined peaks is due to the multitude of low-amplitude dispersive waves constantly emitted by the breathing mechanism and by the fact that each SHG breather has an amplitude which oscillates between |ψp |2 ∈ [0.2, 0.8] W and |ψs|2 ∈ [0.1, 0.5] W for the first and second harmonics, respectively. To make sure that the difference in the PDFs is not simply a consequence of a different number of high-amplitude pulses in the system, one can measure the number of peaks with an amplitude higher than a fixed threshold. Figure 7(c) displays the temporal evolution of the peak count for the first-harmonic dynamics in regions I and II. The peak count is averaged over 25 different realizations characterized by different white noises in the initial condition. A threshold value of 0.04 Wis selected, slightly lower than the peak power of the single-SHG-soliton solution. Figure 7(c) reveals that during the initial stage of the dynamics, both regions exhibit an equal number of peaks. However, due to the crystal instability, a sudden increase in peak count is observed specifically in the dynamics associated with region II. This phenomenon can be attributed to the influence of walk-off, resulting in an asymmetric shape within the solitons forming the crystal. This asymmetry is manifested by a minor amplitude bump [as depicted in Fig. 5(b)] on the right side of each peak within the crystal. Over time, the Hopf instability amplifies this small bump, giving rise to the generation of multiple pulses within the resonator. Consequently, for a brief period, the number of these pulses seems to be twice the count of peaks characterizing the crystal. However, after this transitional phase, both the SHG dissipative breather and soliton gas settle around an average peak count, which becomes comparable.
VI.		TURBULENCE LOCKING
As a final remark, we will focus on the SHG dissipative breather gas, showing how the system reaches a statistically stationary equilibrium. This can be seen by measuring the ensemble average of the rhs of Eq. (5) over many realizations, each characterized by a different initial random noise (see

FIG. 9. (a) Plot of the linear momentum for the first (black line) and second (orange line) harmonics with respect to different soliton velocities in the absence of walk-off. (b) Ratio of linear momenta vs soliton velocity for nonzero walk-off, U /2π = 5 MHz. Red circles indicate the numerically evaluated velocity selection for dissipation ratios κs /κp = 4 and κs /κp = 1. The detuning parameters and cw-laser power are the same as those used in Fig. 1(b).

Fig. 8). For this system, statistical stationarity implies that
^Mp^ =-κs
^M ^	κ ,
sp

(14)

where ^·^ stands for the ensemble average over many realizations. As pointed out in the Appendix, if there is no walk-off, the single-SHG solution is quiescent, V = 0. This means that both Mp and Ms are equal to zero, as shown in Fig. 9(a). The same idea can be applied to the case of statistical stationarity. If U is zero, then the ensemble averages ^Mp ^ and ^Ms ^ are also zero. On the other hand, when U is not zero, both ^Mp^ and ^Ms^ are nonzero, and their ratio is fixed by the dissipation value. In our system, where κp/κs = 4, Fig. 8 shows how the ratio Mp /Ms oscillates around the value -4. This observation suggests that despite the chaotic state, the first and second harmonics remain locked together and move with a nonzero average velocity. We refer to this phenomenon as turbulence locking.
VII.		CONCLUSIONS
We presented a study of the impact of walk-off on the stability of a second-harmonic-generated soliton comb in a realistic quadratic WGM microresonator. We showed that SHG breathers can be obtained in microresonator cavities away from the cascading regime. The soliton-breather transition can, indeed, be achieved away from the cascading regime due to the instability caused by the walk-off between the firstand second-harmonic light fields. Soliton instability was also shown to be responsible for the transition into a turbulent regime where a multitude of breathers coexist with dispersive waves. A study of the PDF for the wave amplitude confirmed the presence of a walk-off-induced SHG dissipative breather gas. The dissipative breather gas was then compared to its SHG dissipative soliton counterpart. Further analysis of the statistical average of the linear momentum showed the presence of locking between the turbulent state in the first and second harmonics. Finally, a semianalytical method was presented in order to estimate the velocity at which second- harmonic-induced combs move in the presence of walk-off. In parallel with this work, experimental observations of the SHG

FIG. 10. (a) First-harmonic components of neutral eigenvectors of matrix (A5) and its transpose. Each eigenvector is normalized to unity. (b) Relative error for analytical prediction for the soliton locking velocity.

(A3)

for the system of equations
^p-iV∂θ - 12d2p∂θ2 Ap-gpAsA∗p+h=0,
^^s - i(V -U)∂θ - 12d2s∂θ2^As-gsA2p=0.	(A2)
We would like to stress that such a procedure can be easily applied in the Hamiltonian limit since for each value of the walk-off parameter soliton solutions admit a continuous range of velocities. From Figs. 9(a) and 9(b) it is possible to see how the linear momenta for each harmonic Mp,s vary for different velocities. Even if we are dealing with a Hamiltonian limit, we can still extract further information concerning the locking velocity in the dissipative case by noticing that Eq. (5) forces moments to have opposite signs:
M (V )	κ
p	-s
Ms(V) κp .

	-kp - gp^As0^i	Lp+gpAs0 r	i gp ^Ap^	r -gp Ap		
M=	⎢⎢ -Lpgp^As0^r
⎢	^ 0^i
⎢	-2gs Ap	-kp + gp^As0^i
-2gs^A0p^r	gp^A0p^r
-ks	i gp Ap
Ls	⎥⎥,	(A5)	
	2gs^A0p^r	-2gs^A0p^i	-Ls	-ks	⎦		

where we used the following definitions: Lp ≡ (^p - equation [28]: 1/2d2p∂θ2) and Ls ≡ (^s - 1/2d2s∂θ2). By noting that η ≡
∂θA0 corresponds to the neutral eigenvector of L, such that
Mη = λη with λ = 0, it is possible to derive the following

V = U ^φ|Pη^, ^φ|η^

(A6)

breathers were reported in microresonators made of thin-film lithium niobate [36].
ACKNOWLEDGMENTS
M.O. and A.V. acknowledge the support provided by the Simon Collaboration on Wave Turbulence (Award No. 651471) and the Ministero dell’Universitá e della Ricerca under the PRIN program (Project No. 2020X4T57A). D.V.S. acknowledges support from the Royal Society (SIF/R2/222029).
APPENDIX: VELOCITY LOCKING
In this Appendix we present an analysis of the velocity at which solitons are locked together. We start by evaluating numerically soliton solutions in the Hamiltonian limit. For numerical reasons it is easier to cast Eqs. (3) and (4) in dimensionless units by performing the following transformations: Ap,s = ψp,s/ W, τ = D1pt, d2p,2s = D2p,2s/D1p, ^p,s = δp,s/D1p, kp,s = κp,s/D1p, gp,s = γp,s√W/D1p, and h = H/ W. In order to find soliton solutions moving with velocity V = V/D1p in the presence of walk-off U = U/D1p we will make use of a root-finding Newton-Raphson method to find a stationary solution in a moving frame of reference,
θ → θ -(D1p-V)t,	(A1)
From Fig. 9(a) it clear that in the case of zero walk-off Eq. (A3) is never satisfied unless both Mp and Ms vanish, which is the case only when V = 0. In contrast, when U ^= 0, Eq. (A3) can be verified for a range of velocities, such as the one highlighted in yellow in Fig. 9(b). From Eq. (A3) we can see how the relation between dissipation coefficients allows us to select a specific value of velocity locking. For example, Fig. 9(b) shows the velocity at which the soliton will move considering values of κs/κp = 1 and 4. The estimated velocity is remarkably close to the one obtained numerically by making use of the velocity-selective Newton-Raphson method [see red circles in Fig. 9(b)].
We now reintroduce dissipation in the system and present a semianalytical method to estimate the value of the velocity locking. Such a method requires only knowledge of the dissipative soliton solution A0p,s = Ap,s (V = 0, U = 0), which is, in the case of zero walk-off, easy to evaluate numerically since, as shown already, the locking velocity vanishes too. Starting from the case when U = 0 and V = 0, we can assume that by introducing a small walk-off in the system, U ∼ ^ , the solutions will change accordingly, Ap,s = A0p,s + ^ (Arp,s + iAIp,s), resulting in the first and second harmonics moving at nonzero locking velocity V ∼ ^ . Considering Eqs. (A2) under such conditions and separating the real and imaginary parts, we can write the first order in ^ as
MA+V∂θA0 = U ∂θ A(0),	(A4)
where each column vector contains real and imaginary parts of the first- and second-harmonic components such that A = (Arp, Aip, Ars, Ais)T. Note that the matrix M is given by
where Pη = (0, 0,∂Ars,∂Ais)T, ^·|·^ defines the scalar product, and φ corresponds to the neutral eigenvector of M†. Figure 10(a) shows the first-harmonic components of η and φ. Despite the fact that Eq. (A6) was derived in the
limit of small walk-off, our semianalytical formula provides excellent results for all possible values of walk-off where solitons exists; see Fig. 10(b), where the percent error is plotted as a function of the value of the walk-off parameter.

[1]		T. Herr, V. Brasch, J. D. Jost, C. Y. Wang, N. M. Kondratiev, M. L. Gorodetsky, and T. J. Kippenberg, Nature Photon. 8, 145 (2014).
[2]		A. Fülöp, M. Mazur, A. Lorences-Riesgo, Ó. B. Helgason, P. H. Wang, Y. Xuan, D. E. Leaird, M. Qi, P. A. Andrekson, A. M. Weiner, and V. Torres-Company, Nat. Commun. 9, 1598 (2018).
[3]		S.-W. Huang, J. Yang, S.-H. Yang, M. Yu, D.-L. Kwong, T. Zelevinsky, M. Jarrahi, and C. W. Wong, Phys. Rev. X 7, 041002 (2017).
[4]		W. Weng, R. Bouchand, E. Lucas, E. Obrzud, T. Herr, and T. J. Kippenberg, Nat. Commun. 11, 2402 (2020).
[5]		M. Karpov, M. H. Pfeiffer, H. Guo, W. Weng, J. Liu, and T. J. Kippenberg, Nat. Phys. 15, 1071 (2019).
[6]		E. Obrzud, I. Breunig, J. Szabados, K. Buse, S. J. Herr, S. Lecomte, T. Herr, V. Brasch, and Y. Jia, Opt. Lett. 43, 5745 (2018).
[7]		J. Szabados, D. N. Puzyrev, Y. Minet, L. Reis, K. Buse, A. Villois, D. V. Skryabin, and I. Breunig, Phys. Rev. Lett. 124, 203902 (2020).
[8]		A. Villois, N. Kondratiev, I. Breunig, D. N. Puzyrev, and D. V. Skryabin, Opt. Lett. 44, 4443 (2019).
[9]		N. Amiune, Z. Fan, Z. Fan, V. V. Pankratov, V. V. Pankratov, D. N. Puzyrev, D. N. Puzyrev, D. V. Skryabin, D. V. Skryabin, D. V. Skryabin, K. T. Zawilski, P. G. Schunemann, I. Breunig, I. Breunig, and I. Breunig, Opt. Express 31, 907 (2023).
[10]		D. N. Puzyrev, V. V. Pankratov, A. Villois, and D. V. Skryabin, Phys. Rev. A 104, 013520 (2021).
[11]		K. Beloy et al., Nature (London) 591, 564 (2021).
[12]		M. G. Suh, X. Yi, Y. H. Lai, S. Leifer, I. S. Grudinin, G. Vasisht, E. C. Martin, M. P. Fitzgerald, G. Doppmann, J. Wang, D. Mawet, S. B. Papp, S. A. Diddams, C. Beichman, and K. Vahala, Nature Photon. 13, 25 (2019).
[13]		T. J. Kippenberg, A. L. Gaeta, M. Lipson, and M. L. Gorodetsky, Science 361, eaan8083 (2018).
[14]		H. Guo, E. Lucas, M. H. P. Pfeiffer, M. Karpov, M. Anderson, J. Liu, M. Geiselmann, J. D. Jost, and T. J. Kippenberg, Phys. Rev. X7, 041055 (2017).
[15]		A. A. Afridi, H. Weng, H. Weng, J. Li, J. Liu, M. McDermott, Q. Lu, W. Guo, W. Guo, and J. F. Donegan, Opt. Continuum 1, 42 (2022).
[16]		M. Yu, J. K. Jang, Y. Okawachi, A. G. Griffith, K. Luke, S. A. Miller, X. Ji, M. Lipson, and A. L. Gaeta, Nat. Commun. 8, 14569 (2017).
[17]		M. H. Anderson, A. Tikan, A. Tusnin, J. Riemensberger, A. Davydova, R. N. Wang, and T. J. Kippenberg, Phys. Rev. X 13, 011040 (2023).
[18]		H. J. Chen, Q. X. Ji, H. Wang, Q. F. Yang, Q. T. Cao, Q. Gong, X. Yi, and Y. F. Xiao, Nat. Commun. 11, 2336 (2020).
[19]		S. Coulibaly, M. Taki, A. Bendahmane, G. Millot, B. Kibler, and M. G. Clerc, Phys. Rev. X9, 011054 (2019).
[20]		D. V. Skryabin, Phys. Rev. E 60, R3508(R) (1999).
[21]		P. Parra-Rivas, C. Mas Arabí, and F. Leo, Phys. Rev. Res. 4, 013044 (2022).
[22]		F. Baronio, Opt. Lett. 42, 1756 (2017).
[23]		E. A. Kuznetsov, Akad. Nauk SSSR Doklady 236, 575 (1977).
[24]		Y.-C. Ma, Stud. Appl. Math. 60, 43 (1979).
[25]		N. N. Akhmediev and V. I. Korneev, Theor. Math. Phys. 69, 1089 (1986).
[26]		F. Leo, T. Hansson, I. Ricciardi, M. De Rosa, S. Coen, S. Wabnitz, and M. Erkintalo, Phys. Rev. Lett. 116, 033901 (2016).
[27]		D. V. Skryabin, J. Opt. Soc. Am. B 37, 2604 (2020).
[28]		D. V. Skryabin and A. R. Champneys, Phys. Rev. E 63, 066610 (2001).
[29]		A. Villois and D. Skryabin, Opt. Express 27, 7098 (2019).
[30]		L. A. Lugiato and R. Lefever, Phys. Rev. Lett. 58, 2209 (1987).
[31]		M. Johansson, V. E. Lobanov, and D. V. Skryabin, Phys. Rev. Res. 1, 033196 (2019).
[32]		D. V. Skryabin and W. J. Firth, Opt. Lett. 24, 1056 (1999).
[33]		D. N. Puzyrev and D. V. Skryabin, Opt. Express 30, 39396 (2022).
[34]		G. El and A. Tovbis, Phys. Rev. E 101, 052207 (2020).
[35]		S. Wang, Z. Qi, T. F. Carruthers, A. M. Weiner, C. R. Menyuk, J. Jaramillo-Villegas, G. D’Aguanno, and M. Qi, Optica 6, 1220 (2019).
[36]		J. Lu, D. N. Puzyrev, V. V. Pankratov, D. V. Skryabin, F. Yang, Z. Gong, J. B. Surya, and H. X. Tang, Nat. Commun. 14, 2798 (2023).