Spatial solitons have fascinated researchers since the early days of nonlinear optics [1,2]. In combination with self-phase modulation, the self-focusing effect offers the possibility for three-dimensional contraction of an optical wave packet and a concomitant intensity increase. While there exist numerous reports on such light bullets [3–5], this intriguing nonlinear mechanism found little application, probably because of the limiting action of a spatial modulation instability [6,7]. Starting from small imperfection in the beam profile, this process induces a rapid small-scale breakup of the beam profile into filaments when the critical power $P_{\mathrm{crit}}$ is exceeded [8,9], thus limiting the obtainable nonlinear interaction length. Given this severe constraint, high-power pulse compression and nonlinear conversion techniques have resorted to hollow capillaries for extended nonlinear interaction length [10–13]. While several other competing techniques [14,15] have been discussed for the compression of pulses with gigawatt peak powers, the hollow fiber is currently the most established compression technique and has found widespread application in attosecond pulse generation and other high-field experiments; see, e.g., Refs. [16,17]. Utilizing the advanced stretched-fiber technique [18,19] recently enabled record-breaking continuous powers above 300 W [20] and peak powers exceeding 1 TW for the first time [21]. Despite the widespread utility of this technique, however, there exist only relatively few analytical approaches [22,23] for modeling the nonlinear broadening processes inside the hollow fiber. Numerical simulations often resorted to a simplified one-dimensional approximation, as full modal expansions [24–26] are numerically cumbersome. In the following, we present a completely analytical approach for determining spatial soliton solutions in nonlinear multimode fiber geometries. Similar spatiotemporal solitons have previously been observed in numerical simulations [27,28]. Moreover, our approach is mathematically similar to the cage soliton solutions of the Haus master equation of mode locking [29,30], and it is also applicable for the thriving field of multimode fiber nonlinear optics [31–36]. Assuming adiabatic pulse shaping, the results of this analysis enable the derivation of universal scaling laws for the design of nonlinear multimode waveguides.

The linear optical properties of a cylindrical hollow dielectric waveguide with radius $a$ were first modeled by Marcatili and Schmeltzer [10]. Assuming linear polarization, one finds hybrid solutions of the wave equation, which were originally designated as ${\mathrm{EH}}_{mn}$ modes. In solid core multimode (SCM) fibers, widely similar solutions are referred to as ${\mathrm{LP}}_{(m-1)n}$ modes. Strictly speaking, both designations are not identical because of the Goos–Hänchen effect [37], which causes the ${\mathrm{LP}}_{(m-1)n}$ modes to extend into the cladding. In contrast, the ${\mathrm{EH}}_{mn}$ modes exhibit a node of the electric field at the dielectric interface. Assuming $m=1$, i.e., azimuthal homogeneity, the radial field profile of the ${\mathrm{EH}}_{1n}$ mode is given by $$E_n(r)\propto J_0\left(\frac{u_{n}r}{a}\right),$$(1)with the zero-order Bessel function $J_0$ and its $n$th zero $u_n$ ($n>0$); see Fig. 1(a). The complex-valued propagation constant of these modes is given by $${\kappa }_n={\beta }_n+i{\alpha }_n=k_0-\frac{u_n^2}{2k_0{a}^2}(1-i\frac{ϵ+1}{k_{0}a\sqrt{ϵ-1}}),$$(2)where $k_0=2\pi /\lambda$ is the wavenumber and $\sqrt{ϵ}=n_{\mathrm{clad}}$ the refractive index of the cladding material. Here ${\beta }_n$ denotes the propagation constant of the $n$th mode. The propagation constant varies with the eigenmode, which causes intermodal dispersion and is attributed to diffractive effects. ${\alpha }_n$ describes mode-dependent losses in a hollow-fiber geometry. Exploiting the trigonometric relationship between the propagation constant and the wavenumber displayed in Figs. 1(b) and 1(c), one finds an approximate parabolic dependence [38,39] of the propagation constant, $${\beta }_n=k_0\cos \left(\frac{u_n}{k_{0}a}\right)\approx {\beta }_1-\frac{{\pi }^2{(n-1)}^2}{2k_0{a}^2}={\beta }_1+\mathcal{B}{(n-1)}^2.$$(3)Here ${\beta }_1$ is the propagation constant of the fundamental mode, and $\mathcal{B}$ describes intermodal dispersion due to diffraction. A negative value of $\mathcal{B}$ is formally equivalent to anomalous group-velocity dispersion in the formation of time-domain solitons. The quality of the approximation is depicted in Fig. 1(d), clearly indicating the close analogy between diffraction and dispersion for the formation of spatial and temporal solitons, respectively. In order to show that small deviations from a perfect parabolic dependence in Eq. (3) play no major role, we used the exact relationship in computing all of the following example cases. The identical parabolic approximation can be made for SCM fibers [40], yet with positive curvature $\mathcal{B}$. Exploiting the general relation, $${\beta }^2+k_{\perp }^2=k_0^2{\sin }^2\theta +k_0^2{\cos }^2\theta =k_0^2,$$(4)depicted in Fig. 1(c) and reinserting discretized values ${\beta }_n$ and ${\theta }_n=u_n/k_{0}a$, one also finds an approximate linear relationship for the transverse wavenumber, $$k_{\perp ,n}=k_0\sin \left(\frac{u_n}{k_{0}a}\right)\approx (n-\frac{1}{4})\frac{\pi }{a}.$$(5)

