Resonant Kushi-comb-like multi-frequency radiation of oscillating two-color soliton molecules

Emission of radiation by physical objects is omnipresent in physical systems. For instance, in classical electrodynamics, oscillating and linear accelerating point charges generate electric and magnetic fields. Charged particles confined to circular orbits emit synchrotron radiation, decelerating charged particles produce Bremsstrahlung, and oscillating electric dipoles cause dipole-radiation characteristic of radiative transitions in atoms and molecules [1].

In the field of nonlinear optics many processes can be described by wave equations of the type of the nonlinear Schrödinger equation (NSE) [2, 3], which, in the integrable case, supports localized field pulses given by solitary waves [4]. Such localized field pulses resemble particles, and can, under the right circumstances, sustain periodic variations of their amplitude and width. The emission of resonance radiation caused by such oscillations has been studied in different nonlinear-optics settings. For example, considering a NSE in which a varying dispersion coefficient perturbs transmitted pulses periodically, dispersion-managed solitons were found to radiate through a resonance mechanism. The underlying parametric resonance and the radiative decay of dispersion-managed solitons were studied for rapidly varying dispersion maps [5], and for weakly varying dispersion maps [6]. In dispersion oscillating fibers, such parametric resonances were investigated numerically and also experimentally for solitons (in case of anomalous dispersion) and shock fronts (in case of normal dispersion) [7], demonstrating the generation of resonant radiation (RR) at multiple resonant frequencies. In a standard NSE with added third-order dispersion, higher-order solitons were shown to support a quite similar radiation mechanism [8]. Therein, the periodic variation of the solitons peak intensity enabled the coupling of the localized state to continuum modes on the opposing side of the zero-dispersion point. As a result, RR in form of multi-frequency comb-like spectral bands is emitted. The generation of such RR was shown to also proceed in presence of the Raman-effect and pulse self-steepening [8]. In the case of the Lugiato–Lefever equation with third-order dispersion, i.e. a variant of the NSE with added damping and driving, a Hopf-bifurcation can result in oscillating dissipative solitons that shed RR in a series of well resolved spectral peaks [9].

In analogy to quantum mechanics, where molecules are formed as bound-states of multiple atoms, nonlinear optics supports various types of molecule-like bound states of multiple field pulses. For instance, usual soliton molecules which consist of two identical subpulses separated by a fixed time-delay, can be realized in a dispersion managed NSE [10, 11]. Mutually bound solitons also arise in models of coupled NSEs [12–19]. Specifically, for twin-core fibers with higher order dispersion, bound solitons near the zero-dispersion point were studied and shown to sustain periodic center-of mass and amplitude variations with different oscillation periods, giving rise to RR with different frequencies [20]. Further, dissipative optical soliton molecules can be generated in mode-locked fiber lasers [21, 22]. Quite recently, a different kind of molecule-like pulse compound has been demonstrated theoretically [23, 24]. It consists of two subpulses at distinctly different center frequencies, which superimpose to form a single complex in the time-domain. An attractive, cross-phase modulation induced potential holds these two-color pulses compounds together [23]. The binding mechanism that enables their stable propagation is different from that of usual soliton molecules. It requires a propagation constant with (at least) two separate domains of anomalous dispersion within which group-velocity matched co-propagation of pulses is possible. A strong incoherent attractive interaction between the subpulses of such a two-color soliton, mediated across (at least) two zero-dispersion points, is ensured if the group-velocity mismatch between the pulses is negligible [25]. The emphasis on group-velocity matching is similar as for the repulsion between a soliton and a dispersive wave mediated across a single zero-dispersion point [26–30], where quasi group-velocity matching is vital to ensure a strong and efficient interaction mechanism. In terms of a NSE with added negative fourth-order dispersion, two-color soliton molecules were identified as members of a large family of generalized dispersion Kerr solitons [24]. Such two-color solitons were recently verified experimentally in mode-locked laser cavities [31, 32]. Two-color soliton microcomb states of similar structure also exist in the framework of the Lugiato–Lefever equation [33, 34]. The strong mutual interaction of the subpulses of such molecule-like bound states [35] renders them ideal for studying the effects of oscillations exhibited by coupled localized states in nonlinear optics systems.

Here, we discuss the rich dynamical behavior of such two-color soliton molecules, consisting of group-velocity matched pulses located in distinct domains of anomalous dispersion, separated by a vast frequency gap. We show that nonstationary dynamics of the subpulses results in emission of RR. The location of the newly generated multi-peaked spectral bands can be understood by extending existing approaches [7, 9, 20, 36, 37] to two-frequency pulse compounds [38]. Due to the manifold of two-frequency pulse compounds with differing substructure, their emission spectra manifest in various complex forms. Viewing a two-color soliton molecule as a superposition of two distinct subpulses immediately suggests the existence of two different types of dynamics: first, oscillations of the subpulses with respect to a common center in time, and second, pure amplitude oscillations of the subpulses. In previous studies of two-color soliton molecules [23, 25, 35], a combination of both types of oscillations has generally been found. Therein, the first type of oscillation was observed to be heavily damped [25, 35, 38]. In the presented work we focus on pure amplitude oscillations, as they represent the persistent oscillation mode of two-frequency soliton molecules.

The article is organized as follows. In section 2 we detail the propagation model used for our theoretical investigations and we specify initial conditions that directly generate special soliton-like two-color soliton pairs [39], which we below also refer to as two-color soliton molecules (or soliton molecules). This approach surpasses previous 'seeding' procedures, which result in a largely uncontrolled generation of localized two-color pulse compounds. In section 3 we will perturb such a two-color soliton molecule by systematically increasing its initial amplitude. In analogy to usual single-pulse nonlinear Schröinger solitons, where increasing the amplitude results in an oscillating pulse, this triggers amplitude oscillations of the entire pulse compound, yielding a vibrating soliton molecule. We numerically demonstrate that the mutual driving of the subpulses results in the generation of resonant multi-frequency radiation and we compare the observed location of the resonances to analytic predictions obtained from a phase-matching analysis. We demonstrate that symmetries, inherent in the propagation setting, can lead to degenerate resonance conditions, and we show how existing degeneracies can be lifted by a perturbation of the propagation constant. Section 4 concludes with a summary.

