Self-Starting Soliton-Comb Regimes in χ(2) Microresonators

1. Introduction

Generation of broad frequency combs, i.e. of a large sequences of equidistant coherent light lines, in monochromatically pumped high-Q Kerr microresonators [1,2,3,4,5] can be regarded as one of the main discoveries in nonlinear optics during the last two decades. Remarkably, the spectral positions of these lines do not generally coincide with (are not close to) the discrete resonator modes because of dispersion of the latter. This is why interpretation of the combs as a sequential parametric excitation of new resonator modes is rather fruitless. The right physical picture of the frequency comb in question is proven to be generation of a spatially narrow soliton propagating along the resonator rim with a constant velocity $v_0$ [5]. Fourier expansion of the light field employing its periodicity over the circumference ensures equidistance of the light lines with the frequency spacing $v_0/R$, where R is the major resonator radius. Thus, the comb and soliton are two physical terms expressing different properties of the same nonlinear light state. Most probably, the soliton nature of broad frequency combs has no alternative.

The soliton-comb concept is nowadays generally accepted for ${\chi }^{\left(3\right)}$ resonators. It is applicable not only to glass based Kerr microresonators, but also to resonators made of ${\chi }^{\left(3\right)}$ crystals [5,6,7,8,9,10]. One can speak of the corresponding comb technology. Amount of publications on the ${\chi }^{\left(3\right)}$ combs can be estimated as $\sim {10}^2$. Note that the dissipative solitons in question are essentially different from (and more complicated then) the so-called conservative solitons associated typically with the nonlinear Schrödinger equation [11,12]. The dissipative solitons balance not only the dispersion broadening and nonlinear compression, but also the gain owing to external pumping and internal resonator losses.

Transfer of the soliton-comb concept to the ${\chi }^{\left(2\right)}$ microresonators, possessing the lower second-order nonlinearity, promises further lowering of the pump power, entering new spectral ranges, and new attractive operation regimes, see also below. This transfer cannot be straightforward because of specific and important features inherent in the ${\chi }^{\left(2\right)}$ nonlinearity [13]. In contrast to the ${\chi }^{\left(3\right)}$ case, elementary nonlinear processes involve not four but three light quanta. These processes are the thresholdless generation of the second harmonic (or sum frequencies) and the inverse threshold process of parametric oscillation. The ${\chi }^{\left(2\right)}$ comb generation can be viewed as a cascade of these processes. Running of them requires employment of the phase matching (PM) for the first harmonic (FH) and second harmonic (SH) resonator modes. The PM conditions in ${\chi }^{\left(2\right)}$ microresonators can in turn be either natural, relevant to certain light wavelengths, or artificial, using the so-called radial poling [14,15,16,17,18], see also below.

As soon as mutually coupled FH and SH modes are distinguished, we must deal with duality of the ${\chi }^{\left(2\right)}$ solitons and combs – each of them must consist of FH and SH components, see also Figure 1. The presence of sub-combs in the FH and SH frequency ranges can definitely be regarded as a positive feature. Also, two different pumping schemes (pumping into FH or SH resonator modes) have to be envisaged. These FH and SH pumping cases are substantially different. The presence of temporal walk-off between the FH and SH modes, which can be viewed as the difference of the corresponding group velocities, is also an important feature of the ${\chi }^{\left(2\right)}$ case. The effect of the walk-off on the soliton-comb generation is generally negative – it has to be minimized when possible. Furthermore, the presence of two dispersion coefficients, relevant to the FH and SH spectral ranges, has to be taken into account. Generally, the ${\chi }^{\left(2\right)}$ case provides much more room to the search for soliton-comb regimes as compared to the ${\chi }^{\left(3\right)}$ case. Conservative solitons in ${\chi }^{\left(2\right)}$ media, that can be regarded as a basis for dissipative cavity solitons, are reviewed in [19].

Axially symmetric ${\chi }^{\left(2\right)}$ microresonators are wide spread [17,18,20,21,22,23]. Their linear properties are largely common with those of the ${\chi }^{\left(3\right)}$ resonators [17,24,25,26]. Each resonator can be characterized by a major radius R, see Figure 1. The modes (whispering gallery modes) are localized near the resonator rim; the shape of the rim can be different. Each mode can be labeled by an integer azimuth number m and, additionally, two transverse modal numbers. All modal functions are proportional to $exp\left(im\phi \right)$, where $\phi$ is the azimuth angle. While exact analytical solutions exist only for the spherical case, approximate methods to calculate the modal frequencies, functions, and nonlinear interactions are well developed [24,27,28,29,30,31]. This is also relevant to the methods for coupling light into and out of resonator. Owing to the presence of a coupler (tapered fiber, prism, etc.), the inverse total (loaded) quality factor $Q^{-1}$ consists of internal and coupling contributions.

The most common lithium niobate (LN) based ${\chi }^{\left(2\right)}$ resonators are discs with the major radius $R\sim 1$ mm, the minor radius $r\lesssim R$, and quality factors $Q=({10}^7-{10}^8)$ [17,18]. The modes can be viewed as quasi-plane waves with frequencies $\omega \equiv 2\pi c/\lambda =k_{m}c/n$, discrete wavenumbers $k_m=m/R$, and an effective refractive index n depended on the modal numbers, R, and r. For $R>1$ mm, n is close to the bulk index $n_b\left(\lambda \right)$; the difference $n-n_b$ characterizes the effects of geometric dispersion. In the optical range, the azimuth number $m=2\pi nR/\lambda$ is $\sim {10}^4$. The frequency distance $\delta \omega =c/nR$ between neighboring modes with $\delta m=1$ is much larger than the line width. As LN is a uniaxial optical material, modes with ordinary (o) and extraordinary (e) polarizations have to be distinguished.

Investigations of the SH generation and parametric oscillation in ${\chi }^{\left(2\right)}$ microresonators are also known [17,18]. In LN resonators, the PM condition for the SH generation, ${\omega }_{2m}=2{\omega }_m$, can be fulfilled naturally at the FH wavelength ${\lambda }_1\simeq 1064$ nm with the use of o- and e-polarized FH and SH modes, respectively [14]. The parametric oscillation for this natural (birefringent) PM corresponds to the SH wavelength ${\lambda }_2={\lambda }_1/2\simeq 532$ nm. Here and later subscripts 1 and 2 refer to the FH and SH ranges respectively. Employment of the radial poling of LN resonators [15,16,17] enables one to fulfill the quasi-PM conditions at practically any wavelength on demand via a controllable breaking the axial symmetry. It is also necessary to mention that, because of very small resonator line widths, fine tuning means (like the temperature or geometric tuning) have to be used to control the linear and nonlinear phenomena [17,18,25,32].

Both ${\chi }^{\left(3\right)}$ and ${\chi }^{\left(2\right)}$ combs originate from parametric instabilities inside the resonators. This is why continuous-wave monochromatic pump is able to generate temporal solitons. The threshold pump power for the instability ${\mathcal{P}}_{\mathrm{th}}$ is largely controlled by the relevant modal Q factors. Roughly speaking, in the ${\chi }^{\left(3\right)}$ and ${\chi }^{\left(2\right)}$ cases ${\mathcal{P}}_{\mathrm{th}}$ is proportional to $Q^{-2}$ and $Q^{-3}$, respectively. For $Q\gtrsim {10}^7$, the ${\chi }^{\left(2\right)}$ nonlinearity is typically dominating, and ${\mathcal{P}}_{\mathrm{th}}$ belongs to the $\mu$W range [17,18].

Initial experimental attempts to get ${\chi }^{\left(2\right)}$ combs dealt with bulk PPLN systems [33,34,35,36] that are quite specific in physics. While some tens of comb lines were observed, neither solitons nor broad comb spectra were reported. Recently, first experimental efforts to generate the ${\chi }^{\left(2\right)}$ combs in microresonators have been made [37,38,39], no soliton-comb regimes are observed so far to the best of our knowledge. Furthermore, we can claim that this failure is in agreement with our theoretical results on the necessary conditions for realization of such regimes. Therefore, we see the main tasks of this mini-review in (i) determination of the conditions for the ${\chi }^{\left(2\right)}$ soliton-comb generation and ($ii$) in predictions of the main properties of the soliton-comb states. Also, we will indicate where and why the necessary conditions were violated in the known experiments. Lastly, we will compare our findings with the known early predictions for the ${\chi }^{\left(2\right)}$ soliton-comb states [40,41,42,43]. As far as possible, the cases considered will be applied to experimentally realizable and important situations.

A considerable part of the literature on nonlinear phenomena in microresonators, see, e.g., [4,17,20,39], extensively employs the quantum-mechanical terminology and notation deeply within the classical domain, where applicability of the classical Maxwell equations is beyond any doubt. Being well developed and attracting, the quantum-mechanical formalism looks here out of place and confusing. We claim, that it can easily be replaced by the classical Hamiltonian formalism [44,45]. Keeping all merits of the quantum formalism, it is fully adequate to the situation.

2. Basic Equations, Parameters, and Definitions

2.1. Periodic and Antiperiodic States

Let first the pump frequency ${\omega }_p$ be close to a modal frequency and the FH-SH phase matching be natural, i.e. the azimuth symmetry takes place. Then, two variants have to be considered, namely pumping into a FH and into a SH mode. In the FH pumping case, the PM condition $\omega \left(2m_1^0\right)\simeq 2\omega \left(m_1^0\right)$ can be satisfied with a high accuracy only for a single azimuth number $m_1^0$. It is implied that the frequency detunings ${\omega }_p-\omega \left(m_1^0\right)$ and $\omega \left(2m_1^0\right)-2\omega \left(m_1^0\right)$ are much smaller in the values than the intermodal distance $c/nR$. The light electric field $\mathit{E}(\phi ,t)$ can be represented as

$$\mathit{E}={\mathit{E}}_1{e}^{{\displaystyle \mathrm{i}(m_1^0\phi -{\omega }_{p}t)}}+{\mathit{E}}_2{e}^{{\displaystyle 2\mathrm{i}(m_1^0\phi -{\omega }_{p}t)}}+c.c.,$$