Assuming propagation of the ${\tilde{E}}_j$ at identical phase velocity, i.e., as a solitonic wave packet, we find a wavenumber offset $\psi$ relative to the fundamental mode ${\beta }_1$. We further renormalize real-valued electric field amplitudes $a_n=|{\tilde{E}}_n/{\tilde{E}}_1|$ to yield $a_1=1$ and redefine an effective nonlinearity $\gamma =\mathrm{\Gamma }{\mathcal{B}}^{-1}{|{\tilde{E}}_1|}^{-2}$. Using these simplifications, we can extract an algebraic discriminant for the resulting spatial soliton, $$\begin{gathered} \psi =\gamma (1+{\eta }_{1234}{a}_2{a}_3{a}_4+{\eta }_{123}{a}_2^2{a}_3+{\eta }_{12}{a}_2^2+{\eta }_{13}{a}_3^2+{\eta }_{14}{a}_4^2) \\ =\gamma [{\eta }_{1234}{a}_3{a}_4+{\eta }_{123}{a}_2{a}_3+{\eta }_{234}{a}_3^2{a}_4+({\eta }_{12}+a_2^2+{\eta }_{23}{a}_3^2+{\eta }_{24}{a}_4^2)a_2]+a_2 \\ =\gamma [{\eta }_{1234}{a}_2{a}_4+{\eta }_{123}{a}_2^2+{\eta }_{234}{a}_2{a}_3{a}_4+({\eta }_{13}+{\eta }_{23}{a}_2^2+a_3^2+{\eta }_{34}{a}_4^2)a_3]+4a_3 \\ =\gamma [{\eta }_{1234}{a}_2{a}_3+{\eta }_{234}{a}_2{a}_3^2+({\eta }_{14}+{\eta }_{24}{a}_2^2+{\eta }_{34}{a}_3^2+a_4^2)a_4]+9a_4. \end{gathered}$$(10)One can now retrieve all real-valued roots of Eq. (10) and compute beam diameter $w_{\mathrm{eff}}$, loss ${\alpha }_{\mathrm{eff}}$, and, most importantly, the resulting “soliton phase” $\psi$ [47], which is more correctly defined as the total propagation constant involving both linear and nonlinear effects relative to the phase velocity frame. $\psi$ vanishes for $\gamma \to 0$, as we already accounted for linear propagation effects by subtracting the propagation constant ${\beta }_1$ in Eq. (10). This effective subtraction of the linear propagation phase now gives immediate access to the nonlinear phase via $\partial {\phi }_{\mathrm{nl}}/\partial z=\psi$. Using the relation ${\phi }_{\mathrm{nl}}=k_{0}z{n}_{2}P{A}_{\mathrm{eff}}^{-1}$ with the nonlinear refractive index $n_2$ and the effective area [41] of the fundamental ${\mathrm{EH}}_{11}$ mode, $$A_{\mathrm{eff}}=\frac{1}{2\pi }\frac{{\left(\int |E(r{)|}^{2}r\mathrm{d}r\right)}^2}{\int |E(r{)|}^{4}r\mathrm{d}r}\approx 1.50a^2,$$(11)one can then relate the soliton phase to the power $P$ in the nonlinear waveguide. Inserting the definition of the critical power for self-focusing $P_{\mathrm{cr}}=0.147{\lambda }^2/n_2$ [8,9,23], we then finally yield $$P\approx 10.2P_{\mathrm{cr}}\frac{\psi a^2}{k_0{\lambda }^2}.$$(12)