2.1. Propagation model

For studying the propagation dynamics of oscillating two-color pulse compounds, we consider the modified NSE

$$\begin{equation} i \partial_z A = \left(\frac{\beta_2}{2} \partial_t^2 - \frac{\beta_4}{24}\partial_t^4\right) A - \gamma |A|^2 A, \end{equation} \tag{ 1 }$$

for a complex-valued field $A\equiv A(z,t)$, on a periodic time-domain of extent T with boundary conditions $A(z,-T/2) = A(z,T/2)$. In equation (1), $\beta_2\gt0$ (in units of $\mathrm{fs^2\,\mu m^{-1}}$) is a positive-valued group-velocity dispersion coefficient, $\beta_4\lt0$ ($\mathrm{fs^4\,\mu m^{-1}}$) is a negative-valued fourth-order dispersion coefficient, and γ ($\mathrm{W^{-1}\,\mu m^{-1}}$) is the nonlinear coefficient. For the discrete set of angular frequency detunings $\Omega\in \frac{2\pi}{T}\mathbb{Z}$, the expressions

$$\begin{align} A_\Omega(z)&=\mathsf{F}[A(z,t)] \equiv \frac{1}{T} \int_{-T/2}^{T/2} A(z,t)\,\mathrm{e}^{i\Omega t}~\mathrm{{d}}t, \end{align} \tag{ 2 }$$

$$\begin{align} A(z,t) &= \mathsf{F}^{-1}[A_\Omega(z)] \equiv \sum_{\Omega} A_\Omega(z)\,\mathrm{e}^{-i\Omega t}, \end{align} \tag{ 3 }$$

specify forward (equation (2)), and inverse (equation (3)) Fourier transforms that relate the field envelope $A(z,t)$ to the spectral envelope $A_\Omega(z)$. Using Parseval's identity for these transforms, the total energy in both domains is given by

$$\begin{equation} E(z) = \int_{-T/2}^{T/2} |A(z,t)|^2~\mathrm{{d}}t = T \sum_\Omega |A_\Omega(z)|^2, \end{equation} \tag{ 4 }$$

wherein $|A(z,t)|^2$ ($\mathrm{W = J/s}$) specifies the instantaneous power, and $|A_\Omega(z)|^2$ ($\mathrm{W}$) the power spectrum.

2.2. Propagation constant

Using the identity $\partial_t^n\,\mathrm{e}^{-i \Omega t} = (-i\Omega)^n\,\mathrm{e}^{-i\Omega t}$, the frequency-domain representation of the propagation constant is given by the polynomial expression $\beta(\Omega) = \frac{\beta_2}{2}\Omega^2 + \frac{\beta_4}{24} \Omega^4$. The inverse group-velocity of a mode at detuning Ω is given by $1/v_g(\Omega) = \beta_1(\Omega)\equiv \partial_\Omega \beta(\Omega)$, and the group-velocity dispersion is $\beta_2(\Omega)\equiv \partial_\Omega^2 \beta(\Omega)$. Assuming $\beta_2\gt0$, and $\beta_4\lt0$, the group-velocity dispersion exhibits a downward parabolic shape with zero-dispersion points defined by $\beta_2(\Omega) = 0$ at $\Omega_\mathrm{{Z1},\mathrm{Z2}} = \mp \sqrt{2\beta_2/|\beta_4|}$. Subsequently we set $\beta_2 = 1\,\mathrm{fs^2\,\mu m^{-1}}$, and $\beta_4 = -1\,\mathrm{fs^4\,\mu m^{-1}}$, resulting in the dispersion characteristics shown in figures 1(a)–(c). We further set $\gamma = 1\,\mathrm{W^{-1}\,\mu m^{-1}}$. Let us note that the phenomena reported below are not limited to the particular choice of the above parameters.

Zoom In Zoom Out Reset image size

2.3. Propagation algorithm

The numerical simulations reported in section 3 are performed using an algorithm with adaptive stepsize control, implementing the 'Conservation quantity error' method (CQE) [40, 41], which employs a conservation law of the underlying propagation equation to guide selection of the propagation stepsize h. Specifically, we use the relative error

$$\begin{equation} \delta_{E}(z) = \frac{|E(z+h) - E(z)|}{E(z)}, \end{equation} \tag{ 5 }$$

where E is the energy (4), conserved by the NSE (1). The CQE method keeps δ_E within the goal error range $[0.1\,\delta_\textrm{{G}}, \delta_\textrm{{G}}]$, specified by the local goal error $\delta_\textrm{{G}}$, by decreasing h when necessary while increasing h when possible. For our numerical experiments we use $\delta_\textrm{{G}} = 10^{-10}$. To advance the field from position z to z + h, we use the fourth-order 'Runge-Kutta in the interaction picture' method [42].

2.4. Spectrograms

To assess the time-frequency interrelations within the field A at a selected propagation distance z, we use the spectrogram [43, 44]

$$\begin{equation} P_{S}(t,\Omega) = \frac{1}{2 \pi} \left|\int_{-T/2}^{T/2} A(z,t^{\prime})h(t^{\prime}-t) \mathrm{e}^{-i \Omega t^{\prime}}~\mathrm{ d}t^{\prime}\right|^2. \end{equation} \tag{ 6 }$$

Therein, $h(x) = \exp(-x^2/2\sigma^2)$ specifies a Gaussian window function with root-mean-square width σ, used to localize the field in time.

2.5. Initial conditions

Given the modified NSE (1), two-color pulse compounds can be 'seeded' by using initial conditions consisting of two superimposed fundamental solitons in the form