(1)

where ${\mathit{E}}_1(\phi ,t)$ and ${\mathit{E}}_2(\phi ,t)$ are slowly varying (and generally vectorial) FH and SH amplitudes (envelopes). Obviously, these amplitudes are $2\pi$-periodic in $\phi$, and the SH azimuth number $m_2^0=2m_1^0$ is even.

In the SH pumping case, the situation is different, as illustrated by Figure 2. The above PM condition can be fulfilled only if the azimuth number of the pumped mode $m_2^0$ is even (subfigure a). If $m_2^0$ is odd (subfigure b), two FH modes (with even and odd azimuth numbers) must be excited. Here, the simplest PM conditions are ${\omega }_{m_2^0}={\omega }_{m_1^0}+{\omega }_{m_1^0+1}$, $m_2^0=2m_1^0+1$. Let us use (regardless of parity of $m_2^0$) the ansatz

$$\mathit{E}={\mathit{E}}_1{e}^{{\displaystyle \mathrm{i}(m_2^0\phi -{\omega }_{p}t)/2}}+{\mathit{E}}_2{e}^{{\displaystyle \mathrm{i}(m_2^0\phi -{\omega }_{p}t)}}+c.c.$$

(2)

for the light field. Then, the first exponential factor is $2\pi$ periodic for even $m_2^0$ and $2\pi$ antiperiodic for odd $m_2^0$. As the true light field $\mathit{E}(\phi ,t)$ is $2\pi$ periodic, we must accept that the FH amplitude ${\mathit{E}}_1(\phi ,t)$ is antiperiodic for odd $m_2^0$, i.e. ${\mathit{E}}_1(\phi ,t)=-{\mathit{E}}_1(\phi \pm 2\pi ,t)$. The SH amplitude is still $2\pi$ periodic, ${\mathit{E}}_2(\phi ,t)={\mathit{E}}_2(\phi \pm 2\pi ,t)$. A dual FH-SH state with antiperiodic amplitude ${\mathit{E}}_1$ and periodic amplitude ${\mathit{E}}_2$ is named the antiperiodic state [46,47]. The antiperiodic sates (solutions) are topologically different from the conventional periodic states. No temporal evolution can make conversion between these two states.

We stress that the antiperiodic states are relevant to pumping into an odd SH mode. Otherwise, we have a periodic dual state, such that ${\mathit{E}}_1(\phi ,t)={\mathit{E}}_1(\phi \pm 2\pi ,t)$ and ${\mathit{E}}_2(\phi ,t)={\mathit{E}}_2(\phi \pm 2\pi ,t)$.

Using above the vectorial notation, we wanted to stress the generality of the periodic and antiperiodic light states for the ${\chi }^{\left(2\right)}$ nonlinearity. While the FH and SH resonator modes can be different in polarization, one can always deal with scalar FH and SH amplitudes denoted in what follows as F and S.

2.2. Equations for the Modal Amplitudes

Now, we write down nonlinear equations for the complex modal amplitudes $F_{m_1}$ and $S_{m_2}$ relevant to the FH and SH frequency domains. It is assumed, as in all known theoretical studies, that the modal polarization and the transverse modal numbers stay the same within each of these domains. We start from the phase matching case relevant to an even azimuth number $m_2^0=2m_1^0$ in the absence of radial poling (full azimuth symmetry). This takes place for the FH pumping, and is also applicable to the SH pumping into an even mode, see Figure 2a. It is useful to count the FH and SH frequencies from ${\omega }_1$ and ${\omega }_2$. In accordance with Equations (1, 2), ${\omega }_1={\omega }_p$, ${\omega }_2=2{\omega }_p$ in the FH pumping case, while ${\omega }_1={\omega }_p/2$, ${\omega }_2={\omega }_p$ in the case of SH pumping. Furthermore, we introduce instead of $m_1$ and $m_2$ more convenient FH and SH modal numbers $j=m_1-m_1^0=0,\pm 1,\dots$ and $l=m_2-m_2^0=0,\pm 1,\dots$. With these preliminaries, we have

$$\begin{array}{ccc}\hfill & \mathrm{i}& {\dot{F}}_j-Z_{1,j}{F}_j=2\mu \sum _l{S}_l{F}_{l-j}^{*}+\mathrm{i}{h}_1{\delta }_{j0}\hfill \\ & \mathrm{i}& {\dot{S}}_l-Z_{2,l}{S}_l=\mu \sum _j{F}_j{F}_{l-j}+\mathrm{i}{h}_2{\delta }_{l0}.\hfill \end{array}$$

(3)

Here the dots indicate differentiation in time, $Z_{1,j}={\omega }_j-{\omega }_1-\mathrm{i}{\gamma }_j$, $Z_{2,l}={\omega }_l-{\omega }_2-\mathrm{i}{\gamma }_l$, ${\gamma }_{j,l}={\omega }_{j,l}/2Q_{j,l}$ are modal decay rates, $\mu$ and $h_{1,2}$ are the nonlinear coupling coefficient and the internal pump amplitudes, and ${\delta }_{j0,l0}$ are Kronecker deltas. For the FH and SH pumping, we have $h_2=0$ and $h_1=0$, respectively. The modal amplitudes $F_j$ and $S_l$ are defined such that ${\omega }_j{\left|F_j\right|}^2$ and ${\omega }_l{\left|S_l\right|}^2$ are the modal energies. This definition corresponds to the classical Hamiltonian formalism [44]. The dimensionless quantities $|F_j{|}^2/\hslash \gg 1$ and $|S_l{|}^2/\hslash \gg 1$ are nothing but the number of modal FH and SH light quanta in the classical limit. The presence of a factor of 2 in set (3) reflects the Hamiltonian nature of the nonlinear interaction – the total energy is conserving for $h_{1,2}=0$ and ${\gamma }_{j,l}=0$.

The coefficients $\mu$ and $h_{1,2}$ can be determined in a quite general manner via modeling of the modes and their coupling in and out of the resonator [29,31]:

$$\mu \simeq \frac{d{\omega }_p^{3/2}}{2n^3{R}^{1/2}{\sigma }_{\mathrm{eff}}^{1/2}},h_{1,2}\simeq {\left(\frac{{\mathcal{P}}_0}{{\tilde{Q}}_{1,2}}\right)}^{1/2},$$

(4)

where d is the relevant component of the conventional susceptibility tensor, ${\sigma }_{\mathrm{eff}}$ is the effective interaction cross-section, ${\mathcal{P}}_0$ is the pump power, ${\tilde{Q}}^{-1}$ is the coupling contribution to $Q^{-1}$, and $n\simeq n_{\mathrm{o}}\simeq n_{\mathrm{e}}$ is a representative refraction index. Further details can be found in Sect. 1 (S1) of the Supplementary Materials.

The rate coefficients $Z_{1,j}$ and $Z_{2,l}$ can be expressed by the group velocities $v_{1,2}$ and dispersion coefficients $v_{1,2}^{\prime }$ if we use the expansion [48] $\omega \left(k\right)=\omega \left(k_0\right)+v\delta k+v^{\prime }\delta k^2/2$ with $\delta k=k-k_0$ near the points $k_0=m_1^0/R$ and $m_2^0/R$:

$$\begin{array}{ccc}\hfill Z_{1,j}& =& {\mathsf{\Delta }}_1-\mathrm{i}{\gamma }_1+v_{1}j/R+v_1^{\prime }{j}^2/2R^2\hfill \\ \hfill Z_{2,l}& =& {\mathsf{\Delta }}_2-\mathrm{i}{\gamma }_2+v_{2}j/R+v_2^{\prime }{l}^2/2R^2.\hfill \end{array}$$

(5)

Parameters v and $v^{\prime }$ are expressed by the refractive index $n\left(\lambda \right)$: $v=c/(n-\lambda n^{\prime })$ and $v^{\prime }=-c{\lambda }^3{n}^{\prime \prime }/2\pi {(n-\lambda n^{\prime })}^2$, where the prime indicates differentiation in $\lambda$. Knowledge of the dependences $n_{\mathrm{o},\mathrm{e}}\left(\lambda \right)$, with or without the geometric dispersion, enables one to calculate $v_{1,2}$ and $v_{1,2}^{\prime }$ as functions of ${\lambda }_2$ and ${\lambda }_1=2{\lambda }_2$. Further details can be found in S2 of supplementary materials. Expressions for detunings ${\mathsf{\Delta }}_{1,2}$ depend on the pumping scheme. For the FH pumping, ${\mathsf{\Delta }}_1={\omega }_{m_1^0}-{\omega }_p$ and ${\mathsf{\Delta }}_2={\omega }_{m_2^0}-2{\omega }_p$. In the SH pumping case, ${\mathsf{\Delta }}_2={\omega }_{m_2^0}-{\omega }_p$ and ${\mathsf{\Delta }}_1=(2{\omega }_{m_1^0}-{\omega }_p)/2$.

Consider now the remaining case of SH pumping into a mode with odd azimuth number $m_2^0$, when two FH modes with even and odd azimuth numbers enter the phase matching condition, see also Figure 2b. The SH modes can be numerated again with integer $l=m_2-m_2^0=0,\pm 1,\dots$. The only way to numerate symmetrically the FH modes is to introduce semi-integer number $j=\pm 1/2,\pm 3/2,\dots$, such that $m_1=m_2^0/2+j$. The form of sets (3) and (5) remains the same, and the only difference is in the expression for detuning ${\mathsf{\Delta }}_1$: ${\mathsf{\Delta }}_1=({\omega }_{m_1^0}+{\omega }_{m_1^0+1}-{\omega }_p)/2$.

2.3. The Impact of Radial Poling