The resulting peak power versus effective nonlinearity is depicted in Fig. 2(a). In the anomalous regime (hollow fibers), the higher-order modes ${\mathrm{EH}}_{1n}$, $n>1$ of the spatial soliton solutions carry the opposite sign of the fundamental ${\mathrm{EH}}_{11}$ mode. This can be understood by a partial cancellation of nonlinearity and modal dispersion, similar to the formation of Schrödinger solitons [2]. In turn, a depression of the peak intensity and an increased effective area $A_{\mathrm{eff}}$ of the beam profile result [Fig. 2(b)]. Consequently, the hollow waveguide can host beams with a peak power of $\approx 1.4P_{\mathrm{cr}}$, yet at a characteristic donut profile; cf. Fig. 3. In the normal dispersion regime (SCM fibers), peak intensities are enhanced as nonlinearity and modal dispersion add up, leading to cusp-like beam profiles for $\kappa >0$. Comparing both waveguide dispersion regimes, it is striking that the asymmetry in the $P(\gamma )$ relation is only caused by beam profile variations. For the SCM case, this effect leads to a maximum peak power hosting of 0.5$P_{\mathrm{cr}}$. Apart from these two fundamental soliton branches, our investigation identifies a second solution branch in the anomalous dispersion regime (red curves in Fig. 2). In this branch, the limit-value donut solution converges toward a “Mexican hat” solution upon subsequent reduction of power. Because of this positive feedback on the losses, the increasing amount of energy in the lossy higher-order modes leads to a runaway, i.e., this solution branch is unstable. This situation can therefore be understood as a blowup of the beam profile in $k_{\perp }$ space, i.e., the energy is quickly transferred into increasingly lossy higher-order modes. In contrast, the fundamental soliton branches (blue and green in Fig. 2) self-stabilize upon propagation as their higher-order contents reduce with decreasing peak power, i.e., there is a negative feedback mechanism.

To this end, it appears illustrative to compute the relevant interaction lengths [41] inside a hollow fiber. Given the rather low group-velocity dispersion of noble gases like argon, which may be additionally cancelled out by the waveguide (group-velocity) dispersion of the hollow fiber, dispersion lengths $L_D\approx {\tau }^2/|{\beta }_2|$ for pulse durations ${\tau }_0>20\mathrm{fs}$ generally exceed the fiber length by a large factor for experimental conditions for typical hollow-fiber compressors in Refs. [19–21]. The absorption length is commonly also chosen longer than the actual fiber length. This is contrasted by the nonlinear length $L_{\mathrm{NL}}={\mathrm{\Gamma }}^{-1}{P}^{-1}$, which amounts to only a few millimeters for powers $P>0.1P_{\mathrm{cr}}$ at identical experimental conditions. One can therefore conclude that spatial soliton effects strongly dominate the nonlinear dynamics inside the hollow fiber, causing an adiabatic reshaping of the solitons as a reaction to the comparatively slow waveguide losses. This adiabatic scenario is contrasted by the use of hollow fibers for high-harmonic generation [48], where much shorter hollow fibers of smaller diameter are typically employed. Apart from the significantly increased modal dispersion, the already compressed input duration also requires consideration of shock terms [42], which we have left out in our analysis. In case of adiabatic reshaping, however, nonsolitonic contents are stripped off into linearly propagating higher-order modes, which travel at reduced group velocities. For example, for $a=200\mathrm{\mu m}$ and $\lambda =1\mathrm{\mu m}$, the linear group velocity difference of the ${\mathrm{EH}}_{12}$ relative to the fundamental mode amounts to about 100 fs/m, i.e., the nonsolitonic contents will lead to the formation of a temporal continuum background after recompression. Ignoring reshaping effects due to group-velocity dispersion, one can now compute the structure of the emerging spatiotemporal light bullets; see Fig. 3. Here we have chosen the highest possible peak power in Fig. 2, i.e., $P=1.4P_{\mathrm{cr}}$, which leads to the formation of a donut spatial structure at the pulse center. At lower intensities, this structure goes over into the more common ellipsoidal shape of conventional light bullets. It needs to be emphasized that temporal broadening effects will accelerate the decrease of peak power. Yet, this effect of group-velocity dispersion can be easily compensated for by chirped mirrors. Uncompensated modal dispersion effects, however, cannot be fixed after nonlinear propagation in the fiber anymore. The severity of uncompensated intermodal dispersion becomes clear from reinspecting Fig. 1(d). An intermodal dispersion of 8 rad/cm between fundamental and ${\mathrm{EH}}_{14}$ modes results in a group delay of 140 fs in a 1 m long linear fiber, which needs to be compared with the typical sub-5 fs pulse duration of the compressed pulses observed after dispersion compensation.

