Spatiotemporal solitons in dispersion-managed multimode fibers
Thawatchai Mayteevarunyoo, Boris A. Malomed, Dmitry V. Skryabin
SSpatiotemporal solitons in dispersion-managed multimode fibers
Thawatchai Mayteevarunyoo ∗ , Boris A. Malomed , , and Dmitry V. Skryabin Department of Electrical and Computer Engineering,Faculty of Engineering Naresuan University, Phitsanulok 65000, Thailand School of Electrical Engineering, Faculty of Engineering, Tel Aviv University, Tel Aviv 69978, Israel Instituto de Alta Investigaci´on, Universidad de Tarapac´a, Casilla 7D, Arica, Chile and Department of Physics, University of Bath, Bath, BA2 7AY, UK
We develop the scheme of dispersion management (DM) for three-dimensional (3D) solitons in amultimode optical fiber. It is modeled by the parabolic confining potential acting in the transverseplane in combination with the cubic self-focusing. The DM map is adopted in the form of alternatingsegments with anomalous and normal group-velocity dispersion. Previously, temporal DM solitonswere studied in detail in single-mode fibers, and some solutions for 2D spatiotemporal “ light bullets”,stabilized by DM, were found in the model of a planar waveguide. By means of numerical methods,we demonstrate that stability of the 3D spatiotemporal solitons is determined by the usual DM-strength parameter, S : they are quasi-stable at S < S ≈ .
93, and completely stable at
S > S .Stable vortex solitons are constructed too. We also consider collisions between the 3D solitons, inboth axial and transverse directions. The interactions are quasi-elastic, including periodic collisionsbetween solitons which perform shuttle motion in the transverse plane. I. INTRODUCTION
Multidimensional solitons represent a vast research area comprising optics, Bose-Einstein condensates (BECs) inultracold gases, plasmas, liquid crystals, and other areas [1–11]. A fundamental problem is that the ubiquitous self-focusing cubic nonlinearity, represented by the Kerr term in optics [12] or the collisional one in the Gross-Pitaevskiiequation for self-attractive BEC [13], creates two- and three-dimensional (2D and 3D) solitons which are unstablebecause the same cubic terms give rise to the critical and supercritical collapse in 2D and 3D, respectively [14]. Onepossibility for the stabilization of 3D matter-wave solitons in BEC, including ones with embedded vorticity [15], isthe use of the trapping parabolic, alias harmonic-oscillator (HO), potential [16], which is, in any case, a necessaryingredient in the experimental realization of BEC [13]. A qualitatively similar mechanism helps to stabilize 3D “optical bullets” [17] (spatiotemporal solitons) in multimode optical fibers, i.e., ones defined by the radial patternof the graded (refractive) index (GRIN) in the transverse cross-section plane. Such fibers have been a subject offundamental and applied research since long ago [18–21], due to their potential for the use in optical sensors [22, 23],high-speed interconnects [24, 25], and space-division multiplexing [26–28]. The GRIN structure, supporting manytransverse modes (see, e.g., Ref. [29]), makes it possible to consider 3D solitons as nonlinear superpositions of suchmodes, self-trapped in the temporal dimension, i.e., along the fiber’s axis [30–34]. Recently, this approach to thestudy of spatiotemporal solitons has drawn much interest [35–44].Another method for stabilization of both one- and multidimensional solitons, in the form of oscillating breathers,with an intrinsic phase chirp [12], is provided by ”management” techniques. They are represented by periodicmodulation of basic parameters of the medium between positive and negative values, along the propagation distance,in terms of optics, or in time, in terms of BEC [45]. A well-known example is the dispersion management (DM),i.e., transmission of temporal solitons through a composite line built as a concatenation of single-mode fibers withanomalous and normal group-velocity dispersions (GVD) [46]. DM provides the remarkable stabilization of the solitonsin communication lines against various perturbations, such as the Gordon-Haus jitter (induced by the interactionof solitons with random optical radiation) [47]. Although the DM soliton runs through the chain of segments withopposite values of the GVD coefficient, which drives strong intrinsic oscillations in it, extremely accurate simulations ofthe respective nonlinear Schr¨odinger equation (NLSE) have demonstrated that the ensuing oscillations of the soliton’sshape do not destabilize it, even if the propagation extends over thousands of DM periods [48–51]. Furthermore, it waspredicted that 2D spatiotemporal “ bullets” in a planar waveguide, composed of alternating segments with oppositesigns of GVD, also propagate in the form of robust breathers with strong intrinsic oscillations [52, 53]. In particular,stable 2D breathers may feature periodically recurring fission in two fragments and recombination into a single soliton[53]. In the bulk waveguide, 3D dispersion-managed “ bullets” are unstable, but they may be stabilized by inclusionof the defocusing quintic nonlinearity, which accounts for saturation of the cubic self-focusing nonlinearity [54]. Inaddition to the stabilization of solitons, the DM technique finds other applications to nonlinear optical media, suchas enhancement of supercontinuum generation [55].While DM applies to optical media, the technique of “ nonlinearity management”, i.e., periodic alternation of self-focusing and defocusing, is chiefly relevant to BEC, where it may be implemented by periodically switching the sign a r X i v : . [ phy s i c s . op ti c s ] N ov of the nonlinearity with the help of the Feshbach resonance controlled by an external magnetic field [56, 57], whichis made periodically time-dependent, for that purpose. The analysis, performed in various forms, has predicted veryefficient stabilization of 2D ground-state breathers, while states with embedded vorticity and all 3D solitons remainunstable under the action of the nonlinearity management [58–64]. Matter-wave solitons may be made stable in 3D ifthe time-periodic nonlinearity management is combined with a quasi-1D spatially periodic potential (optical lattice)[65].The fact that DM helps to create and strongly stabilize oscillating solitons in composite single-mode fibers [45,46, 50], and 2D oscillatory “ bullets” in composite planar waveguides [52, 53], suggests to consider a possibility ofthe creation of robust spatiotemporal solitons in dispersion-managed multimode waveguides, composed of alternatingpieces of GRIN fibers with opposite signs of GVD. In particular, the trend of the management to suppress instabilityagainst the collapse [45] offers a possibility to create stable high-power solitons, that may be useful for applications.Although splicing of segments of multimode fibers in the composite system is a technological challenge, it has beenimplemented in various forms, [66–70], including large-transverse-area fibers with the same parabolic profile of thelocal refractive index as considered in the present work [71].Alternative options are to replace one fiber species by adispersive grating [72], or control the effective GVD by means of off-axis light propagation [73].The objective of the present work is to identify stable and quasi-stable spatiotemporal solitons in the DM multimodesystem, with the GRIN structure represented by the HO trapping potential. We also address solitons with embeddedvorticity, as well as collisions between solitons.The evolution of the complex amplitude A ( X, Y, Z, T ) of the electromagnetic field in the multimode fiber is governedby NLSE with propagation distance Z , transverse coordinates ( X, Y ), and reduced time T in the coordinate systemtraveling at the group velocity of the carrier wave [30, 32, 35, 74, 75]: i ∂A∂Z = − k (cid:18) ∂ A∂X + ∂ A∂Y (cid:19) + 12 β (cid:48)(cid:48) ( Z ) ∂ A∂T + k ∆ R (cid:0) X + Y (cid:1) A − γ | A | A, (1)where k is the propagation constant of the carrier wave, and β (cid:48)(cid:48) ( Z ) is the GVD coefficient which takes oppositevalues in alternating segments of the multimode fibers. Further, the coefficient in front of the HO potential, ∆ ≡ (cid:0) n − n (cid:1) / (cid:0) n (cid:1) , is the relative difference of the refractive index, n , between the fibers’ core and cladding, R is the core’s radius, and the nonlinearity coefficient is γ ≡ k n / ( n core A eff ), where n is the Kerr coefficient, andA eff the effective area of the fiber’s cross section. The present model disregards the difference in γ between differentfiber segments, as it is known that the DM is a much stronger factor than the difference between different values ofthe self-focusing coefficient [45].Equation (1) is cast in a dimensionless form by the substitution: x = X/w , y = Y /w , z = Z/ (cid:0) k w (cid:1) ,τ = T /T , u ( x, y, τ ; z ) = (cid:112) γk w A ( X, Y, Z, T ) , (2)where T and w are the temporal and transverse scales (the longitudinal one being Z = k w ). This leads to thenormalized NLSE, i ∂u∂z = − (cid:18) ∂ ∂x + ∂ ∂y (cid:19) u − D ( z ) ∂ u∂τ +( x + y ) u − | u | u, (3)where the transverse scale is chosen as w = (cid:0) R /k ∆ (cid:1) / , to make the coefficient in front of the HO potential equalto 1, and D ( z ) ≡ − (cid:16) R/T √ ∆ (cid:17) β (cid:48)(cid:48) ( Z ).The DM map , i.e., the scheme of the periodic alternation of the GVD coefficient in Eq. (3), is defined as follows: D ( z ) = D anom ,D norm ,D anom , < z < z anom , z anom < z < z anom + z norm , z anom + z norm < z < z map = 0 . . (4)Here, D anom = ∆ D + D av > D norm = − ∆ D + D av < D av is the average GVD value, and ∆ D (cid:29) D av is the DM amplitude. The size of the DM period is fixedin Eq. (4) to be z map ≡ z anom + z norm = 0 . z anom = z norm = 0 . u ( x, y, τ ; z ) into a truncated superposition of commonly known eigenmodes of the isotropic HO potential in the( x, y ) plane, while expansion amplitudes are considered as functions of τ and z , governed by an approximate systemof coupled 1D equations. However, we prefer to develop a “ holistic” approach, relying upon numerical simulations ofthe fully three-dimensional NLSE (3). While the number of actually excited transverse modes is effectively finite forany 3D solution, the advantage of using the full 3D equation is that this number is not restricted beforehand. Theapplicability of the 3D NLSE for the description of the multimode propagation in fibers with an internal transversestructure has been demonstrated, in other contexts, both theopretically and experimentally – see, e.g., recent work[76] and references therein.As concerns physically relevant scales, the characteristic propagation length for DM schemes in singe-mode fibersis measured in many kilometers, while for the multimode waveguides it is limited to a few meters [18]-[28], [35]-[44].This fact makes it possible to neglect losses in Eq. (1). On the other hand, it may be relevant to add, to Eq. (3),higher-order terms, representing, in particular, the third-order GVD and the intra-pulse stimulated Raman scattering.In this work, we focus on the model based on the basic NLSE (1), as it was shown previously that the additionalterms, although affecting the shape of DM solitons, do not lead to dramatic changes in their dynamics [77–79].Basic results for the existence and stability of 3D spatiotemporal solitons (including ones with intrinsic vorticity),under the action of DM, which are produced by numerical analysis of the model based on Eqs. (3) and (4), arereported in Section II. Collisions between 3D solitons in the axial and transverse directions are addressed in SectionIII, in the absence and presence of the DM. The paper is concluded by Section IV. II. BASIC NUMERICAL RESULTS: STABLE SPATIOTEMPORAL SOLITONS UNDER THE ACTIONOF DMA. Families of spatiotemporal solitons in the absence of DM
First, we produce a family of 3D solitons in the absence of the DM, i.e., setting ∆ D = 0 in Eq. (4). Stationarysoliton solutions to Eq. (3), with real propagation constant k , are looked for as u ( x, y, τ ; z ) = e ikz φ ( x, y, τ ) , (5)where the real function φ ( x, y, τ ) is a localized solution of equation kφ − (cid:18) ∂ φ∂x + ∂ φ∂y (cid:19) − D av ∂ φ∂τ +( x + y ) φ − φ = 0 . (6)This equation was solved by means of the Newton conjugate-gradient method [80]. The so obtained soliton solutionsare characterized by the total energy, E = (cid:90) (cid:90) (cid:90) | u ( x, y, τ ) | dxdydτ, (7)and the temporal and spatial FWHM widths, W τ, FWHM and W s, FWHM , which are extracted from the numerical dataat x = y = 0 and τ = 0, respectively. Then, stability of the solitons was investigated by means of the linearizationof Eq. (3) for small perturbations, added to the stationary soliton solutions, and computation of the respectiveeigenvalues. This was done by means of a numerical method borrowed from Ref. [30]. The so predicted stability wasthen verified by direct simulations of Eq. (3).Note that, in the linear limit and in the absence of GVD, D av = 0, eigenvalue of k in Eq. (6) corresponds tothe commonly known ground-state energy of the 2D HO potential: k = −√
2. As seen in Eq. (6), the contributionfrom the anomalous GVD ( D av >
0) of temporally localized states shifts k towards more negative values, while thecontribution to k produced by the self-focusing cubic term is positive, therefore spatiotemporal solitons may exist atboth k < k > k = 0), see Fig. 1(a).The results for the shape and stability of the spatiotemporal solitons are summarized in Fig. 1, by means of theenergy and width curves, E ( k ) and W s,τ, FWHM ( E ), for the soliton family. In particular, the E ( k ) dependence in Fig.1(a) and the stability boundary demonstrated by it are close to those produced in Ref. [16] for a strongly elongated3D HO trapping potential, which is close to the 2D potential in Eq. (3) (in the same work, it was demonstrated that (a) (b) FIG. 1: (a) The integral energy (7) for the family of stationary spatiotemporal solitons, obtained in the absence of the DM[∆ D = 0, D av = 1 in Eq. (4)], versus the propagation constant, k . The critical point, with dE/dk = 0, is located at k = − . W τ, FWHM , and spatial, W s, FWHM , FWHM widths of the stationary solitons (solid and dashed lines, respectively)versus E . The colors have the same meaning as in (a). (a)(b) FIG. 2: (a) The evolution of the stable spatiotemporal soliton with k = − . E = 6 . D = 0, D av = 1 in Eq. (4)], corresponding to the pink dot in Fig. 2. The evolution is shown by means of the set ofisosurfaces of local intensity, | u ( x, y, τ ) | = 0 .
8, and its temporal profiles. (b) The evolution of the temporal (left) and spatial(right) local intensity. the E ( k ) curve can be accurately predicted by means of the variational approximation). In particular, there are twodifferent solutions for a given energy, their stability obeying the necessary condition, dE/dk >
0, i.e., the celebratedVakhitov-Kolokolov criterion [81, 82]. The stability, as it is shown in Fig. 1, was identified through the computationof eigenvalues for small perturbations. In addition to stable and unstable stationary 3D solitons, a narrow interval ofrobust breathers, replacing unstable solitons, is found close to the critical point, dE/dk = 0 [the short green segmentin Fig. 1(a)]. The coexistence of stable and unstable families of 3D solitons, demonstrated by Fig. 1, is a genericfeature of the 3D NLSE with the cubic self-focusing and HO trapping potential [16].Results for the stability, displayed in Fig. 1, were verified by direct simulations of Eq. (3) for the evolutionof the spatiotemporal solitons, performed by means of the split-step fast-Fourier-transform algorithm. An example,presented in Fig. 2 for the soliton with propagation constant k = − E = 6 . ∼
200 dispersion/diffractionlengths of the soliton. In the course of propagation, the soliton energy is conserved with relative accuracy ∼ − ,demonstrating that the soliton is an absolutely robust solution of the underlying DM model. In the absence of DM,the phase chirp of stable solitons remains equal to zero, in the course of the simulations. B. Dispersion-managed spatiotemporal solitons
The main parameter which controls the action of DM in the temporal domain is the
DM strength [45, 46, 50], S ≡ ( D anom z anom + | D norm | z norm ) / (cid:0) W τ, FWHM (cid:1) min (8)where subscript min refers to the smallest value of the temporal width of the periodically oscillating DM soliton.To identify the temporal and spatial chirps of the soliton oscillating under the action of the DM, C τ and C s , itis represented in the Madelung’s form, as u ( x, y, τ, z ) ≡ | u | exp ( iχ ). Then, the chirps were computed from thenumerically identified phase, χ ( x, y, τ, z ), as C τ = ∂ ∂τ χ ( x = 0 , y = 0 , τ, z ) | τ =0 ,C s = ∂ ∂x χ ( x, y = 0 , τ = 0 , z ) | x =0 , (9)The simulations of Eq. (3) with the DM map (4) were initiated with an input taken as stable stationary solitonsnumerically produced in the system without the DM (see above). It was found that the outcome of long-distancesimulations is adequately determined by the value of the DM strength (8). First, for S < S ≈ .
93 (10)(relatively weak DM), the simulations produce quasi-stable 3D solitons, as shown in Fig. 3 for S ≈ .
62. While thesoliton keeps its overall integrity and does not decay in the course of evolution (Fig. 3(d) demonstrates that, havingpassed 300 DM periods, the soliton has lost <
1% of the initial energy, through emission of small-amplitude radiation),it develops quasi-random oscillations of its characteristic parameters, although with a relatively small amplitude. Thisis seen in Fig. 3(c), that demonstrates a random walk of the soliton’s trajectory, which is trapped in a small domainof the plane of relevant dynamical parameters, viz ., the temporal width and chirp.It is relevant to mention that the spatiotemporal soliton displayed in Fig. 3(b) has the dispersion and diffractionlengths < ∼
1. This is comparable to the underlying DM period, z map = 0 .
5, which corroborates that the DM is anessential ingredient of the system under consideration. The same pertains to the completely stable spatiotemporalsoliton displayed below in Fig. 4(b).The boundary value (10) of the region of quasi-stable DM solitons corresponds, e.g., to the DM map with ∆ D ≈ D av = 1), and temporal width ( W τ, FWHM ) min ≈ .
84. At
S > S , the propagation of the spatiotemporalsolitons is completely stable under the action of moderate or strong DM. A typical example is displayed in Fig. 4for S ≈ .
13. The most essential manifestation of the full stability is that, in Fig. 4(d), the input loses (cid:39) .
5% ofits total energy at the initial stage of the evolution ( z <
50 DM periods), adjusting itself to the propagating state,and then, at z >
50, the emission of radiation completely ceases. The fully regular dynamics of the soliton is alsodemonstrated by the evolution of its spatial and temporal parameters displayed in Fig. 4(a), cf. Fig. 3(a) for thequasi-stable spatiotemporal soliton.Numerically exact DM solitons can be produced by means of the averaging method, which was previously elaboratedfor temporal solitons in dispersion-managed single-mode fibers [49, 50]. The method is based on collecting a set ofshapes of an oscillating soliton, produced by the straightforward simulations, at points where the shape is narrowest,and computing an average of the set. The result is the pulse which propagates in a strictly periodic form. A conclusionof further numerical analysis is that the direct simulations displayed in Fig. 4, as well as in other cases of the completelystable propagation, converge precisely to the DM solitons produced by the averaging method, extended to the present3D setting. It is relevant to stress that the averaging procedure needs to be applied only in the temporal direction,while in the plane of ( x, y ) the solution readily converges by itself to the one predicted by the averaging method. Asan example, the DM soliton, to which the evolution displayed in Fig. 4 converges, is shown in Fig. 5. In Fig. 5(a),the isosurface plot at the top displays the spatiotemporal evolution of the DM soliton while it passes three DM maps,0 < z < .
5, while other panels display the variation of the soliton’s temporal and spatial characteristics. In Fig.5(b), the isosurface plot corroborates the full stability of the soliton over the propagation distance equivalent to 100maps.
C. Three-dimensional solitons with embedded vorticity
The creation of stable spatiotemporal DM solitons with embedded vorticity is a challenging objective. To the bestof our knowledge, 3D solitons of such a type have not been reported before. To this end, it is necessary, first, toconstruct self-trapped vortex states as solutions to Eq. (3) in the absence of DM (∆ D = 0). In terms of the polarcoordinates, ( r, θ ) in the plane of ( x, y ), stationary vortex solitons are sought for as u = exp ( imθ ) φ ( r, τ ) , (11)with integer vorticity m ≥ φ obeying equation kφ − (cid:18) ∂ ∂r + 1 r ∂∂r − m r − D av ∂ ∂τ (cid:19) φ +( x + y ) φ − φ = 0 , (12) (a)(b)(c) (d) FIG. 3: (a) The left column shows the evolution of the peak intensity, temporal width and chirp of a quasi-stable spatiotemporalsoliton under the action of weak DM with strength S = 0 .
62, corresponding to parameters ∆ D = 5 and D av = 1 in Eq. (4),over 300 map periods (which corresponds to z = 150; here and in other figures, values z map refer to the propagation distancemeasured in units of the DM map). Black lines represent values of the same variables, taken at the beginning of each DMmap (4). The input is the same stationary soliton which is shown in Fig. 2 in the absence of DM, with temporal width W τ, FWHM ≈ .
65. The right column displays the evolution of the soliton’s spatial width and chirp. Magenta lines representtheir values at the beginning of each DM map. (b) Profiles of the oscillatory solitons in the temporal (left) and spatial (right)cross sections, plotted at the beginning of each DM map. (c) The trajectory of the soliton in the plane of the temporal width( W τ, FWHM ) and chirp ( C τ ). Black and blue markers denote the initial and final points of the trajectory. (d) The total energy(7) of the soliton versus the propagation distance. with boundary condition φ ∼ r m at r →
0. In the linear limit and for D av = 0, k is determined by energy eigenvaluesof the 2D HO potential, i.e., k linear = − (1 + m ) √ m = 1, similar to that reported in Ref. [16]for the NLSE with a strongly anisotropic HO trapping potential (i.e., the Gross-Pitaevskii equation), is displayed inFig. 6. As seen in the figure, the Vakhitov-Kolokolov criterion, dE/dk >
0, is necessary but not sufficient for thestability of the vortex solitons. An example of the stable vortex, shown by means of the isosurface of | u ( x, y, τ ) | ,clearly shows the inner hole, maintained in the soliton by the embedded vorticity. Unstable vortices follow the usualscenario, spontaneously splitting in a pair of fragments [15], as shown in Fig. 7. Eventually, the fragments merge intoa quasi-turbulent state.In the presence of the DM, quasi-stable vortex solitons were found by direct simulations. A typical example isshown in Fig. 8: using a stable vortex from Fig. 6, which was found for ∆ D = 0, as the input, direct simulationsdemonstrate the formation of a robust soliton which keeps the intrinsic vorticity and features slow shape oscillationswith a large period, ∆ z (cid:39)
80 DM periods, under the action of the relatively strong DM, with ∆ D = 30. Theshape oscillations may be removed by accurately tuning the input, and a stability area for vortex solitons may bethus identified in the parameter space. Here, we do not aim to report such results in a comprehensive form, as it (a)(b)(c) (d) FIG. 4: The same as in Fig. 3, and with the same input, but for the propagation of a stable spatiotemporal soliton under theaction of moderately strong DM, with strength S ≈ .
13, corresponding to ∆ D = 30 and D av = 1 in Eq. (4). is extremely time-consuming to accumulate the necessary amount of numerical data. We did not consider vortexsolitons with m > III. COLLISIONS BETWEEN THREE-DIMENSIONAL SOLITONSA. Collisions in the longitudinal direction
Once stable DM solitons have been found, it is relevant to consider collisions between them. The DM soliton canbe set in longitudinal motion by the application of a kick (i.e., longitudinal boost) to it, with arbitrary frequency shiftΩ. Indeed, Eq. (3) is invariant with respect to the Galilean transformation, which generates a new solution from agiven one, u ( x, y, τ ; z ): u ( x, y, τ ; z ) = exp (cid:18) − i Ω τ − i (cid:90) D ( z ) dz (cid:19) × u (cid:18) x, y, τ − τ + Ω (cid:90) D ( z ) dz ; z (cid:19) , (13)where τ is an arbitrary constant shift of the solution as a whole. Note that, in addition to the progressive motionwith average speed speed = − D av Ω , (14) (a)(b) FIG. 5: (a) The top panel displays the spatiotemporal dynamics of the numerically exact DM soliton, with strength S ≈ . W τ, FWHM ) min ≈ .
65, to which converges the soliton presented in Fig. 4. The evolution ofthe spatiotemporal DM soliton, produced by means of the averaging method (see the main text), is displayed by means ofan isosurface of the local intensity, | u ( x, τ, z ) | = 0 .
25, in the cross section of y = 0, as the soliton passes three DM maps,0 < z < .
5. The bottom plots show the variation of the soliton’s temporal and spatial characteristics. (b) The isosurface of | u ( x, , τ, z ) | = 0 .
25, shown at the beginning of each DM map, corroborates the full stability of the soliton as it passes 100DM periods, 0 < z < the substitution of DM map (4) in Eq. (13) gives rise to oscillatory motion with spatial period z map = 0 . τ = (1 / D | Ω | z map . (15)Thus, it is possible to create initial conditions for simulating collisions between spatiotemporal DM solitons movingin opposite temporal directions. In studies of temporal DM solitons in single-mode optical fibers, collisions werestudied in detail for solitons carried by different wavelengths in the wavelength-division-multiplexed (WDM) system,which is an important practical problem, as the collision-induced jitter is a source of errors in data-transmissionschemes [83, 84].Numerical simulations demonstrate elastic collisions between boosted solitons, created as per Eq. (13), with fre-quency shifts ±| Ω | and centers initially set at points τ = ± τ , see an example in Fig. 9 for Ω = ± .
3. This “ slow”collision is strongly affected by the DM, because, with the smallest value of the solitons’ width ( W τ, FWHM ) min ≈ .
65 and relative speed 2 D av | Ω | = 0 . Z coll ≈ W τ, FWHM ) min / (2 | Ω | ) (cid:39) z map . This estimateimplies that multiple collisions take place between the two solitons, while they stay strongly overlapping. Indeed. Eq.(15) yields an estimate for the characteristic overlapping degree: ( W τ, FWHM ) min / ∆ τ ≈ .
62 for the present values ofparameters. This circumstance explains a complex intermediate structure observed at z = 16 in Fig. 9(a), and in theleft bottom plot of Fig. 9(b).In addition to the initial temporal separation 2 τ , a pair of colliding solitons is characterized by a phase shiftbetween them, ∆ χ . However, additional numerical results clearly demonstrate that, similar to what is well known inmany other systems, results of the collisions do not depend on ∆ χ . Actually, the collisions are driven by ∆ χ in thecase when the initial pair is taken with zero relative velocity and relatively small separation [85], which is not the casehere. Nevertheless, collisions displayed below in Fig. 11 are sensitive to the fact they start with ∆ χ = 0. FIG. 6: The energy of numerically generated vortex-soliton solutions of Eq. (3) with D av = 1 and m = 1 [see Eq. (11)], in theabsence of DM (∆ D = 0), vs. the propagation constant. Black and red segments denote stable and unstable subfamilies, asidentified through the computation of eigenvalues for small perturbations. The top inset shows the isosurface plot, | u ( x, y, τ ) | =0 .
35, of a typical stable vortex soliton corresponding to the dot on the E ( k ) curve, with E = 15 .
97 and k = − .
3. The bottominset: the stability of this soliton is illustrated by simulations of its perturbed evolution, displayed on the line of y = τ = 0.An example of unstable evolution of the vortex soliton, corresponding to k = 0, is presented in Fig. 7.
25 30 35 40 z = 0 10 15 20 FIG. 7: The isosurface evolution of an unstable vortex soliton with m = 1, simulated in the framework of Eq. (3) with ∆ D = 0and D av = 1. The input is a stationary solution of Eq. (12) for k = 0, which is definitely unstable, as seen in Fig. 6. It is worthy to note that the collision produces, in a short interval of the propagation distance, a peak with alarge amplitude, as seen in the left bottom panel of Fig. 9(b). This picture somewhat resembles the formation ofoptical rogue waves, see, e.g., recent work [86] and references therein. However, the analogy is rather superficial, as,unlike traditional rogue waves, here the peak appears not spontaneously, but as a result of the collision, it is notfed by a continuous-wave background, and the configuration is not intrinsically unstable (as it does not include amodulationally unstable background).
B. Collisions in the transverse direction
In the present context, it is relevant to consider, as a new option, collisions between 3D solitons moving in thetransverse direction in the 3D setting. For this purpose, solitons may be set in motion by initially placing them atoff-axis positions, with nonzero initial coordinates of the soliton’s centers, ( x , y ), and letting them roll down towards x = y = 0 in the HO potential (a similar setup was employed to initiate collisions between matter-wave solitonstrapped in the 3D isotropic HO potential [85]).0 (a) (b) FIG. 8: (a) Quasi-stable evolution of a vortex soliton with m = 1, simulated in the framework of Eq. (3) with relatively strongDM, viz ., ∆ D = 30 and D av = 1. The input is the stationary vortex soliton with ∆ D = 0 and D av = 1 shown in Fig. 6. (b)The total energy of the vortex soliton versus the propagation distance. z = 0 5.0 10.013.0 16.0 18.020.0 25.0 30.0 (a)(b) FIG. 9: Simulations of slow collision, in the axial (temporal) direction, between two DM solitons, constructed as per Eq. (13)from the stable DM mode produced in Fig. 4 [with ∆ D = 30, D av = 1, ( W τ, FWHM ) min ≈ . τ = ±
5, are boosted by frequency shifts Ω = ± .
3. (a) A set of 3D snapshots, taken in the course of the collision,are depicted by means of isosurfaces of the local intensity, | u ( x, y, τ ) | . (b) The collision shown by means of the isosurface inthe cross section of y = 0 (the top plot), and by the evolution of the intensity on the lines of ( x = y = 0) and ( y = τ = 0) (theleft and right bottom plots, respectively). Thus, it is possible to consider the collision between two identical solitons initially shifted to positions ± ( x , y ).It makes sense to address this possibility, first, in the absence of DM (∆ D = 0), as it was not addressed in previousworks. As shown by a typical simulation displayed in Fig. 10, the colliding solitons, originally created (at z = 0)at positions x = y = ±
3, pass through each other quasi-elastically (QE) at z ≈ .
1, separate, then return underthe action of the transverse HO potential, and again feature the QE collision at z ≈ .
3. The extension of thesimulation demonstrates that the chain of QE collisions between the solitons, which perform the shuttle motion inthe HO potential, continues indefinitely long, with intervals∆ z coll ≈ . ω ⊥ = √ z oscill = π/ω ⊥ ≈ .
22, which readily explains the size ofthe interval between the collisions, see Eq. (16). Similar results were produced by simulations of head-on collisionsbetween a soliton released from the initial position with coordinates x = y > x = y = 0). Similar to what is mentioned above for collisions in the longitudinal direction, head-oncollisions in the transverse plane are not sensitive to ∆ χ either.1 z = 0 0.7 1.1 1.52.2 3.3 3.7 4.3 FIG. 10: Collisions, in the transverse plane, between two 3D solitons in the absence of DM, i.e., with ∆ = 0 and D av = 1 in Eqs.(12) and (4). At z = 0, the solitons, identical to the one presented in Fig. 2, are created at off-axis positions x = y = ± | u ( x, y, τ ) | . The quasi-elasticcollisions between the solitons repeat periodically in the course of indefinitely long simulations.FIG. 11: The same as in Fig. 10, but in the case when centers of the two in-phase solitons are initially separated in the axialdirection by shifts τ = ± It is also relevant to consider interactions between two solitons initially created with centers placed at points withcoordinates ( x, y ) = ± ( x , y ), additionally separated in the axial direction, i.e., with temporal coordinates τ = ± τ .In this case, the solitons do not collide head-on; nevertheless, the simulation displayed in Fig. 11 demonstrates thatthe attractive interaction between the in-phase 3D solitons, mediated by their tails [87], leads to the collision betweenthem. As mentioned above, these collisions, unlike those displayed in Figs. 9 and 10, are sensitive to the fact thatthe solitons were created with zero initial phase difference, ∆ χ = 0; in the case of ∆ χ = π , the collision does nottake place, as the solitons interact repulsively (not shown here in detail). In this case, simulations also demonstratea periodic chain of QE collisions with the same interval as predicted by Eq. (16). However, the corresponding fullperiod of the collisional dynamics is double (including two collisions), ∆ z (full)coll ≈ .
4, because each collision leads torotation of the line connecting centers of the solitons, as seen in Fig. 11.In the presence of the DM, collisions in the transverse direction are similar to those shown in Figs. 10 and 11.In particular, Fig. 12 demonstrates simulations of the head-on collision between the DM soliton, released from theposition with x = y = 3, and a quiescent one, placed at x = y = 0. Note that, at the respective values of theparameters, the transverse speed of the soliton rolling down under the action of the HO potential at the collisionmoment is (speed) ⊥ = − ω ⊥ (cid:112) x + y = − √
2. Then, the propagation distance necessary for the completion of thecollision in the transverse direction can be estimated as(∆ z coll ) ⊥ (cid:39) W s, FWHM / | (sp) ⊥ | (cid:39) . , (17)where estimate W s, FWHM (cid:39) . FIG. 12: Collisions, in the transverse plane, between two 3D solitons in the presence of DM, with ∆ D = 30 and D av = 1 inEqs. (12) and (4). At z = 0, a DM soliton, the same as the one presented in Figs. 4 and 5, is placed at the position with x = y = 5, and an identical soliton is placed at x = y = 0. The dynamics of the collision is shown by means of the set ofisosurfaces of the local intensity, | u ( x, y, τ ) | . Quasi-elastic collisions between the solitons recur periodically in the course ofthe evolution.FIG. 13: Recurring collisions between the same DM solitons as in Fig. 12, but initially placed at positions with coordinates x = y = τ = 5 and x = y = τ = 0, i.e., with the additional initial separation ∆ τ = 5 in the axial (temporal) direction. z map = 0 . z (full)coll ≈ .
4, including two collisions, as the line connecting centers of the interacting solitons rotates, asa result of each collision. The comparison of Figs. 11 and (13) suggests that, as well as in the case of the head-oncollision, DM produces a moderate effect on the interaction of the spatiotemporal solitons.3
IV. CONCLUSION
The objective of this work is to extend the concept of DM (dispersion management) for 3D spatiotemporal solitonsin multimode nonlinear optical fibers. In previous works, the DM concept was elaborated in detail, theoretically andexperimentally, for 1D temporal solitons in single-mode fibers, as well as, in a theoretical form, for 2D spatiotemporal“ light bullets” in planar waveguides. We have here produced a family of 3D solitons in the model combining theGRIN structure of the refractive index in the transverse plane, approximated by the HO (harmonic-oscillator, i.e.,quadratic) profile, Kerr self-focusing nonlinearity, and the usual DM map, based on periodic alternation of anomalous-and normal-GVD segments. It is found that the stability of the spatiotemporal DM solitons is determined by the DMstrength parameter, S : periodically oscillating DM solitons self-trap from localized inputs, as fully stable modes, at S > S (cid:39) .
93. After a transient rearrangement of the input, stable DM solitons keep constant energy. At
S < S ,the simulations demonstrate quasi-stability: the self-trapping gives rise to spatiotemporal solitons with persistentsmall-amplitude random intrinsic vibrations, which very slowly lose their energy, remaining robust over propagationdistances corresponding to hundreds of DM periods. Stable three-dimensional DM solitons with embedded vorticitywere constructed too. Collisions between DM solitons were considered by boosting them in opposite axial directions,with a conclusion that the collisions are quasi-elastic. Collisions in the transverse plane were also addressed, byinitially placing one or both solitons at off-axis positions, and letting them roll down under the action of the HOpotential. In this case, the pair of solitons, which perform the shuttle motion in the confining potential, demonstratea periodic sequence of quasi-elastic collisions, in the absence or presence of DM.The present work can be extended by considering the system including DM combined with the nonlinearity manage-ment. This extension may be relevant because different segments of multimode fibers in the composite waveguide mayhave different values of the nonlinearity parameter. As mentioned above, the realistic DM model may also include theintra-pulse stimulated-Raman-scattering effect and third-order dispersion, which deserves the consideration. Anotherdirection for the development of the analysis may be the systems of the WDM type, with two or several distinctcarrier wavelengths, modeled by a system of nonlinearly coupled NLSEs. V. ACKNOWLEDGMENTS
We appreciate valuable discussions with T. Birks and L. G. Wright. This work was supported, in part, by theThailand Research Fund (grant No. BRG 6080017), Israel Science Foundation (grant No. 1286/17), and RussianFoundation for Basic Research (grant No. 17-02-00081). TM acknowledges the support from Faculty of Engineering,Naresuan University, Thailand.
VI. REFERENCES [1] Malomed B A, Mihalache D, Wise D, and Torner L 2005
Journal of Optics B: Quantum and Semiclassical Optics R53–R72[2] Malomed B A, Mihalache D, Wise D, and Torner L 2016
Journal of Physics B: Atomic, Molecular and Optical Physics Chaos Physics Reports
Advances In Atomic, Molecular, and Optical Physics Reports on Progress in Physics The European Physical Journal Special Topics
Romanian Reports in Physics Physical Review Letters
Nature Reviews Physics Nature Communications Optical Solitons: From Fibers to Photonic Crystals (Academic Press)[13] Pitaevskii L P and Stringari S 2003
Bose–Einstein Condensation (Oxford University Press)[14] Fibich G 2015
The Nonlinear Schr¨odinger Equation: Singular Solutions and Optical Collapse (Heidelberg: Springer)[15] Malomed B A (INVITED) 2019
Physica D
Physics Letters A [17] Silberberg Y 1990 Optics Letters Bell System Technical Journal IEEE Journal of Quantum Electronics
QE-10
AEU - International Journal of Electronics and Communications Applied Physics Letters Optics Letters Optics Express Journal of Lightwave Technology Journal of Lightwave Technology Nature Photonics Advances in Optics and Photonics IEEE Communications: Surveys and Tutorials Optics Express Optics Communications
Optics Letters Optics Communication
Nature Photonics Optics Letters Nature Communications Optica Physical Review Letters
Optics Express Light:Science & Applications e16208[40] Guenard R, Krupa K, Dupiol R, Fabert M, Bendahmane A, Kermene V, Desfarges-Berthelemot A, Auguste J L, TonelloA, Barth´el´emy A, Millot G, Wabnitz S, and Couderc V, 2017 Optics Express IEEE Photonics Technology Letters Nature Photonics Optica Optics Express Soliton Management in Periodic Systems (Springer: New York).[46] Turitsyn S K, Bal B G, and Fedoruk M P 2012
Physics Reports
Optics Letters Japanese Journal of Applied Physics L889-L891[49] Smith N J, Knox F M, Doran N J, Blow K J, and Bennion I 1996
Electronics Letters Electronics Letters Optics Letters Physical Review E Physical Review E IEEE Journal of Selected Topics in Quantum Electronics Optics Express Physical ReviewLetters Physical Review. Letters
Journal of the Optical Society of America B Physical Review A Physical Review Letters Physical Review Letters Physical Review A Journal of the Optical Society of America B Physical Review A Physical Review Letters Applied Optics Applied Optics Applied Optics Chinese Optics Letters Optics & Laser Technology
Journal of Lightwave Technology Journal of Lightwave Technology Optics Communications , 145-157[74] Karlsson M, Anderson D, and Desaix M 1992
Optics Letters Optics Letters and Wise F W, 2020 Nature Physics Pure and Applied Optics L57–L62[78] Lakoba T I and Agrawal G P 1999
Journal of the Optical Society of America B Optics Letters , 808-810[80] Yang J 2009 Journal of Computational Physics
Radiophysics and Quantum Electronics Physics Reports
Optics Letters Journal of the Optical Society of America B Nature Physics Photonics Research Physical Review E58