$$\begin{equation} A_0(t) = a_1 \,\mathrm{{sech}}(t/t_1)\, \mathrm{e}^{-i \Omega_1 t} + a_2\, \mathrm{{sech}}(t/t_2)\,\mathrm{e}^{-i\Omega_2 t}, \end{equation} \tag{ 7 }$$

where $a_n^2 = |\beta_2(\Omega_n)|/(\gamma \,t_n^2)$, with $n\in(1,2)$, are the peak intensities of the two subpulses which are subject to the group-velocity matching condition $\beta_1(\Omega_1) = \beta_1(\Omega_2)$. Upon propagation, such an initial condition will evolve into a narrow localized pulse which can exhibit internal dynamics, reminiscent of molecular vibrations, and is accompanied by free radiation [23]. Recently we demonstrated that if the spectral extent of each subpulse is small compared to the frequency gap separating them, the mutual interaction of both pulses is effectively incoherent, allowing to construct special solutions in the form of two-color soliton pairs [39]. Considering the propagation constant used in equation (1), an entirely symmetric two-color soliton pair can be constructed by choosing both subpulse loci at the two group-velocity matched frequencies $\Omega_{1,2} = \mp\sqrt{6\beta_2/|\beta_4|}$, highlighted in figures 1(a)–(c). They exhibit $\beta_1(\Omega_{1}) = \beta_1(\Omega_2) = \beta_1(0) = 0$ (figure 1(b)), and $\beta_2(\Omega_1) = \beta_2(\Omega_2) = -2 \beta_2 \lt 0$ (figure 1(c)). Assuming two identical subpulse envelopes in equation (7), realized by $t_1 = t_2\equiv t_0$ and $a_1^2 = a_2^2 = 2\beta_2/(3 \gamma t_0^2)$ [39], then yields

$$\begin{equation} A_0(t)= F(0,t)\,\cos\left( \sqrt{2 \beta_2/|\beta_4|}t\right), \end{equation} \tag{ 8 }$$

wherein

$$\begin{equation} F(z,t)= \sqrt{P_0}\,\mathrm{{sech}}(t/t_0)\,\mathrm{e}^{i \kappa z} \end{equation} \tag{ 9 }$$

specifies a pulse-enclosing envelope with peak power $P_0 = 8 \beta_2/(3 \gamma t_0^2)$ and nonlinear wavenumber $\kappa = \beta_2/t_0^2$ [39]. These two-color soliton pairs correspond to the fundamental metasoliton describing generalized dispersion Kerr solitons, recently obtained in terms of a multi-scales analysis [24].

2.6. Fundamental two-color soliton molecule

To investigate the radiation that emanates from an oscillating two-color soliton molecule, we need to realize stable z-periodic dynamics. Specifically, in section 3.1 we demonstrate how to induce amplitude and width oscillations of soliton molecules and we study their characteristics. In section 3.2 we demonstrate that oscillating soliton molecules shed RR in the form of multi-peaked spectral bands, and we compare our numerical results to predictions based on a phase matching analysis. In section 3.3 we show that possible degeneracies in the resulting radiation can be lifted by perturbations.

3.1. Oscillating two-color soliton molecules

3.1.1. Amplitude oscillations

So as to excite amplitude and width oscillations of two-color soliton molecules, we consider initial conditions of the form of equations (8) and (9), and include an additional dimensionless factor N through $P_0\to P_0^{\prime} \equiv N^2 P_0$, allowing to increase their initial amplitude. The choice N = 1 yields the basic soliton molecule of equation (8) (figure 1(f)), for which the pulse-enclosing envelope F is coined the fundamental metasoliton in reference [24]. In analogy with the nomenclature of usual single-pulse nonlinear Schrödinger solitons, choosing N > 1 define higher-order metasolitons [24], resulting in oscillating soliton molecules. For brevity, we will refer to N simply as the order of the soliton molecule. Figure 2 demonstrates such oscillating soliton molecules for N = 1.1 (figure 2(a)), N = 1.4 (figure 2(b)), and N = 1.65 (figure 2(c)). If the initial perturbation of the two-color soliton molecule is small (N = 1.1), the variation of the scaled peak intensity is largely harmonic (figure 2(d)), i.e. a single spatial Fourier mode describes the observed intensity variation (figure 2(g)). As evident from figures 2(e), (h) and (g), (i), increasing the molecule order results in nonlinear peak intensity oscillations, and, consequently, an increasing number of spatial Fourier modes. Thereby, the variation of the molecule peak intensity is characterized by an oscillation period Λ that decreases with increasing order N.

Zoom In Zoom Out Reset image size

3.1.2. Amplitude oscillation wavelength

To quantify these dynamic amplitude and width variations we performed a series of numerical simulations for soliton molecules with different initial duration and order in the range $[1.05,2]$, see figure 3. For N → 1, the observed oscillation period agrees well with that of standard nonlinear Schrödinger solitons of order N governed by a standard NSE, indicated by the solid line (labeled C) in figure 3. This finding is reasonable, since, for N ≈ 1, the coupling between the molecules subpulses remains effectively incoherent so that a given subpulse can be described by a standard NSE with modified coefficients [39]. For N → 2, the oscillation period of standard NSE solitons approaches the value $z_\textrm{{sol}} = (\pi/2) z_0$, i.e. the soliton period of a higher-order (single) soliton with dispersion length z₀, indicated by the horizontal line in figure 3. We find that the oscillation period of soliton molecules with large duration $t_0 = 12\,\mathrm{fs}$, labeled B in figure 3, are in excellent agreement with this behavior. This finding is consistent with results obtained for a N = 2 metasoliton of duration $t_0 = \sqrt{2 |\beta_4|/(3 \beta_2 \epsilon^2)}\approx 16.3\,\mathrm{fs}$ ( = 0.05), studied in reference [24]. In contrast, soliton molecules with smaller duration $t_0 = 6\,\mathrm{fs}$, labeled A in figure 3, deviate from this behavior for $N\gtrapprox 1.6$. From our numerical experiments we found that for soliton molecules with small duration, even small amplitude oscillations suffice to cause a molecule to adiabatically change its shape, affecting its oscillation period and explaining the observed trend of the data. Overall, the close correspondence between the oscillatory behavior of soliton molecules and standard NSE solitons allows to directly transfer approximate descriptions of the evolution of near-fundamental standard NSE solitons, obtained in terms of variational approaches [45, 46], to the present case. A simple qualitative model of the observed subpulse oscillations can be derived on basis of a coupled system of variational equations, governing the evolution of the parameters of a hyperbolic-secant trial pulse [46]. From this, an equation of motion can be obtained that captures the variation of the pulse width about some preferred, stationary value (see appendix A). In lowest order approximation, it describes a harmonic oscillation. The corresponding wavelength