Utilizing the spatial soliton solutions of Eq. (8), one can now derive a few design rules for hollow fiber compressors. As nonlinearly broadened spectra typically exhibit near-perfect spectral symmetry, one can employ the simple relation [41] derived for unchirped Gaussian input pulses, $$\chi (z)=\frac{\mathrm{\Delta }{\nu }_{\mathrm{rms}}}{\mathrm{\Delta }{\nu }_0}=\sqrt{1+\frac{4}{3\sqrt{3}}{\phi }_{\mathrm{nl}}^2(z)},$$(13)to estimate the compressibility of the input pulse. Here $\mathrm{\Delta }{\nu }_0$ and $\mathrm{\Delta }{\nu }_{\mathrm{rms}}$ are the input and output root mean square width of the pulse, respectively. In view of applications, increase of peak power is typically considered more important than ultimate shortness of the pulse. Let us therefore define the figure of merit, $$M(a,\lambda )={\max [\chi (z)\exp (-2{\alpha }_{\mathrm{tot}}z)]|}_z,$$(14)as the criterion for the maximum beneficial propagation length $z_{\max }$ inside the hollow fiber. This limit is reached when detrimental losses start to outweigh the higher compressibility that would result from additional spectral broadening. Consequently, neglecting any possible loss in the subsequent compression process (e.g., by chirped mirrors), peak powers will not further increase upon additional propagation. Pertinent computations are shown in Fig. 4, indicating the dependence of the total (linear and nonlinear) loss ${\alpha }_{\mathrm{tot}}$ as a function of $a^3/{\lambda }^2$; see Fig. 4(a). Compared with experimental data [blue symbols in Fig. 4(a)], it appears striking that early work with relatively short fibers [19] reported losses that deviate from predictions of our model, whereas more recent reports with longer fibers [20,21] appear to completely agree with our model. Moreover, we also compare the maximum beneficial length and the observed spectral broadening between theory and previous experimental findings in Fig. 4(b). Here it is not overly surprising that our projections are too optimistic, with maximum compressibilities and beneficial fiber lengths that are about two-thirds of our predictions. Nevertheless, the relations derived from our spatial soliton model clearly explain the trends observed in previous experiments, confirming that long hollow fibers promise superior performance compared to the single-meter-long segments of early experimentation.

In conclusion, the field of multimode nonlinear optics bears a number of appealing applications, which are mostly ruled by an $a^3/{\lambda }^2$ relationship. In particular, in hollow fibers, where group-velocity dispersion plays an inferior role, spatial soliton formation appears to take a previously unrecognized lead role at large diameters and lengths. Given that dispersive and absorptive lengths are orders of magnitude larger than the soliton length, adiabatic reshaping dominates the nonlinear dynamics of pulse propagation through the hollow fiber. As the critical power $P_{\mathrm{cr}}$ plays a decisive limiting role for the performance of a nonlinear multimode waveguide, further upscaling of compressible peak powers requires usage of lower pressures or gases with lower refractive index than the commonly used argon. As the currently demonstrated highest peak powers already used fibers with several-meter lengths, hosting even higher powers requires the use of significantly longer fibers of tens or even hundreds of meters’ length to accumulate sufficient nonlinear phase for the broadening process. Such dimensions appear to be out of range for universities but could certainly be implemented in large-scale facilities, in particular linear accelerators. Using SCM fibers near their zero-dispersion wavelength instead, spectral broadening can be accomplished at much higher pulse energies than in single-mode fibers. In contrast to previous demonstrations of nonlinear multimode optics, our theoretical investigations suggest that the exact refractive index profile plays only a minor rule, enabling the use of simple step-index architectures rather than relying only on parabolic profiles.

Acknowledgment. G. S. gratefully acknowledges fruitful discussions with Pavel Sidorenko and Frank Wise (Cornell University), Boris Malomed (Tel Aviv University), as well as with Howard Milchberg (UMD).