2.1.

Experimental Setup

Figure 1 illustrates the laser configuration and detection system employed in investigating the nonlinear response of soliton molecules. A typical dispersion-managed passively mode-locked fiber laser is used to generate multipulse states, utilizing a mode-locking mechanism based on nonlinear polarization rotation (NPR). A 980/1550 nm wavelength division multiplexer connects the laser diode and the gain fiber (Liekki Er 110-4/125), efficiently coupling the pump light into the laser cavity. The free-space components consist of three wave plates, an isolator, and a polarizing beam splitter, utilized, respectively, for polarization adjustment, ensuring unidirectional light transmission, and outputting the laser signals. The laser operates in a near-zero, slightly negative dispersion regime, with a net dispersion of approximately $-0.07\mathrm{p}{\mathrm{s}}^2$. Due to the intensity-dependent loss introduced by NPR, adjusting the gain strength and polarization state offers a convenient and flexible method for tuning the laser into the desired mode-locking regime. Soliton molecules demonstrate robust performance and stability in this regime, even under external injection. The external driving signal is injected through an optical coupler. A continuous wave (Santec TSL-550) signal is initially amplified using a custom-built erbium-doped fiber amplifier. The amplified signal is then modulated by an electro-optic modulator (Conquer KG-AM-15-2.5G) and inscribed with a sweeping driving signal using a function generator. This signal is subsequently coupled into the cavity through the coupler.

Fig. 1

Concept and demonstration of resonant excitation of nonlinear dynamics of soliton molecules. (a) Experimental setup of the mode-locked fiber laser: LD, laser diode; WDM, wavelength division multiplexer; EDF, Er-doped fiber; OC, output coupler; Col., collimator; PBS, polarization beam splitter; QWP, quarter-wave plate; HWP, half-wave plate; ISO, isolator. Driving module: CW, continuous wave; EDFA, erbium-doped fiber amplifier; F.G., function generator; and EOM, electro-optic modulator. (b) Mode-locking status monitoring module: PD, photodiode; OSC, oscilloscope; ESA, electrical spectrum analyzer; OSA, optical spectrum analyzer. (c) BOC system used to detect variations in intrapulse separation. (d) The resonant excitation induces forced oscillations in the intrapulse separation, implying that the system periodically traverses different positions on the anharmonic potential curve under external driving. (e) Schematic of resonance frequency shifts induced by Duffing-type nonlinearity. (f) Schematic of harmonic/subharmonic responses and chaotic dynamics. ${\omega }_d$: driving frequency; ${\omega }_r$: eigenfrequency.

The mode-locked pulse train, optical spectrum, and radio-frequency (RF) spectrum of the laser repetition rate are characterized by an oscilloscope (Agilent Infinium), optical spectrum analyzer (Yokogawa AQ63700), and electrical spectrum analyzer (Rigol DSA815), respectively. A flippable mirror is utilized to direct the laser output into the BOC system [see Fig. 1(c)], enabling the detection of intramolecular dynamics. Scanning the movable mirror in one arm of the Michelson interferometer results in an S-shaped curve, with the linear region (marked in green) corresponding to the zero crossing. For details on the principle and adjustment of the BOC system, refer to our paper,⁴⁸ in which the subfemtosecond resolution of the BOC technology was experimentally demonstrated. By implementing precise control techniques and high-resolution real-time observation methods, we have experimentally observed the resonant response of pulse separations within soliton molecules. The susceptibility of this response exhibits a pronounced shift with increasing resonance strength, providing deep insights into the anharmonic potential governing the interactions between optical solitons [Fig. 1(d)]. We achieved an excellent fit of the resonance frequency shift using the Duffing equation, indicating that Duffing-type nonlinearity plays a crucial role in the nonlinear resonant response within the molecule when the soliton molecule is analogized to a driven damped oscillator [Fig. 1(e)]. Furthermore, we successfully trigger a subharmonic response away from the fundamental resonance frequency. Under appropriate driving strengths, the subharmonic and chaotic response can dominate intramolecular dynamics, reflecting the multiple thresholds required for different oscillating modes in soliton interactions [Fig. 1(f)].