$$\begin{equation} \Lambda^\mathrm{{HO}}(N) = \frac{\pi^2}{\sqrt{2} N^2} z_0, \end{equation} \tag{ 10 }$$

derived in appendix A, is in qualitative agreement with our numerical results, see the dashed line (labeled D) in figure 3. A less intuitive, but more accurate estimate of the oscillation period can be obtained by using Gaussian trial functions to approximate NSE solitons in a variational analysis of the standard NSE, yielding the approximate oscillation period [45]

$$\begin{equation} \Lambda^\mathrm{{var}}(N) =\frac{\pi \sqrt{2} N^2}{(\sqrt{2}N^2-1)^{3/2}} z_0, \end{equation} \tag{ 11 }$$

which is in good agreement with our numerical results, as shown by the dash-dotted line (labeled E) in figure 3. A fit of the molecule data to a polynomial model of the form

$$\begin{equation} \Lambda^\mathrm{{fit}}(N) = z_0 \sum_{n=0}^3 c_n N^{-n}, \end{equation} \tag{ 12 }$$

yields the parameters $(c_0,c_1,c_2,c_3) = (-15.1, 64.3, -93.2, 56.2)$ for soliton molecules with initial duration $t_0 = 6\,\mathrm{fs}$, shown as short dashed line (labeled F) in figure 3. For molecules with $t_0 = 12\,\mathrm{fs}$, a fit yields $(c_0,c_1,c_2,c_3) = (-3.8, 24.7, -47, 39)$.

Zoom In Zoom Out Reset image size

3.2. Degenerated multi-frequency radiation

Below we show that the nonstationary dynamics of soliton molecules, which is reminiscent of molecular vibrations, results in emission of resonant multi-frequency radiation.

3.2.1. Example for N = 1.65

As a first example, we consider a soliton molecule with duration $t_0 = 6\,\mathrm{fs}$, and order N = 1.65. Figure 4(a) shows the time-domain intensity on a logarithmic scale (cf figure 2(c)), indicating that, in response to the initial perturbation, the localized molecule sheds radiation along either t-direction. The periodic variation of the spectrum, accompanying the amplitude and width oscillations of the molecule with wavelength $\Lambda/z_0 = 2.23$, is shown in figure 4(b). Let us note that the oscillation of the molecule is clearly visible in the close-up view in figure 4(a). Distinct, comb-like bands of frequencies can be seen to develop in the vicinity of the subpulse loci Ω₁ and Ω₂. These newly generated spectral bands appear particularly pronounced when the molecule assumes its maximum time-domain width during an oscillation period, so that the resonances are not drowned by the spectrum of the molecule, demonstrated for the specific choice $z/z_0 \approx 28.6$ in figures 4(b) and (d). The spectral location of this RR can be predicted in terms of wave number matching conditions, derived for the dynamically evolving, oscillating molecule state. For clarity, below we summarize the derivation of the resonance conditions detailed in reference [38]. Considering an oscillating soliton molecule, consisting of two subpulses U₁ and U₂, we can account for amplitude and width oscillations by representing both in terms of a Fourier series [38]

$$\begin{equation} U_n(z,t) = \sum_{m} C_{nm}(t)\,\exp\left(i\kappa_n z + im\frac{2\pi}{\Lambda}z\right), \quad n\in(1,2),~m\in\mathbb{Z} \end{equation} \tag{ 13 }$$

wherein C_nm are expansions coefficients, κ_n are the wavenumbers that govern the z-propagation of the subpulses, and Λ is the wavelength of the periodic z-variation of the molecule. Below, we abbreviate the wavenumbers associated with the harmonics of the molecules z-oscillation wavelength Λ by $K_m = 2\pi m/\Lambda$. Under the assumption that these subpulses satisfy a pair of coupled NSEs [25, 39], and by assessing how dispersive radiation is excited by their mutual periodic driving [38], the resonance conditions

$$\begin{align} D_n(\Omega_\mathrm{RR}) - \kappa_n &= K_m, &\qquad n\in(1,2),~m\in\mathbb{Z} \end{align} \tag{ 14 }$$

$$\begin{align} D_n(\Omega_\mathrm{RR}) -2 \kappa_n + \kappa_{n^{\prime}} &= K_m, &\qquad n, n^{\prime} \in (1,2),~n\neq n^{\prime},~m\in\mathbb{Z} \end{align} \tag{ 15 }$$

can be obtained. Therein, $D_n(\Omega)\equiv \beta(\Omega)-\beta(\Omega_n) -\beta_1(\Omega)(\Omega-\Omega_n)$, for $n = 1,2$, defines the dispersion profile in the reference frame of each of the subpulses, and $\Omega_{RR}$ specifies frequencies at which RR is excited. The individual subpulse wavenumbers are related to the subpulse amplitudes a₁ and a₂ through

$$\begin{align} \kappa_1 &= \frac{\gamma}{2} (a_1^2 + 2 a_2^2), \end{align} \tag{ 16 }$$