The radial poling enables one to adjust the phase matching to practically any desirable spectral point with no effect on the linear properties. In polar ${\chi }^{\left(2\right)}$ materials, like LN, the largest components of the quadratic susceptibility tensor, $d_{13}$ and $d_{33}$, change sign under inversion of the spontaneous polarization. The periodic radial poling results in $\phi$ periodic sign changes of any of these components, as illustrated schematically by Figure 3a.

Let it be d with the absolute value $d_0$. The periodic sign-changing function $d\left(\phi \right)$ can be represented by the Fourier series

$$d=\sum _s{d}_{s}exp\left(\mathrm{i}s\phi \right).$$

(6)

If $\mathcal{N}$ is the number of poling periods, the nonzero components $d_s$ correspond to $s=0,\pm \mathcal{N},\pm 2\mathcal{N},\dots$. In the case of ± symmetric radial domain structure ($50\%$ duty cycle ratio), which is the most suitable for quasi-phase matching, only the Fourier harmonics with odd ratios $s/\mathcal{N}=\pm 1,\pm 3,\dots$ are nonzero. For these harmonics $|d_s|=2d_0\mathcal{N}/\pi \left|s\right|$, see also Figure 3b. The reduction factor $|d_s|/d_0$ decreases with increasing $\left|s\right|/\mathcal{N}$, but remains comparable with 1 for $\left|s\right|=\mathcal{N}$. Typically $\mathcal{N}\gg 1$, and the neighboring peaks are well separated from each other.

Employment of the largest first Fourier harmonics of $d\left(\phi \right)$ for quasi-phase matching corresponds to the PM condition $2{\omega }_{m_1^0}={\omega }_{2m_1^0\pm \mathcal{N}}$ for the SH generation. The sign "plus" is relevant to the most typical case of decreasing wavelength dependence of the refraction index $n\left(\lambda \right)$ [49]. We see that the SH azimuth number $m_2^0=2m_1^0\pm \mathcal{N}$ is even for even alternation number $\mathcal{N}$ and odd for odd number $\mathcal{N}$. In the second case, SH pumping leads to the excitation of periodic states. In the case of odd combination $m_2^0\mp \mathcal{N}$, the above PM condition cannot be fulfilled. However, we can fulfill the PM condition ${\omega }_{m_1^0}+{\omega }_{m_1^0+1}={\omega }_{2m_1^0\pm \mathcal{N}}$ leading to the excitation of the antiperiodic states. Thus, the quasi-PM via the radial poling results in a controllable shift of the FH azimuth number $m_1^0$. Set (3) for the FH and SH amplitudes stays unchanged.

Note that microresonators made of commercial PPLN wafers [20] are not appropriate for the spectral adjustment in question. The spectrum $|d_s|$ includes here many significant peaks $|d_{\mathcal{N}\pm 1}|,|d_{\mathcal{N}\pm 2}|,\dots$ in addition to $|d_{\mathcal{N}}|$ leading to unwanted interference of many cascading ${\chi }^{\left(2\right)}$ processes.

The possibility to adjust the phase matching to any spectral range raises the question about specific goals and properties of the radial poling. Figure 4 shows the dependences $v_{1,2}\left(\lambda \right)$ and $v_{1,2}^{\prime }\left(\lambda \right)$ for e-polarized modes in a LN based resonator with $R=1.5$ mm. The effects of geometric dispersion are very weak here. Zero walk-off point corresponds to ${\lambda }_2^c\simeq 1349$ nm, where the dispersion coefficients $v_1^{\prime }$ and $v_2^{\prime }$ are opposite in signs. Vicinity of this point, as will be shown below, is especially attractive for the soliton-comb generation. The number of periods $\mathcal{N}$, necessary to achieve quasi-PM at this point can be estimated as 284. Further details can be found in S3 of supplementary materials.

2.4. Fundamental $z,t$ Representation

The primary set of ordinary differential Equation (3) is convenient for numerical simulations. However, many fundamental features of the soliton-comb states become more evident if we transform it to a set of partial differential equations. This can be done in terms of the FH and SH amplitudes $F(\phi ,t)$ and $S(\phi ,t)$ defined by the Fourier expansions

$$F=\sum _j{F}_{j}exp\left(\mathrm{i}j\phi \right),S=\sum _l{S}_{l}exp\left(\mathrm{i}l\phi \right).$$

(7)

For integer values of $j,l$, this leads automatically to the periodic states, while for integer l and semi-integer j we arrive at a antiperiodic state.

Using Equations (3), (5), and (7) and introducing the rim coordinate $z=R\phi$, one can make sure that the periodic and antiperiodic states obey the same set of nonlinear equations for any phase matching scheme [46,47]:

$$\begin{array}{ccc}\hfill & [& \mathrm{i}({\partial }_t+v_1{\partial }_z)+0.5v_1^{\prime }{\partial }_z^2-{\mathsf{\Omega }}_1]F=2\mu S{F}^{*}+\mathrm{i}{h}_1\hfill \\ \hfill & [& \mathrm{i}({\partial }_t+v_2{\partial }_z)+0.5v_2^{\prime }{\partial }_z^2-{\mathsf{\Omega }}_2]S=\mu F^2+\mathrm{i}{h}_2,\hfill \end{array}$$

(8)

where ${\mathsf{\Omega }}_{1,2}={\mathsf{\Delta }}_{1,2}-\mathrm{i}{\gamma }_{1,2}$. At perfect resonant pumping and phase matching, we have ${\mathsf{\Delta }}_{1,2}=0$. All parameters entering set (8) are known or controllable in experiment. For the ${\chi }^{\left(2\right)}$ solitons, set (8) plays the role similar to the role of the nonlinear Lugiato–Lefever equation for the ${\chi }^{\left(3\right)}$ soliton-comb states [4,5]. However, these two cases are fundamentally different.

Consider some important properties and useful transformations of set (8):

–

This nonlinear set admits dual steady-state solutions $F(z-vt)$, $S(z-vt)$ circulating with a common velocity $v_0$ without shape changes. Solitons belong indeed to this class of solutions. Velocity $v_0$ has to be determined simultaneously with the shape of the envelopes.

–

Set (8) is written for a static coordinate frame. It is practical to rewrite it for a coordinate frame moving with velocity $v_1$. To do this, it is sufficient to drop the term $v_1{\partial }_z$ in the first equation and replace $v_2$ by $-v_{12}=v_2-v_1$ in the second one. It is assumed from now on.

–

In a dual steady state propagating with a constant velocity $v_0$, the modal amplitudes $F_j$ and $S_l$ oscillate in time as $exp(-\mathrm{i}j{v}_{0}t/R)$ and $exp(-\mathrm{i}l{v}_{0}t/R)$. This elementary property is useful to control the establishment of dual steady states in numerical simulations, see Section 1.

–

Attempts to rewrite set (8) for the antiperiodic states using the replacements $F(z,t)=\tilde{F}(z,t)exp(\pm \mathrm{i}z/2R)$ lead to restoration of the periodicity but also to an explicit z-dependence of the right-hand sides (non-autonomous system). This in inappropriate.

The structure of Equation (8) gives a hint about normalization of parameters relevant to the linear properties of the system. Such a normalization enables one to decrease the number of variable parameters. In what follows, we will use the following normalized detunings, walk-off, and dispersion parameters ${\delta }_{1,2}$, $\alpha$, and ${\beta }_{1,2}$:

$${\delta }_{1,2}=\frac{{\mathsf{\Delta }}_{1,2}}{{\gamma }_{1,2}},\alpha =\frac{v_{12}}{{\gamma }_{2}R},{\beta }_{1,2}=\frac{v_{1,2}^{\prime }}{2{\gamma }_{1,2}{R}^2}.$$

(9)

As a rule, the actual values of detunings are restricted by the inequalities $|{\delta }_{1,2}|\lesssim 1$. The actual values of $\alpha$ and ${\beta }_{1,2}$ substantially depend on the group velocity difference $v_{12}=v_1-v_2$, the dispersion coefficients $v_{1,2}^{\prime }$, and the decay rates ${\gamma }_{1,2}$.

For LN resonators, typical values of the decay rates are $\sim {10}^7$, such that $|{\beta }_{1,2}|\ll 1$. On contrary, we have typically $\left|\alpha \right|\gg 1$. Exception from this rule is possible in a close vicinity of zero walk-off point ${\lambda }_2\simeq 1349$ nm. The decay rates ${\gamma }_1$ and ${\gamma }_2$ are typically not strongly different, such that assumption ${\gamma }_{1,2}=\gamma$ can be used sometimes to simplify modeling.

Generally speaking, the values of F, S, and $h_{1,2}$ also need some normalization. It is, however, different in the FH and SH pumping cases. We will provide the relevant information when necessary.

One more important property of set (8) is the presence of spatially uniform steady-state background solutions $\overline{F}\ne 0$, $\overline{S}\ne 0$ [50]. In steady state, $F_0=\overline{F}$ and $S_0=\overline{S}$. Obviously, such a dual background is possible only for the periodic states. The dual background is closely related both to the properties of solitons and to the thresholds of instabilities leading to soliton-comb formation. Omitting all derivatives in (8), one can easily come to a closed algebraic equation for $|\overline{F}{|}^2$ that is, however, different for the FH and SH pumping cases. We have $|\overline{S}|=\mu |\overline{F}{|}^2/\left|{\mathsf{\Omega }}_2\right|$ and $|{\mathsf{\Omega }}_1|/2\mu$ for the FH and SH pumping cases, respectively. Dependence of $|\overline{F}{|}^2$ on the pump amplitudes $h_{1,2}$ is presented graphically in Figure 5.

In both pumping cases, dependences $|\overline{F}{|}^2\left(h_{1,2}\right)$ are controlled by the phase $\mathsf{\Phi }=arg\left({\mathsf{\Omega }}_1{\mathsf{\Omega }}_2\right)$, such that multi-valued solutions exist for sufficiently small $\left|\mathsf{\Phi }\right|$.