2.2.

Experimental Resonant Excitation of Soliton Molecules

Stationary soliton molecules exhibiting constant intrapulse separation are realized by adjusting the pump strength within an appropriate range, specifically $\sim 520\text{to}600\mathrm{mA}$. The optical spectrum corresponding to the mode-locked state at a pump current of 530 mA, shown in Fig. 2(a), reveals pronounced interference fringe patterns indicative of the established bound-soliton structure. The fringe period of $\sim 2.32\mathrm{nm}$ corresponds to an intrapulse separation of about 3.5 ps, with a pulse duration of $\sim 300\mathrm{fs}$. The rapid intramolecular oscillations exceed the mechanical scanning speed of the spectrometer, resulting in an averaging of the fast dynamics in the output spectrum. Nevertheless, the stability of the intrapulse separation can still be inferred from the fringe visibility, which reaches 100%. In contrast, the intrapulse separation varies over time for oscillating or chaotic soliton molecules, significantly reducing the modulation depth of the interference fringes.^31,41 Figure 2(b) illustrates the pulse sequence, where the interval between consecutive pulses is measured to be 22 ns, corresponding to the round-trip time within the laser cavity. This interval aligns with a laser repetition frequency of $\sim 45\mathrm{MHz}$.

Fig. 2

Experimental measurements of the mode-locked state and resonant excitation of soliton molecules. (a) Optical spectrum and (b) pulse sequence. (c) Response of intrapulse separation under external driving detected by the BOC system. Resonant pulse energy for a driving period is also shown as the brown curve. (d) Resonant susceptibility $\chi (f)$ with Fano line shape fitting, where the red curve represents phase evolution. The driving frequency is also presented as the gray curve. (e) The two-dimensional spectra extracted from the short-time Fourier transform of (c) display fundamental response and overtones. (f) RF spectrum of the laser repetition rate under external driving.

Next, we proceed with the resonant excitation of the intrapulse separation in soliton molecules. The stationary state with constant separation depicted in Fig. 2(a) is employed as the initial state. We apply a logarithmic frequency-swept intensity-modulated sinusoidal signal with a central wavelength set within the gain spectrum of the Er:fiber at 1530 nm, and an average power of $\sim 3\mathrm{mW}$. Empirically, when the injected power exceeds 5 mW, the mode-locked state is readily disrupted. The sweep period is set to 50 ms, covering a frequency range from 100 Hz to 10 MHz. Representative responses of the intromolecular separation over several periods are captured in Fig. 2(c), in which the amplitude of the pulse separation response varies with changes in the driving frequency. A clear correlation between the response amplitude and the driving frequency is observed. The response of intramolecular separation peaks at a specific driving frequency, where the maximum amplitude approaches 100 fs. To explore the pulse energy variation during the frequency sweep, a 10 MHz low-pass filter was used to isolate low-frequency energy fluctuations while filtering out the high-frequency pulse sequences. The brown curve in Fig. 2(c) shows that the pulse energy exhibits a consistent response pattern with the intrapulse separation.

The frequency-dependent susceptibility $\chi (f)$ of the change in intrapulse separation can be simply evaluated from the envelope of the BOC data across the entire frequency sweep. The susceptibility function facilitates the fitting of the resonance curve, thereby evaluating the resonance frequency. As depicted in Fig. 2(d), the blue curve shows the Fano fit applied to the experimental data (scatter points).⁴⁹ For further clarity, the driving frequency over a single period is also plotted (gray curve). From this analysis, we estimate the resonance frequency to be $\sim 5.2\mathrm{kHz}$. We observe a slight asymmetry in the susceptibility curve, which can be attributed to the presence of a nonresonant background. The nonresonant background, typically associated with continuum states, has a complex origin that may include instantaneous Kerr nonlinearity or an inhomogeneous noise floor induced by laser gain depletion.⁵⁰ Additionally, the phase extracted from the susceptibility curve reveals a shift of approximately $\pi$ across the entire resonance region (red curve). Figure 2(e) shows the short-time Fourier transform of a single resonance period of the BOC signal depicted in Fig. 2(c). Along with the fundamental resonance mode, the transform reveals prominent overtones, including the second and third harmonics. The observation of these harmonics underscores the nonlinear characteristics of the dynamical system (soliton molecules here), which are also evident in other nonlinear systems, such as optomechanical oscillators.⁵¹ In addition, we assess the mode-locking stability of soliton molecules during the intrapulse resonance. As depicted in Fig. 2(f), while applying the external signal markedly raises the noise baseline of the RF spectrum of the laser repetition rate, the signal-to-noise ratio remains above 80 dB. This indicates that the mode locking is highly stable.