$$\begin{align} \kappa_2 &= \frac{\gamma}{2} (a_2^2 + 2 a_1^2) + \beta(\Omega_2)-\beta(\Omega_1). \end{align} \tag{ 17 }$$

Zoom In Zoom Out Reset image size

Equation (14) yields resonance conditions for the generation of Cherenkov radiation emitted by each of the molecules subpulses, related to the mth harmonic of the driving period Λ. Equation (15) defines further resonance conditions indicative of four-wave mixing (FWM) processes that involve both subpulses.

In reference [38], the subpulse wavenumbers κ₁ and κ₂ had to be extracted a posteriori from the simulated propagation dynamics. In the present case, we can establish a reasonable a priori estimate by assuming the molecule to be close to a two-color soliton pair with subpulse amplitudes $a_1 = a_2 = \sqrt{2 \beta_2/ (3\gamma t_0^2)}$ [39] and $\beta(\Omega_1) = \beta(\Omega_2) = 3\beta_2^2/(2 |\beta_4|)$, giving

$$\begin{equation} \kappa_1 = \kappa_2 = \frac{\beta_2}{t_0^2}. \end{equation} \tag{ 18 }$$

Consequently, the resonance conditions (14) and (15) are degenerate. The wavenumber profile in the reference frame of the subpulse at Ω₁, given by $D_1(\Omega)-\kappa_1$, is shown in figure 4(c). In the reference frame of the subpulse at Ω₂, the wavenumber profile $D_2(\Omega)-\kappa_2$ looks identical. This is a consequence of the symmetry of the propagation constant about the point $\Omega = 0$, and the symmetric choice $-\Omega_1 = \Omega_2 = \sqrt{6}\,\mathrm{rad\,fs^{-1}}$ of the subpulse loci. A graphical solution of the resonance conditions (14) and (15) can be obtained by superimposing the wavenumber profile with the wavenumbers K_m , associated with the harmonics of the molecules z-oscillation wavelength. This is demonstrated in figure 4(c), where K_m for $m = -8\ldots 1$ are included as horizontal dashed lines. The crossing points, emphasized by dots in figure 4(c), indicate frequency loci $\Omega_\mathrm{RR}$ at which the coupling of energy into dispersive waves is expected to be efficient. As shown by the vertical lines in figures 4(c) and (d), we find that these resonant frequencies are in excellent agreement with our numerical results. The large number of resonant frequencies, excited and continuously fed during propagation of the molecule, is related to the large number of modes required by a Fourier decomposition of the nonlinear vibration of the molecule for N = 1.65 (cf figures 2(f) and (i)). Molecules of order N = 1.4 yield fewer and less pronounced resonances. As shown by the spectrogram figure 6(a), all resonantly excited modes leave the molecule at a characteristic group-velocity, characteristic for their location in the vicinity of Ω₁ and Ω₂: all resonantly excited modes in $[0,\Omega_1]$ and $[\Omega_2,\infty]$ propagate towards larger positive times; those in $[0,\Omega_2]$ and $[\Omega_1,\infty]$ propagate towards increasingly negative times. We can simplify the above resonance condition to obtain a simple expression that specifies the location of the first few resonances reasonably well: by making a parabolic approximation of the propagation constant in the vicinity of Ω₁ and Ω₂, and by neglecting the contribution of the subpulse wavenumbers κ₁ and κ₂, we find

$$\begin{equation} \Omega_\mathrm{{RR},\textit{m}}^{(1,2),\pm} \approx \Omega_{1,2} \pm \sqrt{\frac{|K_m|}{\beta_2}}. \end{equation} \tag{ 19 }$$

The approximate positions of the resonances for $m = -1$ are indicated by vertical dashed lines in figure 4(d).

3.2.2. Example for N = 2

According to equation (19) it should be possible to control the separation between consecutive resonances by either decreasing β₂, resulting in a more 'flat' dispersion profile $D_1(\Omega)$ in figure 4(c), or by increasing the wavenumbers K_m , resulting in more widely spaced horizontal dashed lines in figure 4(c). To verifyy this prediction we next consider a soliton molecule of order N = 2. Since such a molecule exhibits $\Lambda(N = 2)\lt\Lambda(N = 1.65)$ (figure 3), we expect to observe resonant frequencies with larger spacing as in the previously considered case. As evident from figures 4(e) and (f), the evolution of the molecule proceeds in two distinct stages: during the initial propagation stage $z/z_0 \lessapprox 10$, a large amount of dispersive radiation emanates from the molecule; then, for $z/z_0\gtrapprox 10$, the molecule seems to have reached a steadily oscillating stage. The initial burst of energy results in an adiabatic reshaping of the molecule in course of which its temporal width decreases to about $t_0\approx 4\,\mathrm{fs}$, and its spectrum broadens. Consequently, the molecule also develops a notable spectral density within the domain of normal dispersion (figure 4(f)). The wavelength of the z-periodic oscillations also decrease slightly from initially $\Lambda/z_0\approx 1.1$ to $\Lambda/z_0\approx 0.92$ for $z/z_0\gtrapprox 10$. This causes the wavenumbers K_m , indicated by horizontal dashed lines in figure 4(g), to shift downwards, and the corresponding resonant frequencies to shift away from the nearest subpulse location. The associated spectral lines then appear broadened (figure 4(h)). A similar observation has been made in a study of the resonance mechanism of higher-order solitons in a NSE with added third order dispersion [8]: adding the Raman-effect, which causes a self-frequency-shift of solitons, was found to wash out the resonantly excited radiation bands. The approximate positions of the resonances for $m = -1$ are again indicated by vertical dashed lines in figure 4(h). Contrasting the dynamics for N = 2 with that for N = 1.65 provides an unexpected observation: since the possible driving wavenumbers K_m are now more widely spaced, see the horizontal dashed lines in figure 4(g), the molecule can support only fewer resonant frequencies; as a result, its amplitude and width oscillations appear more harmonic. As evident from the spectrogram figure 6(b), due to the aforementioned decrease of Λ for increasing z, some resonances, which are present during the initial 'nonlinearly oscillating' propagation stage ($z/z_0\lessapprox 10$), die out and are absent in the later 'harmonically oscillating' propagation stage ($z/z_0\gtrapprox 10$). Let us note that for both cases considered above, the particular resonance for m = 0, indicating the excitation of the usual Cherenkov radiation in absence of amplitude and width oscillations [47], is not excited.