Importantly, the dual background with $\overline{F},\overline{S}\ne 0$ always exists in the FH pumping case, see Figure 5a. At the SH pumping, the situation is different. For sufficiently small values of $h_2$, there is a single background state, such that $\overline{F}=0$ and $|\overline{S}|=h_2/\left|{\mathsf{\Omega }}_2\right|\ne 0$, see Figure 5b. This difference has dramatic consequences for the dual comb-soliton states, see Section 4,Section 5,Section 6. Note also that Figure 5 gives a clear hint about normalization of the amplitudes $F,S$ and pump amplitudes $h_{1,2}$ in the FH and SH pumping cases.

Lastly, we mention two limiting cases for set (8).

–

In the limit $h_{1,2}\to 0$, ${\gamma }_{1,2}\to 0$ (no gain, no modal decay) we proceed to the conservative case reviewed in [19]. A number of exact soliton solutions are known here, but they are far from the subject of our study.

–

In the the dissipativeless driven limit ${\gamma }_{1,2}\to 0$ at ${\mathsf{\Delta }}_{1,2}\ne 0$ and $h_{1,2}\ne 0$, one can also obtain exact soliton solutions [51]. It turns out, however, that they all are unstable, i.e. cannot be realized.

3. Instability Thresholds

As known, generation of second harmonic is thresholdless – it corresponds just to forced oscillation caused by the ${\chi }^{\left(2\right)}$ nonlinearity. This thresholdless process becomes merely more efficient in the presence of FH-SH phase matching. On contrary, excitation of FH harmonics in the presence of SH modes is always a coherent threshold process. This explains a big difference in the properties of parametric excitation for the FH and SH pumping schemes in microresonators [50].

The simplest and most efficient is the SH pumping scheme. At sufficiently small pump powers (small values of $h_2$ at $h_1=0$), we deal simply with the linear excitation of the pumped mode $m_2^0$. Its steady-state amplitude is $S_0=\overline{S}=h_2/({\gamma }_2+i{\mathsf{\Delta }}_2)$; detuning ${\mathsf{\Delta }}_2$ characterizes just the distance to the linear resonance. To see the effect of $S_0$ on the FH harmonics, we consider the first of Equation (3) with $S_l=S_0{\delta }_{l,0}$ and $h_1=0$. The quantities $F_j$ and $F_{-j}^{*}$ with an arbitrary j obey obviously two linear differential equations. Can they have a nonzero steady-state solution? Using the zero determinant solvability condition for linear algebraic equations and Equation (9) for the normalized parameters, we obtain the threshold condition in the form [46,47,50]

$${\eta }_{2,j}^{\mathrm{th}}={\left[1+{({\delta }_1+{\beta }_1{j}^2)}^2\right]}^{1/2}{(1+{\delta }_2^2)}^{1/2},$$

(10)

where ${\eta }_2=2\mu h_2/{\gamma }_1{\gamma }_2$ is the normalized pump amplitude. The same threshold equation can be obtained using a more general approach: Namely, one can set $F_j,F_{-j}^{*}\propto exp(-\mathrm{i}\nu t)$ to determine the complex characteristic exponent $\nu ={\nu }^{\prime }+\mathrm{i}{\nu }^{\prime \prime }$ from an algebraic solvability condition. Then, one can make sure that ${\nu }^{\prime \prime }>0$ for ${\eta }_2>{\eta }_{2,j}^{\mathrm{th}}$, while ${\nu }^{\prime \prime }<0$ for ${\eta }_2<{\eta }_{2,j}^{\mathrm{th}}$. Thus, small amplitudes $F_{\pm j}\left(t\right)$ decay below the threshold and grow exponentially above it. This fully justifies the use of the threshold EquatiEon (10).

Consider in some detail the threshold properties given by Equation (10). Remarkably, it includes neither the walk-off parameter $\alpha$ nor the SH dispersion coefficient ${\beta }_2$. The threshold equation refers to the FH pairs j and $-j$. The value of j can be integer or semi-integer depending on the type of solution in question. The thresholds are generally different for different $\left|j\right|$. However, the minimal in ${\delta }_1$ values of ${\eta }_2^{\mathrm{th}}$ are the same for all $\left|j\right|$. The value of ${\delta }_1$ minimizing the threshold for a pair $\pm j$ is ${\delta }_1^{\left(j\right)}=-{\beta }_1{j}^2$; it is substantially different for different $j^2$. Using this property, one can excite selectively any pair $j,-j$ on demand. The absolute threshold minimum, ${\left[{\eta }_2^{\mathrm{th}}\right]}_{min}=1$ can be achieved indeed at ${\delta }_2=0$, i.e. at the exactly resonant excitation of the pumped SH mode. It is clear now that parameter ${\eta }_2$ can be expressed by the pump power as ${\eta }_2={(\mathcal{P}/{\mathcal{P}}_{min}^{\mathrm{th}})}^{1/2}$.

Now we turn to the FH pumping scheme. It is relevant not only to the comb generation issue, but also to the SH excitation. Complexity of this scheme, compared to the previous one, has fairly simple reasons. The parametric instability requires sufficiently large values of the SH amplitude $S_0=\overline{S}$. They can be achieved, in turn, only for sufficiently large values of the amplitude of the pumped mode $F_0=\overline{F}\left(h_1\right)$. Keeping this fact in mind, we consider what happens if we admit the presence of small perturbations $\delta F(\phi ,t)$ and $\delta S(\phi ,t)$, such that $F=\overline{F}+\delta F$ and $S=\overline{S}+\delta S$. Substituting these expressions in Equation (8) with $h_2=0$, we see that the SH perturbation $\delta S$ becomes linearly linked to the FH perturbation via the term $\overline{F}\delta F$, which was absent earlier. As the result, we must deal with a quartet of perturbations $F_j$, $F_{-j}^{*}$, $S_j$, $S_{-j}^{*}\propto exp(-\mathrm{i}\nu t)$ instead of the former parametric pair $F_j,F_{-j}^{*}$. The characteristic equation for $\nu$ becomes a cumbersome four-degree algebraic equation including all linear parameters ${\delta }_{1,2}$, ${\beta }_{1,2}$, and $\alpha$ [50]. The threshold properties become complicated, they cannot be covered by simple relations similar to Equation (10).

Nevertheless, a general overview and illustrations of the threshold properties, which do not pretend to be complete, are possible. We introduce first the normalized pump amplitude ${\eta }_1=\sqrt{2}\mu h_1/{\gamma }_1\sqrt{{\gamma }_1{\gamma }_2}$, this normalization is different from that relevant to the SH pumping. The presence of parameter $\alpha$, ranging from 0 to very large values, in the threshold equation, raises question about the effect of walk-off. It turns out that the lowest threshold value of ${\eta }_1$ at zero walk-off ($\alpha =0$) corresponds to $j=0$, i.e. to an internal instability of the dual background [50]. For detunings ${\delta }_{1,2}=0$, this value is minimal, ${\eta }_{1,0}^{\mathrm{th}}=3\sqrt{2}$. With increase of $\left|\alpha \right|$, some quartets with $j\ne 0$ become quickly the easiest to excite. However, the optimal values of $\left|j\right|$ depend significantly on the dispersion coefficients. This is illustrated by Figure 6 obtained with the use of numerical methods.

SubFigure 6a is plotted for ${\delta }_{1,2}=0$ and representative dispersion parameters relevant to the quasi-PM in LN. Within a broad range of $\left|\alpha \right|$, the lowest threshold corresponds to harmonics with very large values of $\left|j\right|$. SubFigure 6b corresponds to the natural PM, where the dispersion coefficients ${\beta }_{1,2}$ are negative and relatively large in the absolute values. The effect of $\left|\alpha \right|$ is strongly different here and it strongly depends on $\left|j\right|$. The reason is in competition of the walk-off and dispersion terms in the rate coefficient $Z_{2,l}$ entering Equation (3) for S. Anyhow, selective excitation of modes with $\left|j\right|=1,2$ is practically impossible. The effects of nonzero detunings ${\delta }_{1,2}$ can also be investigated numerically. They also do not facilitate selective excitation of harmonics with small values of $\left|j\right|$.

With increasing $\left|\alpha \right|$, the threshold values ${\eta }_{1,j}^{\mathrm{th}}\left(\alpha \right)$ slowly approach 2 from above. This range is, however, not favorable for the soliton-comb generation because of a strong suppression of the side SH harmonics with $l\ne 0$, see also Section 5 and Section 6.

4. Numerical Methods

We simulated predominantly set (3) in the coordinate frame moving with velocity $v_1$ using the fourth-order Runge-Kutta method. Initially, for verification purposes, we solved also set (8) of partial differential equations using the step-split Fourier method. The results were practically equivalent. The number of Fourier harmonics taken into account (and the mesh points in $\phi$) ranged from 32 to 1024 within each of FH and SH frequency domains. Accuracy of the calculations was controlled by changing the time step and the number of harmonics. Different values of the variable input dimensionless parameters of the system (${\eta }_{1,2},{\delta }_{1,2},{\beta }_{1,2},\alpha$) were used. With harmonics $F_j\left(t\right),S_l\left(t\right)$ calculated, one can investigate temporal evolution of the spatial profiles $F\left(\phi \right),S\left(\phi \right)$ and comb spectra $|F_j{|}^2,{\left|S_l\right|}^2$. Also, it is possible to verify establishment of steady states moving with a constant velocity and, if applicable, determine this velocity $v_0$.

Three main problems were considered numerically:

– We verified numerically the analytical results relevant to the instability thresholds. To do so, we used initial conditions with very small random complex values of $F_j\left(0\right)$ and $S_l\left(0\right)$ to see the initial exponential growth or decay of $|F_j|\left(t\right)$ for different values of ${\eta }_{1,2}$. Coincidence of the analytical and numerical results with a high accuracy was always achieved.

– We verified stability of the dual soliton solutions obtained analytically in [46,51] for the known limiting cases. As initial conditions, we used Fourier transforms of the relevant analytical expressions for $F\left(\phi \right)$ and $S\left(\phi \right)$. Numerical analysis has shown that all the analytical solutions are unstable: Temporal evolution leads inevitably to irregular spatial profiles and comb spectra.