2.3.

Simulation Results

To verify the generality of the intramolecular resonance response and to uncover its deeper physical insights, we conduct numerical simulations using the scalar generalized nonlinear Schrödinger equation. ⁵² This equation incorporates fundamental physical effects within the laser cavity, including dispersion, self-phase modulation, finite gain bandwidth, and saturable gain. The key components of the model include an EDF as the gain fiber, two segments of single-mode fiber, and an output coupler. A fast saturable absorber model is employed to achieve the mode locking. The gain dynamics are controlled via round-trip—varying modulation of the small-signal gain coefficient $g_0=g_{0\mathrm{s}}(1+M_f\sin (2\pi f{R}_{ts}))$, where $R_{ts}$ is the laser round trips and $g_{0\mathrm{s}}$ refers to the small-signal gain coefficient corresponding to the initial stationary soliton molecules. $M_f$ is the modulation depth, and $f$ is the modulating frequency. Given that the simulation model assumes an ideal environment without accounting for the diverse sources of experimental noise, the resonance response can be realized at a low driving strength. For a detailed description of the simulation model for the ultrafast fiber laser, refer to Sec. A in the Supplementary Material.

Figure 3(a) shows the resonant response of the internal dynamics of soliton molecules at a modulation depth of 0.1%. Both the intrapulse separation and pulse energy exhibit clear resonant responses, which are in high agreement with the experimental results. Additionally, we examine the consistency of energy fluctuations for each pulse within the molecule during the resonance process. The inset shows that pulses 1 and 2 always experience synchronized energy variations. This indicates that the resonant response of the intramolecular system is not directly but rather indirectly related to the gain dynamics, reflecting an intrinsic nonlinear dynamical behavior. Perturbations in gain simultaneously affect the energies of both pulses, leading to changes in the underlying interaction potential and ultimately driving the resonance in the intrapulse separation. Figure 3(b) shows the corresponding spectral evolution with the orange-highlighted spectral profile. Away from the resonance frequency, spectral changes are nearly imperceptible. However, the spectrum undergoes clear periodic broadening and narrowing at the resonance region, exhibiting a distinctive breathing pattern. Following the experimental procedure, we apply a short-time Fourier transform to the response of the intrapulse separation, revealing a distinct second-harmonic response in Fig. 3(c). Higher-order harmonics, such as third-harmonic observed in the experiment, do not manifest due to the minimal driving modulation depth of only 0.1%. In line with the experimental results, the susceptibility in Fig. 3(d) also follows a Fano line shape with a phase shift of $\pi$. A prominent resonance peak is evaluated at $\sim 5.8\mathrm{MHz}$. In Fig. 3(e), we stop the frequency sweeping and stabilize the driving at the resonance frequency to monitor the transition of the intrapulse separation from a constant to oscillatory behavior. This transition occurs rapidly, within fewer than 100 round trips, and demonstrates remarkable reproducibility and stability.

Fig. 3

Simulation of resonant excitation of soliton molecules under 0.1% external driving. (a) Resonance response of intrapulse separation, with the pulse energy response shown by the brown curve. The inset illustrates the energy variation of each pulse within the molecule, highlighting a high consistency. (b) Spectral evolution under frequency sweeping, with the white box highlighting an enlarged view at the resonance region. The orange curve represents the single-shot spectra. (c) Two-dimensional spectra extracted from the short-time Fourier transform of intrapulse separation in (a). (d) Resonant susceptibility $\chi (f)$ with Fano line shape fitting, where the red curve represents phase evolution. Inset: the corresponding Nyquist plot in the complex plane. (e) Transient behavior of intrapulse separation transitioning from no driving to stable forced oscillations, with the driving frequency fixed at 5.8 MHz.