In both cases studied above, we find that the amplitude and width oscillations of the soliton molecules are quasistationary and decay only slowly due to the emission of RR. Further, let us note that taking longer pulses will result in more narrow spectra in the vicinity of the center frequencies of both subpulses. This will result in a larger dispersion length and soliton period, see equations (10) and (12). As a consequence, the wavenumber associated with the z-periodic variation of a soliton molecule will decrease. According to equations (14) and (15), this will result in more narrowly spaced resonances in the spectrum.

3.3. Lifting the degeneracies by perturbations

In order to lift the degeneracy of the resonance conditions (14) and (15), which occurs for $\kappa_1 = \kappa_2$, equations (16) and (17) suggest two possibilities: either setting $a_1\neq a_2$ (see appendix B), i.e. assigning different amplitudes for both subpulses, or devising a propagation scenario for which $\beta(\Omega_1)\neq \beta(\Omega_2)$. Subsequently, we opt for the latter and assess the impact of a non-symmetric propagation constant on the dynamics of a soliton molecule. Specifically, we modify the propagation constant to $\beta(\Omega) = \frac{\beta_2}{2}\Omega^2 + \frac{\beta_3}{6} \Omega^3 + \frac{\beta_4}{24} \Omega^4$ using $\beta_3 = -0.2\,\mathrm{fs^3\,\mu m^{-1}}$. This breaks the symmetry inherent in the propagation scenarios studied above. In figure 5 we show how a soliton molecule with N = 1.65, and subpulse center frequencies again at $-\Omega_1 = \Omega_2 = \sqrt{6}\,\mathrm{rad\,fs^{-1}}$, responds to such a perturbation. In this case, the effective values of second order dispersion felt by both subpulses are $\beta_2(\Omega_1)\approx-1.51\,\mathrm{fs^2\,\mu m^{-1}}$, and $\beta_2(\Omega_2)\approx -2.49\,\mathrm{fs^2\,\mu m^{-1}}$. We then start with a (fundamental) two-color soliton pair in the form of equation (7) with $a_1 = \sqrt{(2\alpha-1)|\beta_2(\Omega_1)|/(3\gamma t_0^2)}$ and $a_2 = \sqrt{(2\alpha^{-1}-1)|\beta_2(\Omega_2)|/(3\gamma t_0^2)}$, with $\alpha = \beta_2(\Omega_2)/\beta_2(\Omega_1)\approx 1.65$ [39], and increase its amplitude by the factor N. As a first observation, we note that the molecule acquires an additional shift in its group-velocity due to $\beta_3\neq 0$. In figure 5(a) we compensated for the resulting drift by transforming to a reference frame with velocity $v_0 = -1.877\,\mathrm{\mu m\,fs^{-1}}$, within which the molecule appears stationary. Further, the spectrum in figures 5(b) and (d) is now very different from that in figures 4(b) and (d). The soliton molecule exhibits a non-symmetric spectrum, consisting of a strong pulse, centered at Ω₁, and a weaker pulse, centered at Ω₂. In such a scenario, due to a strong entanglement of the subpulses [35], we expect that the strong pulse drives the weak pulse so that the wavelength of the intensity and width oscillations are equal for both pulses. From our numerical simulation we find $\Lambda/z_0 = 3.2$. Below we perform a phase-matching analysis in the frame of reference of the stronger subpulse (i.e. the one at Ω₁). In this case, we have to resort to measuring the subpulse wavenumbers κ₁ and κ₂ a posteriori from the simulated dynamics. A graphical solution of the resonance conditions (14) for $m = -5\ldots 5$ shows that now there are also multi-peaked spectral bands of resonances excited around m = 0, indicated by the red short-dashed line in figure 5(c). A filtered view of the RR at $\Omega\approx 0.6\,\mathrm{rad\,fs^{-1}}$, enclosed by a box (labeled A) in figure 5(b), is shown in the box (labeled A) in figure 5(a). As evident from this filtered view, the RR is emitted in a pulse-wise fashion and with decreasing intensity for increasing propagation distance [48]. The phase-relationship between theses separated pulses then yields the multi-peaked spectral band in figure 5(d). Additional bands of resonances, associated with Cherenkov resonances, are visible at $\Omega\approx \Omega_1$ and $\Omega\approx 3.5\,\mathrm{rad\,fs^{-1}}$. For clarity, we provide a graphical solution of the phase-matching conditions (14) for $m = -5\ldots 5$, only. Similar effects for perturbation by β₃ have previously been reported for higher-order solitons in the NSE [8], and dissipative solitons in fiber cavities [9]. Additional FWM resonances arising from equation (15) can be identified at $\Omega\approx -3.8\,\mathrm{rad\,fs^{-1}}$. In contrast to an earlier approach, where such resonances were caused by the nonlinear mixing of solitons and dispersive waves [37], we here find that these resonances are caused by a mixing of the subpulses of a soliton molecule.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