– We tried to generate stable steady-state soliton-comb states starting from very small random complex amplitudes $F_j\left(0\right)$ and $S_l\left(0\right)$. As a rule, evolution of our nonlinear systems well above the threshold leads ultimately, after a stage of exponential growth, to highly irregular saturated behavior. For an abrupt switching the pump on, rough features of this irregular behavior depend also on the choice of initial conditions. Nevertheless, for many well defined choices of the input parameters, proper pumping scheme, and not very abrupt switching of the pump on, temporal evolution leads reliably and uniquely to multi-parametric families of comb-soliton states. Below we describe the corresponding adiabatic procedure of growing of soliton-comb states as applied (for definiteness) to the SH pumping case.

Let FH harmonics $\pm j$ possess the lowest threshold. For ${\eta }_2<{\eta }_{2,j}^{\mathrm{th}}$, we have then the steady state with $S_0\ne 0$ and zero other harmonics. Switching to ${\eta }_2$ slightly above ${\eta }_{2,j}^{\mathrm{th}}$, we obtain a steady state including $S_0$ and only few relatively small other essential FH and SH harmonics. After that, the steady-state harmonics are used as new initial conditions at a slightly larger value of ${\eta }_2$. This procedure is repeated then many times to come to a steady state well above the threshold. The number of FH and SH harmonics taken into account has to be chosen properly. The final state is insensitive to the choice of initial random amplitudes. The states achieved in this way can be qualified as self-starting steady states; their stability is doubtless. Let us stress that achievement of the soliton-comb states with the adiabatic procedure is not guaranteed. However, in many cases it takes place.

Restriction on the rise time of the pump, ${\tau }_p$, necessary to achieve the steady state is found to be not very obligating. It is sufficient to fulfill the inequality ${\gamma }_{1,2}{\tau }_p\gg 1$, which corresponds typically to ${\tau }_p\gtrsim 1\mu$s. Note also, that soliton-comb states can be achieved often in a few steps of step-wise increase of the pump power.

Establishment of soliton-comb steady states and the corresponding velocity $v_0$ can be quantified with a high accuracy using a correlation method. We use the fact that Fourier harmonic $A_j\left(t\right)$ ($F_j$ or $S_j$) known in the coordinate frame moving with velocity $v_1$ can be recalculated in the frame moving with an arbitrary velocity $v_{*}$ through multiplication by $exp[\mathrm{i}j(v_1-v_{*})t/R]$. Let us introduce the dimensionless error parameter

$$\epsilon (t,\tau )=\frac{{\sum }_j{|A_j\left(t\right)-A_j(t-\tau )|}^2}{{\sum }_j\left(\right|A_j{\left(t\right)|}^2+|A_j(t-\tau ){|}^2)},$$

(11)

where $A_j$ refers to a coordinate frame moving with velocity $v_{*}$, t is the calculation time and $\tau$ is a variable time shift. Obviously, $\epsilon (t,\tau )$ turns to zero only when we deal with the steady state and, simultaneously, $v_{*}=v$. The error parameter calculated for modestly large evolution times, ${\gamma }_{1,2}t\gtrsim {10}^3$, and minimized over $v_{*}$ shows extremely small values ($\epsilon ={10}^{-14}-{10}^{-15}$) caused by the numerical noise, see also below. For smaller t, i.e., during the transient stage, it is larger by many orders of magnitude. We have thus a numerical tool to control the establishment of steady states and determine precisely velocity $v_0$.

Figure 7 shows three examples of temporal evolution of $\epsilon (t,\tau )$ relevant to our adiabatic procedure applied to the antiperiodic case at ${\gamma }_{1,2}=\gamma ={10}^7$ s⁻¹, $\tau =1/\gamma$, and parameters ${\beta }_1=0.02$, ${\beta }_2=-0.01$, $\alpha =1$ representative for a close vicinity of the zero walk-off point in LN based resonators. The harmonics in question are $F_j$. The normalized pump amplitude ${\eta }_2$ takes sequentially 17 values above the threshold: $1.01$, $1.1$, …, and 30. SubFigure 7a,b,c exhibit evolution of $\epsilon$ after taking the first two and the last values at $t_0$, $t_1$, and $t_{17}$. The most general feature this evolution is always the same: It is an almost exponential decrease of $\epsilon$ from initial values of the order of 1 to extremely low values ${10}^{-14}-{10}^{-15}$. Note also some quantitative details. The longest evolution, subFigure 7a, starts from modal noise. It takes about $(5-7){\gamma }^{-1}$ for the harmonics to grow exponentially and reach the saturation. The corresponding evolution time can be identified with the rise time of near-threshold optical parametric oscillation. The further evolution runs, which start from regular distributions of the FH and SH harmonics achieved during the previous runs, occur substantially faster. Our choice of 17 steps of evolution to reach the state with $\eta =30$ is largely arbitrary. In many cases, the same final state can be achieved in a few steps.

In the next sections we will see in details the spatial and spectral properties of the dual states grown with the use of our adiabatic procedure.

5. Soliton-Comb States for SH Pumping

Below we report the results on soliton-comb generation obtained with our adiabatic procedure in the SH pumping case for LN resonators [47,52]. The $\phi$-dependences of ${\left|F\right|}^2$, ${\left|S\right|}^2$, $arg\left[F\right]$, and $arg\left[S\right]$ and the spectra $|F_j{|}^2$ and $|S_l{|}^2$ are of our prime interest. Also, we are interested in dependences of the soliton-comb characteristics on the pump amplitude ${\eta }_2$, the walk-off parameter $\alpha$, and other variable parameters of the system. Generally, the soliton position in the resonator depends on the choice of initial random amplitudes. Since the azimuth angle $\phi$ can be counted from any point, we can chose its zero point as convenient.

We start our report with the case of quasi-PM via the radial poling. Here, the pump wavelength ${\lambda }_p$ can take practically any value. In particular, it can be made arbitrary close to the zero walk-off point ${\lambda }_2^c\simeq 1349$ nm. This case is especially important. For sufficiently small values of $\delta \lambda ={\lambda }_p-{\lambda }_2^c$, the group velocity difference is given by $v_{12}\left[{10}^5\mathrm{cm}/\mathrm{s}\right]\simeq -5.9\times \delta \lambda \left[\mathrm{nm}\right]$ in accordance with Figure 4a. Correspondingly, for $\left|\delta \lambda \right|\lesssim 10$ nm and typical values of R and ${\gamma }_2$ the walk-off parameter $\alpha =v_{12}/{\gamma }_{2}R$ ranges from 0 to large negative and positive values, while the dispersion coefficients ${\beta }_{1,2}\left(\delta \lambda \right)$ stay practically constant and opposite in signs, ${\beta }_1>0$, ${\beta }_2<0$.

5.1. Antiperiodic solutions near ${\lambda }_2^c\simeq 1349$ nm

Figure 8 corresponds to an antiperiodic steady state moving with velocity $v_{01}=v_0-v_1=0$. It is relevant to zero walk-off ($\delta \lambda =0$, $\alpha =0$), ${\eta }_2=30$, zero detunings (${\delta }_{1,2}=0$), $R=1.5$ mm, and decay rates ${\gamma }_1={10}^7$, ${\gamma }_2={10}^8$ s⁻¹ ($Q_1\simeq 3.3\times {10}^7$, $Q_2\simeq 0.7\times {10}^7$). The total number of harmonics taken into account is 1024. SubFigure 8a,b show the central parts of the normalized FH and SH intensity profiles – far from the center the intensities stay constants. It is evident that we deal with a narrow dual dark-bright soliton. Far from the soliton core, the normalized values of ${\left|F\left(\phi \right)\right|}^2$ and ${\left|S\left(\phi \right)\right|}^2$ ($\simeq 29$ and 1) nicely correspond to the dual background considered in Section 2.

At the first sight, this feature is not consistent with the absence of the dual background for the antiperiodic states. The point is that the amplitude $F\left(\phi \right)$ tends to the opposite left and right constant values far from the center, and these values can be found from Equation (8). The facts that (i) the intensity profiles are even in $\phi$ and (ii) the FH intensity ${\left|F\left(\phi \right)\right|}^2$ turns exactly to zero at the soliton center are exclusively due to our zero walk-off assumption.

SubFigure 8c,d exhibit $\phi$-dependences of the phases $arg\left[F\right]$ and $arg\left[S\right]$. While an overall $\pi$ change of $arg\left[F\right]$ across the resonator is general for the antiperiodic states, the $\pi$ jump at $\phi =0$ is typical only for $\alpha =0$. More specifically, the line $\mathrm{Re}\left[F\right]\left(\phi \right),\mathrm{Im}\left[F\right]\left(\phi \right)$ on the complex plane crosses the zero point at a finite angle because of the symmetry leading to an abrupt $\pi$ jump of the phase. For $\alpha \ne 0$, the function $arg\left[F\right]\left(\phi \right)$ is smooth. As concerned the SH phase $arg\left[S\right]$, it shows only smooth (but large) deviations from zero near the soliton center, see Figure 8d.

The FH and SH comb spectra are shown in subFigure 8e,f. They are symmetric in j and l and very broad. The normalizations of F and S are slightly different and specific for the SH pumping case. On the other hand, they provide the maximum values of $|F_j|$ and $|S_l|$ of the order of 1. The shapes of the FH and SH spectra are notably different. The central peak in Figure 8f) is due to the pump. To quantify the comb spectra, we define the number of significant FH and SH lines ($N_1$ and $N_2$) whose normalized intensities are above the cut-off level ${10}^{-4}$. This characteristic is conditional and rather useful; the lines of this strength are typically well above the light noise. In our case $N_1=122$ and $N_2=161$ corresponding to well developed dual combs.