As a further note, it is inherently interesting to observe that the resonance frequency is significantly higher than the experimentally observed values. In our previous work, we applied semiclassical noise theory to analyze stationary soliton molecules under the influence of amplified spontaneous emission, focusing on their quantum noise limit.³⁶ Although this was not explicitly highlighted, we indeed observe a characteristic eigenfrequency, which was designated as the cutoff frequency. Our prior experimental and simulation results demonstrated a pronounced correlation between this eigenfrequency and the intrapulse separation, which well supports the divergence observed between current experimental and simulation results: a longer separation, ranging from a few picoseconds to several tens of picoseconds, corresponds to a significantly lower cutoff frequency, which can drop to just a few kilohertz or even below 1 kHz. Conversely, a shorter separation (0.9 ps in the simulation) corresponds to a much higher cutoff frequency reaching several megahertz. Here, our simulation results are not intended for quantitative comparison with the experiment but are designed to qualitatively support the general resonance dynamics within soliton molecule dynamics.

After successfully validating the experimentally observed resonance dynamics of intrapulse separation, we further investigate the evolution of the resonance response under strong driving conditions in the simulation model. This scenario is particularly challenging to achieve in experimental settings due to the detrimental effect of experimental noise, which diminishes the robustness of soliton molecules to resist intense external perturbations. As illustrated in Fig. 4(a), the waveforms of the susceptibilities remain stable across varying driving strengths. The peak amplitude increases with driving strengths, with its linear growth trend depicted in Fig. 4(b). Notably, a significant shift in the resonance frequency is observed in Fig. 4(a), which indicates a strong dependence of the resonance characteristics on the driving strength. Such behavior is interpreted within the framework of nonlinear physics and is associated explicitly with the softening effect induced by Duffing-type nonlinearity.³

Fig. 4

Simulation of resonant excitation of soliton molecules under driving strengths ranging from 0.5% to 1.2%. (a) Resonant susceptibility $\chi (f)$. (b) Trend of maximum response amplitude and its linear fit. (c) Backbone curve with linear and Duffing equation fitting. (d) Fano fittings of the resonant susceptibility. Inset: phase evolution. (e) Nyquist plots in the complex plane.

We extract the resonance frequencies and plot them on the vertical axis against the maximum amplitudes on the horizontal axis to construct the backbone curve. The result highlights a distinct critical amplitude, beyond which the system, driven with an amplitude significantly exceeding typical operational conditions, manifests pronounced nonlinear effects, marking the dominance of high amplitude-driven operations. As depicted in Fig. 4(c), the resonance frequency decreases with increasing response amplitude, exhibiting a behavior that significantly deviates from a linear fitting (dashed black line). We model this by solving the amplitude for a driven Duffing oscillator,

Equation (2) is used to fit the backbone curve of the resonant soliton molecules, achieving a high degree of consistency, as indicated by a fitting $R^2$ value of 0.97 [green curve in Fig. 4(c)]. For detailed descriptions of the driven Duffing oscillator and the fitting of the backbone curve of soliton molecules, refer to Sec. B in the Supplementary Material. Additionally, we apply Fano fitting to all susceptibility curves in Fig. 4(a). Although the asymmetry of the curves becomes more pronounced with increased driving, the fits in Fig. 4(d) remain reasonable, accompanied by uniform phase curves [inset in Fig. 4(d)] and Nyquist plots [Fig. 4(e)] consistent with the weak driving scenario in Fig. 3(d), as well as the experimental results in Fig. 2(d).