All of the resonantly generated radiation emanates from the molecule in the vicinity of the Cherenkov, and FWM resonances specified by equations (14) and (15). The pulse-wise emission is clearly evident in the spectrogram figure 6(c). In figure 6(c) we did not compensate for the non-vanishing group velocity of the molecule, so that, at $z/z_0\approx 28$, it has drifted to $t\approx -300\,\mathrm{fs}$. In the non-oscillating case, the observed multi-frequency spectral bands would degenerate into single, sharp spectral lines [36, 37, 47]. Specifically, the points in the spectrum labeled C and D in figure 4(d) correspond to the location of Cherenkov radiation, and the points labeled B and E correspond to additional FWM processes in absence of oscillations. Finally, let us note that in the present propagation scenario we do not observed RR in the vicinity of the point labeled E in figure 5(d). As discussed in references [37, 38], this can be attributed to a low efficiency of the corresponding resonant process. In comparison to the earlier two cases studied in section 3.2, we here find that the amount of RR emanating from the soliton molecule is larger and that the amplitude oscillations decay faster. Let us note that lifting the degeneracies between the resonance conditions (14) and (15) through weak perturbations first leads to a split-up of the resonance lines into multiple lines (see appendix B). Only when increasing the strength of the perturbation, complex spectra as in figure 5 are obtained.

As pointed out above, a viable means to lift existing degeneracies is to consider pulses at a group velocity matched pair of frequencies for which $\beta(\Omega_1)\neq \beta(\Omega_2)$. Above, this was achieved by considering a non-symmetric propagation constant. Similarly, this can be achieved by retaining a symmetric propagation constant but shifting the subpulses to a group-velocity matched pair of frequencies different from the specific pair $\Omega_{1,2}$ considered above. In addition, allowing a frequency dependence of the nonlinear coefficient γ in the nonlinear part of equation (1) directly breaks the symmetry of the subpulses that comprise a molecule. Propagation models of this type were considered in previous studies in references [23, 25, 35, 38].

We studied two-color soliton molecules, i.e. pulse compounds consisting of two tightly bound subpulses at vastly distinct frequencies, supported by nonlinear waveguides with two distinct domains of anomalous dispersion.

First, we analyzed the amplitude and width oscillations of such molecules, induced by increasing the molecule order in a systematic way. The order of a molecule is defined similar to the usual order of solitons in the standard NSE. We contrasted the periodic dynamics exhibited by two-color soliton molecules with that of single solitons in the standard NSE, finding that, in the limit of large molecule duration and despite their complex internal structure, their z-oscillation wavelength matches that of a standard NSE soliton of the same order. By elaborating on a well established variational analysis of the NSE we were led to expect that near-fundamental soliton molecules oscillate harmonically, whereas the amplitude oscillations quickly grow nonlinear with increasing molecule order. This simple and intuitive picture was found to be in good qualitative agreement with our numerical results. Recently it was demonstrated that the oscillation period of breathing solitons in a mode-locked fiber laser, i.e. a dissipative system, depends on the pump power with their relationship following that of a devil's staircase [49]. This relationship is clearly different from that discussed in section 3.1 for our conservative system, which could lead to expect a very different dynamics of two-color soliton molecules in dissipative and conservative systems.

Second, we investigated the emission spectrum of soliton molecules, excited by their periodic amplitude and width oscillations. As soon as a molecule oscillates, it sheds RR in the form of multi-frequency spectral bands, reminiscent of the shape of traditional Japanese Kushi combs. The frequency spacing of resonantly excited modes is determined by the z-oscillation wavelength of the molecule. We demonstrated that the location of the resonant frequencies can be reliably predicted by means of resonance conditions resulting from a phase-matching analysis. Specifically, there are two types of resonance conditions which separately determine the generation of Cherenkov radiation and further FWM processes. We showed that if the propagation scenario is entirely symmetric, these resonance conditions are redundant, so that the generated radiation is degenerate. Breaking the symmetry of the propagation scenario, which we here accomplished by perturbing the propagation constant, the degeneracy can be lifted and distinct multi-frequency spectral bands uniquely associated with Cherenkov radiation, and further FWM processes, are obtained.

The occurrence of such multi-frequency radiation, especially in the degenerate form, constitutes a fundamental phenomenon in nonlinear waveguides. The presented study clarifies the mechanism underlying the generation of RR of two-color soliton pulses, observed previously in terms of theoretical studies [23], and recently also in experiments on multi-wavelength mode-locked solitons in laser cavities [32]. Our results provide further insight into the puzzling propagation dynamics and radiative processes in presence of two-frequency soliton molecules. Due to the manifold of two-frequency pulse compounds with differing substructure, their emission spectra manifest in various complex forms. This demonstrates a peculiar nonlinear optics system for investigating emission of RR.

OM, SW, IB, UM, and AD acknowledge financial support from Deutsche Forschungsgemeinschaft (DFG) under Germany's Excellence Strategy within the Cluster of Excellence PhoenixD (Photonics, Optics, and Engineering–Innovation Across Disciplines) (EXC 2122, Project No. 390833453). AY acknowledges financial support from Priority 2030 Academic Leadership Program and Goszadanie No. 2019-1246. IB also acknowledges support from DFG (Project No. BA4156/4-2). UM also acknowledges support from DFG (Project No. MO 850-20/1).

The data that support the findings of this study are available upon request from the authors.

Considering the standard NSE in nondimensional form

$$\begin{equation} i\partial_z A = -\frac{1}{2}\partial_t^2 A - |A|^2 A, \end{equation} \tag{ A.1 }$$

a variational approach based on the trial function

$$\begin{equation} A(z,t)=\eta\, \mathrm{{sech}}\left(\frac{t}{w}\right) \exp\left(i\theta + ib \frac{t^2}{2w}\right), \end{equation} \tag{ A.2 }$$

yields the set of coupled ordinary differential equations

$$\begin{align} \frac{\mathrm{d}}{\mathrm{d}z}(\eta^2 w) &= 0, \end{align} \tag{ A.3 }$$

$$\begin{align} \frac{\mathrm{d}}{\mathrm{d}z} w &= b, \end{align} \tag{ A.4 }$$

$$\begin{align} \frac{\mathrm{d}}{\mathrm{d}z} b &= \frac{4}{\pi^2 w^3} (1-\eta^2 w^2), \end{align} \tag{ A.5 }$$