With increasing ${\eta }_2$, the dual soliton is progressively narrowing and the comb spectra are getting broader. No deterioration of the soliton-comb states was seen up to the maximal investigated values ${\eta }_2^{max}\approx {10}^2$. Despite an apparent simplicity of the soliton profiles at $\alpha =0$, no analytical solutions for them is known so far.

It is not difficult to estimate the numbers of light quanta $|F_j{|}^2/\hslash$ and $|S_l{|}^2/\hslash$ in the FH and SH modes for the data of Figure 8e,f. Using Equation (4) for $\mu$ and setting ${\sigma }_{\mathrm{eff}}=100\mu$m², we obtain that they are roughly of the order of ${10}^6$ for $\left|j\right|\sim \left|l\right|\sim 1$. Thus, we are deeply within the classical range. This is valid also for the subsequent illustrations. At the same time, far tails of the comb spectra can experience quantum effects.

Next, we consider the effects of nonzero walk-off. Nonzero values of $\alpha$ influence neither the thresholds nor the outcome of the adiabatic procedure – the system still evolves to a single-soliton steady state. The most evident effects are (i) nonzero soliton velocity, $v_{01}\ne 0$, and (ii) asymmetry of the soliton profiles and the comb spectra. The larger $\left|\alpha \right|$, the stronger is this asymmetry. This is illustrated by Figure 9. Compared to Figure 8 we have changed $\delta \lambda$ from 0 to $-1.08$ nm, this corresponds to $\alpha \simeq 0.05$. Asymmetry of the soliton profiles and comb spectra is already evident. Furthermore, the soliton profiles acquire pronounced tails. The value of ${\left|F\left(\phi \right)\right|}_{min}^2$ is now far from zero, and ${\left|S\left(\phi \right)\right|}_{max}^2$ is about two times smaller than earlier. At the same time, the background values of the FH and SH intensities are the same, and the total number of the comb lines $N=N_1+N_2$ is even slightly larger than it was at $\alpha =0$. Changing sign of $\delta \lambda$ results in the mirror reflection of the graphs of Figure 9 about the central vertical line and in changing sign of the soliton velocity $v_{01}$.

To gain new insights, we consider color maps of velocity parameter $|v_{01}|/2\pi R$ and of the total comb line number $N=N_1+N_2$ on the $\delta \lambda ,{\eta }_2$ plane for a few representative combinations of the decay rates ${\gamma }_1$ and ${\gamma }_2$, see Figure 10. Three first subfigures of the upper row represent the maps of $|v_{01}|/2\pi R$ for the combinations ${\gamma }_{1,2}=\gamma$ (a1), ${\gamma }_1=\gamma$, ${\gamma }_2=10\gamma$ (b1), and ${\gamma }_{1,2}=10\gamma$ (c1) with $\gamma ={10}^7$ s⁻¹. Subfigures a2, b2, and c2 of the lower row represent the corresponding maps of the total comb line number N. Subfigures d1 and d2 give additionally cross-sections $v_{01}\left(\delta \lambda \right)$ and $N\left(\delta \lambda \right)$ at ${\eta }_2=20$ for our combinations of ${\gamma }_1$ and ${\gamma }_2$.

As evident from the maps, both $|v_{01}|$ and N grow monotonously with increasing ${\eta }_2$. As concerned dependences of these parameters on $\left|\delta \lambda \right|$ (or $\left|\alpha \right|$), they are decreasing only outside a vertical central strip, whose width depends on ${\gamma }_{1,2}$. The decrease of $N\left(\right|\delta \lambda \left|\right)$ and $|v_{01}\left|\right(\left|\delta \lambda \right|)$ persists for $\left|\delta \lambda \right|>10$ nm. In short, large values of the walk-off parameter, $\left|\alpha \right|\gg 1$ are harmful for the soliton-comb generation in question.

The range of small $\left|\delta \lambda \right|$ (or $\alpha$), illustrated in some details in the subfigures d1 and d2, is worthy of attention. The dependences $v_{01}\left(\delta \lambda \right)$ and $N\left(\delta \lambda \right)$ at large values of ${\eta }_2$ are not always continuous. This means that small variations of $\delta \lambda$ can lead to notably different results when adiabatically increasing ${\eta }_2$. Continuous changes of $v_{01}\left(\delta \lambda \right)$ and $N\left(\delta \lambda \right)$ at ${\eta }_2=20$ in Figure 10 take place only for ${\gamma }_1=0.1{\gamma }_2=\gamma$. The optimum value of $\left|\delta \lambda \right|$ and the corresponding largest peak values of N take place for ${\gamma }_{1,2}=10\gamma$. In accordance with our definition of ${\eta }_2$ they correspond to the largest $h_2$ (and the pump power). Anyhow, the increase of N via optimization of $\delta \lambda$ is not very large. Note lastly that the wavelength range of 20 nm relevant to Figure 10 strongly exceeds the intermodal wavelength distance of $\simeq 0.08$ nm, so that the indicated spectral features are not available for the standard experimental fine tuning means.

Now we turn to the effects of the frequency detunings. These effects are of two kinds. One can consider the effects of ${\delta }_{1,2}$ on initiation of soliton-comb states via an adiabatic increase of ${\eta }_2$. These strong and important effects are considered in Section 5.2. Alternatively, we can vary slowly ${\delta }_{1,2}$ upon reaching the above considered soliton states at large values of ${\eta }_2$. This issue, which is largely (but not only) about the area of existence of developed soliton-comb states, is considered below.

We have extended our adiabatic procedure for slow independent changes of ${\delta }_1$ and ${\delta }_2$. In fact, two similar scanning procedures were used. Within the first one, we determined initially the soliton velocity $v_{01}$ and the comb line number N by increasing and decreasing ${\delta }_1$ at ${\delta }_2=0$. After that, starting from a point ${\delta }_1,0$, we increased and decreased detuning ${\delta }_2$. The achieved steady-state values of the amplitudes $F_j$ and $S_l$ were used as new initial conditions in each step of changing ${\delta }_1$ or ${\delta }_2$. Within the second procedure, the scanning order was inverted: initially we varied ${\delta }_2$ at ${\delta }_1=0$, and after that increased and decreased ${\delta }_1$ starting from points $0,{\delta }_2$. We have verified whether the functions $v_{01}({\delta }_1,{\delta }_2)$ and $N({\delta }_1,{\delta }_2)$ depend on the scanning order. Also, we have made sure that each steady state on the ${\delta }_1,{\delta }_2$ plane is still relevant to a single-soliton antiperiodic state.

Consider first the case zero walk-off, $\alpha =0$. In our coordinate frame moving with velocity $v_1$, the positive and negative propagation directions are equivalent, such that the soliton velocity difference $v_{01}=v_0-v_1$ is expected to be zero. Surprizingly, we have found that this is not always the case. Figure 11a represents color map of $|v_{01}|/2\pi R$ for $\delta \lambda =0$, ${\gamma }_1={10}^7$, and ${\gamma }_2={10}^8$ s⁻¹ (as in Figure 8). The uniformly dark-blue colored lower and upper parts, which are separated from the rest by two slanted borders (bifurcation lines), correspond to $v_{01}/2\pi R=0$ within an accuracy of $\sim {10}^{-9}$ MHz. Between these borders we have a smooth distribution $|v_{01}|({\delta }_1,{\delta }_2)$ whose scale is comparable with that of Figure 10. This distribution practically does not depend on the scanning order. But what is the sign of $v_{01}$? We have found that it randomly depends on fine particularities of our calculation procedure, such as the presence of very weak remnant

velocity perturbations. In essence, we have a clear example of the spontaneous symmetry breaking, when the state with a high symmetry ($v_{01}=0$) becomes unstable against the excitation of one of two equivalent states of lower symmetry with $|v_{01}|\ne 0$. One of the most known examples of such a symmetry breaking is the ferroelectric second-order transition below the Curie temperature.

Let us now take a look at the corresponding color map of N, Figure 11b. The distribution $N({\delta }_1,{\delta }_2)$ is rather smooth. The lower bifurcation line of Figure 11a) is barely seen here, while the upper line is pronounced. Remarkably, increase of ${\delta }_2$ at ${\delta }_1\approx 0$ is favorable for the comb in spite of the increasing instability threshold, see Equation (10). The characteristic scale in ${\delta }_1$ is here about 1, which strongly exceeds the scale ${\beta }_1=0.02$ relevant to the changes in ${\eta }_2^{\mathrm{th}}\left({\delta }_1\right)$. This assertion is general for ${\eta }_2\gg 1$; it is justified by plotting maps of $N({\delta }_1,{\delta }_2)$ for different values of $\delta \lambda$ and different combinations of ${\gamma }_1$ and ${\gamma }_2$. In short, we have a vast family of single-soliton antiperiodic states continuously depending on the input parameters.

Our attempts to generate periodic soliton-comb solutions at ${\delta }_{1,2}=0$ failed. The reason and the possibilities to get new solutions are indicated in Section 5.2.

5.2. Selective Generation of Periodic and Antiperiodic Milti-Soliton States

The idea for generation of periodic and antiperiodic multi-soliton states roots in the threshold conditions of Section 3 and a simple observation. The point is that the FH profile of Figure 8 with a sharp minimum emerges with increasing ${\eta }_2$ as a continuous (without bifurcations) development of the near-threshold profile ${\left|F\right|}^2\propto {sin}^2(\phi /2)$ corresponding to the lowest threshold for FH harmonics with $j=\pm 1/2$. One can suggest thus that the spatial structure of developed solitons can be determined by the threshold conditions.

To develop this idea, we consider the excitation diagrams of Figure 12.

Diagrams a) and b) correspond to the antiperiodic and periodic cases. The cyan circles on the horizontal ${\delta }_1/{\beta }_1$ axis separate the regions where the lowest thresholds correspond to different values of the FH harmonic number $\left|j\right|$. Within a broad range in diagram a), ${\delta }_1>-{\beta }_1/2$, including the zero point, an adiabatic increase of ${\eta }_2$ is expected to lead to single-soliton antiperiodic states. This is in line with the above considered results. For $-9{\beta }_1/4<{\delta }_1<-{\beta }_1/4$, we can expect three-soliton antiperiodic states, etc.