Next, we shift our focus to the frequency range far from the fundamental resonance mode. This approach allows for significantly stronger driving without destabilizing the binding structure of molecules due to excessive resonance effects. We set the driving frequency to 10 MHz, approximately twice the eigenfrequency, and increase the driving strength from 5% to 13%. The resulting separation response is shown in Fig. 5(a), with its short-time Fourier transforms presented in Fig. 5(b). Within the driving strength range of 5% to 8%, the response amplitude increases smoothly with increasing driving strength. However, when the driving strength surpasses 8%, a distinct critical behavior is observed, characterized by a rather sharply defined order-of-magnitude increase in the amplitude and the emergence of a subharmonic frequency at 5 MHz. This behavior indicates that the system has crossed a critical subharmonic threshold, where second-order nonlinear interactions become sufficiently strong to initiate and sustain oscillations at half the driving frequency (5 MHz). Below this threshold, the driving force at 10 MHz primarily excites the fundamental mode, with energy transfer to other modes, including the subharmonic, being too weak to overcome critical damping. However, as the driving amplitude surpasses the threshold, the nonlinear coupling between the fundamental and subharmonic modes intensifies, leading to more efficient energy transfer. This results in a bifurcation in the separation response, where the subharmonic mode not only emerges but stabilizes, increasing in amplitude until it can coexist with or even dominate the oscillation.⁵⁴

Fig. 5

Subharmonic response and deterministic route to intramolecular chaos obtained by sweeping driving strength far from the fundamental resonance region. (a) Intrapulse separation response. Inset: an enlarged view of the boxed region shows a distinct amplitude threshold, above which the subharmonic response dominates. The transition from the fundamental to subharmonic modes manifests period-doubling bifurcation characteristics. (b) Two-dimensional spectra extracted from the short-time Fourier transform of intrapulse separation in (a). The dashed lines spanning (a) and (b) indicate the system’s operation across different modes. (c) Chaotic separation response at a fixed driving strength of 11%. (d) Lyapunov exponent analysis of intrapulse separation in (c). Inset: enlarged view of the linear region. (e) Chaotic spectral evolution.

We observe extreme sensitivity at the bifurcation point in the full-record spectrogram, where pronounced transient instabilities emerge and dissipate. These transient processes are frequently accompanied by spontaneous self-modulation and occur during transitions between various oscillation modes.⁵⁵ As we further increase the driving strength to 9%, beyond the subharmonic threshold, an additional frequency component $f_i$ at $\sim 1.1\mathrm{MHz}$ emerges, incommensurate with the driving frequency. This signals the onset of quasi-periodic dynamics, characterized by a two-dimensional torus in phase space.⁵⁶ In nonlinear systems, quasi-periodic dynamics are intrinsically linked to the emergence of chaos, manifested through irregular and unpredictable oscillations, as demonstrated in our observations of soliton molecules. This chaotic regime emerges when the driving strength surpasses $\sim 11\%$, where the intramolecular system exhibits increasingly complex and erratic behavior. As shown in the spectrogram, the distinct oscillation frequency peaks resulting from nonlinear mode coupling disappear, and the noise floor significantly increases, providing qualitative evidence for the presence of chaos.

To provide quantitative evidence of chaos, we perform a Lyapunov exponent analysis on the response signal of intrapulse separation obtained with a fixed driving strength of 11%. Lyapunov exponents quantify the rate of exponential divergence between nearby trajectories in phase space, which is essential for identifying chaotic dynamics. As shown in Fig. 5(d), the average initial logarithmic divergence of neighboring trajectories follows a linear trend, confirming the presence of chaotic dynamics through positive maximal Lyapunov exponents.⁵⁷ The enlarged view of the linear region is shown in the inset. Finally, the spectral evolution corresponding to Fig. 5(c) is presented in Fig. 5(e), revealing pronounced spectral fluctuations accompanied by oscillations in the interference fringes. These oscillations deviate from the regular breathing periodicity characteristic of breather molecules,³⁹ instead exhibiting erratic and unpredictable oscillatory behavior. We note that the emergence of chaos here does not follow the classic, well-defined routes, but rather involves a new modulated subharmonic route characterized by a “subharmonic–quasi-periodic–chaotic” sequence, highlighting an intrinsic complexity analogous to that of low-dimensional nonlinear systems, as observed in dissipative solitons.^58,59