$$\begin{align} \frac{\mathrm{d}}{\mathrm{d}z} \theta &= \frac{5 \eta^2}{6} - \frac{1}{3 w^3}, \end{align} \tag{ A.6 }$$

describing the z-variation of the amplitude (η), pulse width (w), pulse phase (θ), and, chirp (b) of the trial function (A.2) [46]. Equation (A.3) implies the constant of motion $\eta^2 w = \mathrm{const.}$, linking intensity and width variations of the trial function (A.2). For an approximate analysis of the pulse width variation along z, we consider the z-derivative of equation (A.4), which, together with equation (A.5), defines the nonlinear equation of motion

$$\begin{equation} \frac{\mathrm{d}^2}{\mathrm{d}z^2} w = \frac{4}{\pi^2 w^3} (1-\eta^2 w^2)\equiv f(w). \end{equation} \tag{ A.7 }$$

Assuming a pulse width $w(z) = w_0 + \epsilon(z)$ that performs periodic oscillations with amplitude $\epsilon(z)$ about some preferred value w₀, we may expand equation (A.7) up to first order in to find $\frac{\mathrm{d}^2}{\mathrm{d}z^2} \epsilon \approx f(w_0) + f^{\prime}(w_0)\, \epsilon$ with $f(w_0) = \alpha w_0^{-3} + \beta w_0^{-1}$, and $f^{\prime}(w_0)\equiv \frac{\partial}{\partial w} f(w)|_{w = w_0} = -3 \alpha w_0^{-4} - \beta w_0^{-2}$, with $\alpha = 4 \pi^{-2}$ and $\beta = -\alpha \eta^2$. Imposing the condition $f(w_0) = 0$ fixes $w_0 = \eta^{-1}$, i.e. a reasonable choice for w₀ given the trial function (A.2), and brings equation (A.7) to the form of a harmonic oscillator

$$\begin{eqnarray} &&\frac{\mathrm{d}^2}{\mathrm{d}z^2} \epsilon(z) + K^2(\eta)\, \epsilon(z) \approx 0, \end{eqnarray} \tag{ A.8 }$$

with angular wavenumber $K(\eta) = 2\sqrt{2}\eta^2/\pi$. The wavelength of this harmonic oscillation then reads

$$\begin{equation} \Lambda^\mathrm{{HO}}(\eta) = \frac{2 \pi}{K(\eta)}= \frac{\pi^2}{\sqrt{2} \eta^2}. \end{equation} \tag{ A.9 }$$

We can readily transfer this result to the case of amplitude and width oscillations of two-color soliton molecules. In the limit where the spectral extent of a molecules subpulses is small in comparison to the frequency gap separating them, i.e. under the condition $\pi t_0 \sqrt{6 \beta_2/|\beta_4|} \gg 1$, each subpulse can effectively be modeled in terms of a standard NSE [39]. Then we are led to expect that the oscillation period of each soliton molecule subpulse, and thus the entire soliton molecule, is reasonably well described by the approximate result equation (A.9). We plot $\Lambda^\textrm{{HO}}(\eta)/z_0$, with $z_0 = 1$ for equation (A.1), as dashed line, labeled D, in figure 3. Let us note that expanding equation (A.7) up to second order in yields a correction term $\propto [f^{^{\prime\prime}}(w_0)/2]\epsilon^2$ with $f^{^{\prime\prime}}(w_0)\equiv \frac{\partial^2}{\partial w^2}f(w)|_{w = w_0} = 40\eta^5/\pi^2$, indicating that the anharmonicity of the oscillation increases rapidly with increasing amplitude. This finding is in qualitative agreement with the numerical results in section 3.1.

To complement the strong symmetry breaking perturbations induced by a non-symmetric propagation constant, discussed in section 3.3, we here demonstrate degeneracy-lifted multifrequency radiation resulting from a small difference in the molecules subpulse amplitudes. Specifically, we consider the symmetric propagation constant $\beta(\Omega) = \frac{\beta_2}{2}\Omega^2 + \frac{\beta_4}{24}\Omega^4$ for $\beta_2 = 1\,\mathrm{fs^2\,\mu m^{-1}}$ and $\beta_4 = -1\,\mathrm{fs^4\,\mu m^{-1}}$, and the initial condition

$$\begin{equation} A_0(t) = a_1 \,\mathrm{{sech}}(t/t_1)\, \mathrm{e}^{-i \Omega_1 t} + a_2\, \mathrm{{sech}}(t/t_2)\,\mathrm{e}^{-i\Omega_2 t}, \end{equation} \tag{ B.1 }$$

where $\Omega_{1,2} = \mp \sqrt{6}\,\mathrm{rad\,fs^{-1} }$, $t_0 = 6\,\mathrm{fs}$, and $a_{1,2}^2 = N_{1,2}^2 2 \beta_2/(\gamma \,t_0^2)$, with $N_1 = 1.57$ and $N_2 = 1.66$. Correspondingly, $a_2\gtrapprox a_1$ in equations (16) and (17). A spectrogram of the resulting vibrating soliton molecule at $z/z_0 \approx 80$ [$z_0 = t_0^2/(2\beta_2)$] is shown in figure A1. As evident from the spectrogram range $[\Omega = -2.3 \ldots -1.5\,\mathrm{rad\,fs^{-1}}] \times [t = -900\ldots -300\,\mathrm{fs}]$, highlighted in close-up view A, the degeneracy of the resonance conditions (14) and (15) is lifted. The primary resonances, resulting from solutions to equation (14), are now accompanied by pairs of secondary resonances, resulting from solutions to equation (15) for $n = 1,2$, which appear on either side of a primary resonance. To facilitate comparison to the degenerate scenario, close-up view B shows the corresponding numerical results for $N_1 = N_2 = 1.59$.

Zoom In Zoom Out Reset image size