The periodic case, diagram b), is expected to be different. In the range ${\delta }_1>-{\beta }_1/2$, including the zero point, the lowest threshold corresponds to $j=0$, i.e. to the generation of the dual background. As this background is stable, no periodic solitons are expected. However, in the range $-5{\beta }_1/2<{\delta }_1<-{\beta }_1/2$, where the lowest threshold corresponds to $\left|j\right|=1$, two-soliton periodic state is expected. Similarly, we can hope to generate four-soliton periodic states for ${\delta }_1<-5{\beta }_1/2$. As concerned detuning ${\delta }_2$, it plays a passive role by increasing non-selectively all thresholds according to Equation (10).

Now we consider the results of the corresponding numerical experiments on selective growing of multi-soliton states [53]. Figure 13 shows the periodic two-soliton state grown via adiabatic increase of ${\eta }_2$ for $\delta \lambda =0$, ${\gamma }_{1,2}=\gamma ={10}^7$ s⁻¹, ${\delta }_2=0$, and ${\delta }_1=-{\beta }_1$. Remarkably, the phase $arg\left[F\right]\left(\phi \right)$ shows two opposite $\pi$-steps when crossing each of two narrow soliton areas. Thus, our periodic two-soliton solution consists of two equally spaced mutually repelling antiperiodic solitons. This structure emerges every time spontaneously without imposing the antiperiodic conditions. Velocity of the state is $v_1=v_2$. The FH and SH comb spectra, shown in Figure 13c,f, are symmetric and structurally similar to those of Figure 8e,f. Because of $\pi$-periodicity in $\phi$, the comb line spacing is $2v_1/R$, it is doubled compared to the single-soliton case. Correspondingly, the number of significant comb lines (above ${10}^{-4}$) is about two times smaller than it would be in the single-soliton antiperiodic case at ${\eta }_2=20$.

Changing ${\delta }_1$ sequentially to ${\delta }_1^{\left(j\right)}=-4{\beta }_1$ and $-9{\beta }_1$, which corresponds to the lowest threshold for $\left|j\right|=2$ and 3, we have grown periodic four- and six-soliton states. They all consist of pairs of the same antiperiodic solitons. The comb line spacing increases with $\left|j\right|$ accordingly. Further increase of $\left|j\right|$ does not lead, however, to new multi-soliton states. Temporal development tends to break such states into the former two-, four-, and six-soliton solutions. This feature is not related to insufficiently small time step and/or insufficiently large number of harmonics taken into account. Also, we mention the case ${\delta }_1=0$ relevant to the lowest threshold for $j=0$. Our adiabatic procedure leads here, as expected, to a dual background steady state $\overline{F},\overline{S}$.

Generalization to the case of nonzero walk-off ($\alpha \ne 0$) and different decay rates ${\gamma }_1$ and ${\gamma }_2$ is rather trivial. We still have multi-soliton periodic solutions consisting of pairs of the former antiperiodic solitons moving with velocity $v_0\ne v_{1,2}$.

Now we switch to multi-soliton antiperiodic states relevant to semi-integer values of j, see Figure 12a. Figure 14 shows representatively the FH and SH intensity profiles for the antiperiodic three-soliton state grown for ${\delta }_1^{(3/2)}=-9{\beta }_1/4$ and $\alpha =-1$. The other input parameters are the same as in Figure 13.

Each of three equally spaced solitons has a familiar structure with oscillating tails. Owing to a relatively large negative value of $\alpha$, the tails are strongly pronounced and the spatial asymmetry is inverted as compared to Figure 9. The soliton velocity difference $v_{01}\simeq 2.2\times {10}^5$ cm/s. Interestingly, antiperiodic five- and seven-soliton states can be easily grown for $\left|j\right|=5/2$ and $7/2$ in spite of strong tail overlaps. This contrasts with multi-soliton periodic states.

5.3. Natural PM: The Absence of Solitons

The natural phase matching, that can be realized in LN based resonators at ${\lambda }_2\simeq 532$ nm, is of special interest. Avoiding complicated procedure of the radial poling and availability of tunable laser sources in this range make it very attractive. Furthermore, it was indicated [54] that the value of the walk-off parameter can be substantially decreased for sufficiently small microresonators ($R<1$ mm) using the effects of geometric dispersion. This implies the use of proper Sellmeier equations for $n_{\mathrm{o},\mathrm{e}}\left(\lambda \right)$ [49,55] and relations for the modal frequencies including the size corrections [29], see S2 of supplementary materials for details. It turns out, however, that the main problem with the use of natural PM for the soliton-comb generation is the dispersion issue [56]. In contrast to the above considered cases, the FH and SH dispersion parameters are both negative: $v_1^{\prime }\simeq -0.64$ and $v_2^{\prime }\simeq -1.5$ of ${10}^4$ cm²/s.

Aiming to see the role of the signs of the dispersion parameters, we conducted first some test numerical experiments. For ${\gamma }_{1,2}={10}^7$ s⁻¹, $\alpha =0$, ${\delta }_{1,2}=0$, the absolute values $|{\beta }_1|=0.02$, $|{\beta }_2|=0.01$ and four combinations of the signs of ${\beta }_{1,2}$, we have launched our adiabatic procedure for growing of the antiperiodic soliton states. In two cases with ${\beta }_1{\beta }_2<0$, we got the same known single-soliton solution, see Figure 8. In two other cases with ${\beta }_1{\beta }_2>0$ we obtained no soliton states, and the comb spectra remained undeveloped. Additionally, we varied the values of input parameters $\alpha$, ${\gamma }_{1,2}$, and ${\delta }_{1,2}$ and switched from the antiperiodic to periodic case. No soliton solutions were found as well. The sign of the product ${\beta }_1{\beta }_2$ is thus crucial for the soliton formation.

6. Solutions for FH pumping

The case of FH pumping is important for experiment, because it corresponds to the conditions of SH generation, including the tuning issues. At the first sight, switching to the FH pumping case is trivial – it is sufficient to set $h_2=0$ and $h_1\ne 0$ in sets (3) and (8). The normalized pump amplitude is here ${\eta }_1=\sqrt{2}\mu h_1/{\gamma }_1\sqrt{{\gamma }_1{\gamma }_2}$, and only the periodic states are allowed. While the threshold conditions are more complicated now, they can be evaluated numerically for any set of the input parameters. With these preliminaries, we consider the most important cases [56].

We start from vicinity of the zero walk-off point $2{\lambda }_c\simeq 2698$ nm. The value of the walk-off parameter $\alpha$ is controlled here by the wavelength difference $\delta \lambda ={\lambda }_p-2{\lambda }_c$, and the dispersion coefficients ${\beta }_{1,2}$ can be treated as constants. Setting $\delta \lambda =0$ ($\alpha =0$) and ${\delta }_{1,2}=0$ brings us to the situation where the lowest instability threshold, ${\eta }_1^{\mathrm{th}}=3\sqrt{2}\simeq 4.24$, corresponds to $j=0$, i.e. to auto-oscillations of the dual background [50]. Figure 15 illustrates what happens with increasing ${\eta }_1$ when applying our adiabatic procedure at ${\gamma }_{1,2}=\gamma ={10}^7$ s⁻¹.

Close enough to the threshold we do see auto-oscillations of $|F_0{|}^2$ and $|S_0{|}^2$, see Figure 15a. The oscillation amplitude and period increases and decreases with increasing ${\eta }_1$, respectively, while harmonics with $j,l\ne 0$ remain very small. For ${\eta }_1\gtrsim 10$, these spatial harmonics become notable. Their amplitudes grow gradually, and the number of significant comb lines reaches $N\approx 10$ for ${\eta }_1=11$ with the comb line distance $\delta j=\delta l=1$. The corresponding normalized FH and SH intensity profiles are shown in Figure 15b by lines 1 and 2. Further increase of ${\eta }_1$ results, however, in a sharp transition to an irregular spatial behavior, see subFigure 15c,d. This contrasts strongly with the case of SH pumping. Variations of the input parameters ${\gamma }_{1,2}$ and ${\delta }_{1,2}$ do not change this unfavorable scenario.

With increase of $\left|\alpha \right|$, the threshold values ${\eta }_{1,j}^{\mathrm{th}}$ with $j\ne 0$ quickly decrease and become lower than the internal instability threshold. This is illustrated by Figure 6a relevant to our input parameters. For $\alpha \lesssim 1$, the modes with $\left|j\right|\gg 1$ possess the lowest thresholds. Application of our adiabatic procedure to the case $\left|\alpha \right|=1$, ${\delta }_{1,2}=0$ results, e.g., in the excitation of $\pi /3$-periodic dual state (with six intensity peaks) near the threshold. With increase of ${\eta }_1$, it transforms via a bifurcation to a $\pi /6$ periodic state (at ${\eta }_1\simeq 6.7$). For ${\eta }_1\gtrsim 16$ the system shows again an irregular spatial behavior. The number of significant comb lines N does not exceed 25 for ${\eta }_1\lesssim 16$. Employment of nonzero detunings ${\delta }_{1,2}$ does not help much.

For $\left|\alpha \right|\gtrsim 10$, the threshold values ${\eta }_{1,j}^{\mathrm{th}}\left(\right|\alpha \left|\right)$ approach 2 from above, and selective excitation of harmonics with $\left|j\right|=\left|l\right|=1,2$ becomes possible. However, the thresholds for $\left|j\right|>2$ are only slightly higher. Consider representatively the case $\alpha =11$, ${\delta }_1=1/2$, and ${\delta }_2=0$ relevant to the lowest threshold for $|j,l|=1$. Applying the adiabatic procedure, we have a $\pi$-periodic dual state near the threshold. For ${\eta }_1\approx 7$ and 18 it transforms sequentially via bifurcations to a $\pi /2$- and $\pi /4$-periodic states. Correspondingly, the comb line distance $\delta j=\delta l$ increases from 2 to 4 and, lastly, to 8. Figure 16a,b show complicated steady-state $\pi /4$-periodic asymmetric FH,SH intensity profiles for ${\eta }_1^{max}=20$, while Figure 16c,d exhibit the corresponding comb spectra.

While the steady state found shows a regular spatial behavior, it cannot be regarded either a soliton state or a pronounced comb state.

For $\left|\alpha \right|=50$, the above scenario of decreasing spatial period and increasing the comb line distance with increasing ${\eta }_1$ remains the same. At ${\eta }_1=20$, we have a $\pi /3$ periodic dual state with a modest the comb line number $N=20$. Suppression of SH harmonics with $l\ne 0$ and asymmetry of the intensity profiles and comb spectra are strongly pronounced.

Lastly, we touch the case of natural PM at ${\lambda }_p\simeq 1664$ nm, when small and/or modest values of $\left|\alpha \right|$ are available owing to the geometric dispersion, see S2 of supplementary materials for details. The values of the quality factors $Q_{1,2}$ are expected to be substantially larger than they are at 2698 nm [57] leading to larger absolute values of the dispersion coefficients ${\beta }_{1,2}$. As representative values, we chose $R=1$ mm, ${\beta }_1=-0.04$, and ${\beta }_2=-0.15$, which corresponds to $Q_{1,2}\approx 3\times {10}^8$. The relevant dependences ${\eta }_{1,j}^{\mathrm{th}}\left(\alpha \right)$ are presented in Figure 6b. While they differ from the threshold dependences of Figure 6a, one important feature is the same: It is impossible to excite selectively the modes with $\left|j\right|=1,2$ for $\left|\alpha \right|\lesssim 10$.

A large number of our numerical experiments, conducted at different values of $\alpha$ and ${\delta }_{1,2}$, show persistently the same scenario: Increasing ${\eta }_1$ above the threshold leads (via bifurcations) to dual states with smaller and smaller spatial periods (larger and larger comb line spacing) and, ultimately, to an irregular spatial behavior. The larger $\alpha$, the longer is the interval of ${\eta }_1$ with the regular behavior. On the other hand, the number of significant comb lines N within this interval remains modest; it does not exceed 15. No signs of soliton-comb states were detected.

7. Discussion

Experimental and theoretical studies of dissipative ${\chi }^{\left(2\right)}$ soliton-comb states is a vast, perspective, and almost unexplored research area. It is strongly different from the corresponding Kerr soliton domain, which is, in turn, much broader the studies of conservative Kerr solitons. The ${\chi }^{\left(2\right)}$ nonlinearity is fundamentally different from Kerr’s one resulting in a number of specific qualitative features including duality of the nonlinear states and the necessity to employ the phase matching and to minimize the temporal walk-off. The number of important variable parameters of the system is significantly larger than it is in the Kerr case. Deficit of analytical tools forces researches to rely on numerical methods.

Combining analytical and numerical tools and approaches, we have predicted an unprecedently vast family of dual ${\chi }^{\left(2\right)}$ soliton-comb steady states. Using the term "vast", we mean a continuous dependence of the parameters of the state on all input parameters – the pump power, the modal quality factors, the frequency detunings, and the walk-off. The found states are self-starting, i.e., easily accessible via nor very abrupt switching the pump on and changing the detunings. The soliton-comb states are divided into two topological classes – periodic and antiperiodic depending on the parity of the azimuth number of the pumped SH mode. The steady state solutions can be single- and multi-soliton ones. All our soliton-comb states are dissipative, they cannot be reduced to (obtained from) the known conservative soliton solutions.

Remarkably, the predicted self-starting soliton-comb states are allowed only for the SH pumping. They require opposite signs of the FH and SH dispersion coefficients. This condition is fulfilled near the zero walk-off point ${\lambda }_2\simeq 1349$ nm in LN based resonators. However, the phase matching at this point requires a precise radial poling. Being available, it represents a serious technical problem. For the natural phase matching in LN resonators, the necessary condition is not fulfilled, and soliton-comb states seem to be forbidden. Despite strong efforts, no soliton-comb solutions were found for the FH pumping. This contradicts initial expectations about similarity of the FH and SH pumping cases but has profound reasons behind.

Both intuitive considerations and numerical modeling indicate that substantial differences in the FH-SH group velocities $v_{12}=v_1-v_2$, which are typical of ${\chi }^{\left(2\right)}$ resonators, are unfavorable for the soliton-comb generation. Quantitatively, this means that the combination $\left(\right|v_{12}|/c)\times ({\lambda }_2/R)\times Q_2$ has to be not much larger than 1. This condition is not generally fulfilled in non-structured high-Q resonators. Minimization of the walk-off parameter can be accomplished via the radial poling. Thus, there is one more physical reason to employ the radial poling for the ${\chi }^{\left(2\right)}$ comb generation.

Despite our massive failed efforts to find self-starting soliton-comb steady states for the FH pumping, we cannot claim that all sensible combinations of the input parameters (${\eta }_1$, $\alpha$, ${\beta }_{1,2}$, ${\delta }_{1,2}$, and ${\gamma }_{1,2}$) are exhausted. Also, we cannot claim that the initial assumption of the same transverse structure of the modes within each of FH and SH frequency domains is fully justified. This problem remains still open for both ${\chi }^{\left(2\right)}$ and ${\chi }^{\left(3\right)}$ resonators. It is likely that suppression of the excitation of unwanted transverse modes takes further efforts on the mode management.

A distinctive feature of our approach is employment of a simple classical equivalent of the Hamiltonian quantum formalism of [4,5,17] and many other papers. While quantum electrodynamics aspects of microresonators are important [17], the comb regimes belong indeed to the classical range, where the use of quantum language is redundant. As we have shown, the actual numbers of modal light quanta in the comb spectra are very large.

The effects of ${\chi }^{\left(3\right)}$ nonlinearity in ${\chi }^{\left(2\right)}$ microresonators is one more issue to consider. The possibility to neglect the ${\chi }^{\left(3\right)}$ nonlinearity is largely controlled by the modal quality factors. Let us consider for simplicity the SH pumping case. The pump power thresholds ${\mathcal{P}}_{\mathrm{th}}^{\left(2\right)}$ and ${\mathcal{P}}_{\mathrm{th}}^{\left(3\right)}$ for the parametric instability caused by the ${\chi }^{\left(2\right)}$ and ${\chi }^{\left(3\right)}$ nonlinearities are roughly proportional to $1/Q^3$ and $1/Q^2$, such that for $Q_{1,2}\gtrsim {10}^6$ the effects of Kerr nonlinearity are expected to be negligible. For modest quality factors, $Q\lesssim {10}^5$, which are typical for ring waveguides, the situation can be more complicated. An interplay between the ${\chi }^{\left(2\right)}$ and ${\chi }^{\left(3\right)}$ nonlinear processes becomes possible [42].

The above considered soliton-comb steady states in ${\chi }^{\left(2\right)}$ microresonators are not the first predicted ones and known in the literature. To the best of our knowledge, the other known stable dissipative soliton states [40,41,42,43] belong to the periodic ones. Furthermore, they are not self-starting, i.e., not easily accessible. The single-soliton numerical solutions of [40] are relevant to the FH pumping at the zero walk-off point. Both bright-bright and dark-dark locally stable solitons were exhibited. Similar numerical soliton solutions were found in [41,42]. Additionally, quasi-solitons at very large detunings, $|{\delta }_{1,2}|\gg 1$, were predicted. Remarkably, the soliton-comb states of [40] are allowed for the same and opposite signs of the dispersion coefficients. Generation of flat-top solitons via an additional spatial modulation of the pump amplitude was found numerically in the SH pumping case for opposite signs of the dispersion coefficients [43]. A considerable amount of soliton related models was proposed for bulk parametric oscillators. They are beyond the scope of this review.

The first experimental attempts to get ${\chi }^{\left(2\right)}$ combs in high-Q LN based microresonators dealt with the FH pumping, the natural phase matching, and large walk-off [37,38]. Neither ${\chi }^{\left(2\right)}$ solitons no broad comb spectra similar to the Kerr soliton spectra, were observed. This is in agreement with our theoretical results on the necessary conditions for realization of such regimes. The authors of [39] claim for discovery of Pockels soliton microcomb in an aluminium nitride microring resonator at SH pumping. The parameters used correspond to modest quality factors $Q_{1,2}<{10}^6$, a large walk-off parameter $\left|\alpha \right|\gg 1$, and the threshold power ${\mathcal{P}}_{\mathrm{th}}^{\left(2\right)}$ only slightly below ${\mathcal{P}}_{\mathrm{th}}^{\left(3\right)}$. The observed broad FH comb spectrum was practically symmetric, which is the fingerprint of the absence of dominating contribution of the ${\chi }^{\left(2\right)}$ nonlinearity. The observed and numerically modeled combs are due to the Kerr nonlinearity. The role of quadratic nonlinearity is mostly in initiation of the parametric instability. Thus, net ${\chi }^{\left(2\right)}$ soliton-comb states are still waiting for experimental discovery and exploration.

8. Conclusions

A vast family of stable soliton-comb steady states is accessible in ${\chi }^{\left(2\right)}$ microresonators via not very abrupt switching the monochromatic pump on. The necessary conditions for the soliton-comb generation are (i) pumping into a SH resonator mode and ($ii$) opposite signs of the FH and SH dispersion coefficients. Also, minimization of the group velocity (free spectral range) difference for the phase matched FH and SH modes is highly desirable. For LN based resonators, the latter can be accomplished simultaneously with ($ii$) using proper radial poling of the resonator. Pumping into FH resonator modes is not favorable for realization of the soliton-comb states. Employment of high-Q resonators ($Q>{10}^7$) ensures low generation thresholds (within the $\mu$W range) and negligibly small influence of the Kerr nonlinearity.