Light bullets in coupled nonlinear Schrödinger equations with variable coefficients and a trapping potential

We analyze three-dimensional (3D) vector solitary waves in a system of coupled nonlinear Schrödinger equations with spatially modulated diffraction and nonlinearity, under action of a composite self-consistent trapping potential. Exact vector solitary waves, or light bullets (LBs), are found using the self-similarity method. The stability of vortex 3D LB pairs is examined by direct numerical simulations; the results show that only low-order vortex soliton pairs with the mode parameter values n ≤ 1, l ≤ 1 and m = 0 can be supported by the spatially modulated interaction in the composite trap. Higher-order LBs are found unstable over prolonged distances. © 2017 Optical Society of America OCIS codes: (190.0190) Nonlinear optics; (190.6135) Spatial solitons. References and links 1. B. Malomed, L. Torner, F. Wise, and D. Mihalache, “On multidimensional solitons and their legacy in contemporary Atomic, Molecular and Optical physics,” J. Phys. At. Mol. Opt. Phys. 49(17), 170502 (2016). 2. M. Blaauboer, B. A. Malomed, and G. Kurizki, “Spatiotemporally localized multidimensional solitons in selfinduced transparency media,” Phys. Rev. Lett. 84(9), 1906–1909 (2000). 3. L. Bergé, “Wave collapse in physics: principles and applications to light and plasma waves,” Phys. Rep. 303(56), 259–370 (1998). 4. D. Abdollahpour, S. Suntsov, D. G. Papazoglou, and S. Tzortzakis, “Spatiotemporal Airy light bullets in the linear and nonlinear regimes,” Phys. Rev. Lett. 105(25), 253901 (2010). 5. D. J. Frantzeskakis, H. Leblond, and D. Mihalache, “Nonlinear optics of intense few-cycle pulses: An overview of recent theoretical and experimental developments,” Rom. J. Phys. 59, 767–784 (2014). 6. Y. Silberberg, “Collapse of optical pulses,” Opt. Lett. 15(22), 1282–1284 (1990). 7. L. Bergé, “Wave collapse in physics: principles and applications to light and plasma waves,” Phys. Rep. 303(56), 259–370 (1998). 8. L. C. Crasovan, Y. V. Kartashov, D. Mihalache, L. Torner, Y. S. Kivshar, and V. M. Pérez-García, “Soliton “molecules”: robust clusters of spatiotemporal optical solitons,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 67(4), 046610 (2003). 9. K. Gallo, A. Pasquazi, S. Stivala, and G. Assanto, “Parametric solitons in two-dimensional lattices of purely nonlinear origin,” Phys. Rev. Lett. 100(5), 053901 (2008). 10. X. Liu, K. Beckwitt, and F. Wise, “Generation of optical spatiotemporal solitons,” Opt. Photonics News 82, 4631–4634 (1999). 11. D. E. Edmundson and R. H. Enns, “Particlelike nature of colliding three-dimensional optical solitons,” Phys. Rev. A 51(3), 2491–2498 (1995). 12. Y. F. Chen, K. Beckwitt, F. W. Wise, and B. A. Malomed, “Criteria for the experimental observation of multidimensional optical solitons in saturable media,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 70(4), 046610 (2004). 13. O. Bang, W. Krolikowski, J. Wyller, and J. J. Rasmussen, “Collapse arrest and soliton stabilization in nonlocal nonlinear media,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 66(4), 046619 (2002). 14. W. Królikowski and O. Bang, “Solitons in nonlocal nonlinear media: exact solutions,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 63(1 Pt 2), 016610 (2001). Vol. 25, No. 8 | 17 Apr 2017 | OPTICS EXPRESS 9094 #287083 https://doi.org/10.1364/OE.25.009094 Journal © 2017 Received 24 Feb 2017; revised 23 Mar 2017; accepted 28 Mar 2017; published 10 Apr 2017 15. D. Mihalache, D. Mazilu, L. C. Crasovan, I. Towers, A. V. Buryak, B. A. Malomed, L. Torner, J. P. Torres, and F. Lederer, “Stable spinning optical solitons in three dimensions,” Phys. Rev. Lett. 88(7), 073902 (2002). 16. A. T. Avelar, D. Bazeia, and W. B. Cardoso, “Solitons with cubic and quintic nonlinearities modulated in space and time,” Phys. Rev. E Stat. Nonlin. Soft Matter Phys. 79(2), 025602 (2009). 17. N. K. Efremidis, J. Hudock, D. N. Christodoulides, J. W. Fleischer, O. Cohen, and M. Segev, “Two-dimensional optical lattice solitons,” Phys. Rev. Lett. 91(21), 213906 (2003). 18. Y. V. Kartashov, V. A. Vysloukh, and L. Torner, “Rotary solitons in bessel optical lattices,” Phys. Rev. Lett. 93(9), 093904 (2004). 19. S. Nixon, L. Ge, and J. Yang, “Stability analysis for solitons in PT-symmetric optical lattices,” Phys. Rev. A 85(2), 023822 (2012). 20. D. Mihalache, D. Mazilu, F. Lederer, B. A. Malomed, Y. V. Kartashov, L. C. Crasovan, and L. Torner, “Stable spatiotemporal solitons in bessel optical lattices,” Phys. Rev. Lett. 95(2), 023902 (2005). 21. A. V. Gorbach and D. V. Skryabin, “Spatial solitons in periodic nanostructures,” Phys. Rev. A 79(5), 053812 (2009). 22. A. Szameit, I. L. Garanovich, M. Heinrich, A. Minovich, F. Dreisow, A. A. Sukhorukov, T. Pertsch, D. N. Neshev, S. Nolte, W. Krolikowski, A. Tünnermann, A. Mitchell, and Y. S. Kivshar, “Observation of diffractionmanaged discrete solitons in curved waveguide arrays,” Phys. Rev. A 78(3), 031801 (2008). 23. T. X. Tran and F. Biancalana, “Diffractive resonant radiation emitted by spatial solitons in waveguide arrays,” Phys. Rev. Lett. 110(11), 113903 (2013). 24. C. Weilnau, M. Ahles, J. Petter, D. Träger, J. Schröder, and C. Denz, “Spatial optical (2+1)-dimensional scalarand vector solitons in saturable nonlinear media,” Ann. Phys. 11(8), 573–629 (2002). 25. S. V. Manakov, “Theory of two-dimensional stationary self-focusing of electromagnetic waves,” Sov. Phys. JETP 38, 248–253 (1974). 26. W. Torruellas, Y. S. Kivshar, and G. I. Stegeman, Quadratic Solitons (Springer, 2001). 27. W.-P. Zhong and M. Belic, “Resonance solitons produced by azimuthal modulation in self-focusing and selfdefocusing materials,” Nonlinear Dyn. 73(4), 2091–2102 (2013). 28. V. S. Bagnato, D. J. Frantzeskakis, P. G. Kevrekidis, B. A. Malomed, and D. Mihalache, “Bose-Einstein condensation: Twenty years after,” Rom. Rep. Phys. 67, 5–50 (2015). 29. D. S. Wang, S. W. Song, B. Xiong, and W. M. Liu, “Quantized vortices in a rotating Bose-Einstein condensate with spatiotemporally modulated interaction,” Phys. Rev. A 84(5), 053607 (2011). 30. A. I. Yakimenko, Y. A. Zaliznyak, and V. M. Lashkin, “Two-dimensional nonlinear vector states in BoseEinstein condensates,” Phys. Rev. A 79(4), 043629 (2009). 31. A. I. Yakimenko, K. O. Isaieva, S. I. Vilchinskii, and M. Weyrauch, “Stability of persistent currents in spinor Bose-Einstein condensates,” Phys. Rev. A 88(5), 051602 (2013). 32. A. I. Yakimenko, S. I. Vilchinskii, Y. M. Bidasyuk, Y. I. Kuriatnikov, K. O. Isaieva, and M. Weyrauch, “Generation and decay of persistent currents in a toroidal Bose-Einstein condensate,” Rom. Rep. Phys. 67, 249– 272 (2015). 33. S. L. Xu, M. R. Belić, and W. P. Zhong, “Three-dimensional spatio-temporal vector solitary waves in coupled nonlinear Schrödinger equations with variable coefficients,” J. Opt. Soc. Am. B 30(1), 113–122 (2013). 34. S. L. Xu, G. P. Zhou, N. Z. Petrović, and M. R. Belić, “Nonautonomous vector matter waves in two-component Bose-Einstein condensates with combined time-dependent harmonic-lattice potential,” J. Opt. 17(10), 105605 (2015). 35. S. L. Xu, J. X. Cheng, M. R. Belić, Z. L. Hu, and Y. Zhao, “Dynamics of nonlinear waves in two-dimensional cubic-quintic nonlinear Schrödinger equation with spatially modulated nonlinearities and potentials,” Opt. Express 24(9), 10066–10077 (2016). 36. S. L. Xu, N. Z. Petrovic, M. R. Belić, and W. W. Deng, “Exact solutions for the quintic nonlinear Schrödinger equation with time and space,” Nonlinear Dyn. 84(1), 251–259 (2016). 37. A. S. Desyatnikov and Y. S. Kivshar, “Necklace-ring vector solitons,” Phys. Rev. Lett. 87(3), 033901 (2001). 38. S. L. Xu and M. R. Belic, “Light bullets in coupled nonlinear Schrödinger equation with spatially modulated coefficients and Bessel trapping potential,” J. Mod. Opt. 62(9), 683–692 (2015). 39. W. P. Zhong, M. Belić, G. Assanto, B. A. Malomed, and T. Huang, “Self-trapping of scalar and vector dipole solitary waves in Kerr media,” Phys. Rev. A 83(4), 043833 (2011). 40. R. Hirota, “Exact Solution of the Korteweg-de Vries Equation for Multiple Collisions of Solitons,” Phys. Rev. Lett. 27(18), 1192–1194 (1971). 41. D. Zwillinger, Handbook of Differential Equations, 3rd ed. (Academic, 1997). 42. M. Dehghan and D. Mirzaei, “The dual reciprocity boundary element method (DRBEM) for two-dimensional sine-Gordon equation,” Comput. Methods Appl. Mech. Eng. 197(6-8), 476–486 (2008). 43. B. A. Malomed, D. Mihalache, F. Wise, and L. Torner, “Spatiotemporal optical solitons,” J. Opt. B Quantum Semiclassical Opt. 7(5), R53–R72 (2005). 44. T. I. Lakoba, “Stability analysis of the split-step Fourier method on the background of a soliton of the nonlinear Schrödinger equation,” in Methods for Partial Differential Equations, 4974, 641–649 (2010). Vol. 25, No. 8 | 17 Apr 2017 | OPTICS EXPRESS 9095


Introduction
Three-dimensional optical spatiotemporal solitary waves, also known as light bullets (LBs), are self-sustained wave packets that are localized in both spatial dimensions and time [1][2][3][4]. They may form in dielectric media in which diffraction and dispersion are balanced by medium's nonlinearity. The generation of LBs is a nontrivial task from the analytical and numerical points of view, and even more complicated in real experimental settings [5,6].
It is well known that optical solitons appear in the nonlinear Schrödinger equation (NLSE) with cubic self-focusing optical nonlinearity, which models a wealth of phenomena in nonlinear science in general and physics in particular. Owing to the occurrence of optical beam collapse, optical solitons are unstable in two and three dimensions (2D, 3D) [7,8]. However, in the past decades, the stability of solitons in 2D and 3D has been improved, using different media and methods. One typically uses quadratic and cubic nonlinear optical media that support stable optical solitons regardless the physical dimension in which they form and propagate [8,9]. In 1999, a seminal experimental study was reported by the group of Wise [10], in which, using the achromatic phase matching technique and generating the necessary anomalous GVD, robust (2 + 1)D nonlinear LBs were formed in quadratic nonlinear crystals. We also mention several other physical settings which are adequate for getting robust LBs. These include (a) either saturable [11,12] or nonlocal [13,14] optical materials, (b) optical media with competing quadratic, cubic or cubic-quintic nonlinearities [15,16], (c) confining by two-or three-dimensional optical lattices [17][18][19][20], and (d) periodic waveguide structures with controlled diffraction and/or GVD [21][22][23].
Spatial solitons can be generally divided into scalar solitons (single-component) and vector solitons (multicomponent), according to the number of field components and the nature of the medium [24]. Here, we are interested in the vector type, arising in (3 + 1) dimensions in a Kerr medium with variable coefficients and a trapping potential. It is well known that vector solitons commonly form in the absence of interference between components. In general, there exist many ways to generate vector solitons, such as crossphase modulation (XPM) [25] or considering beams of different wavelengths [26] and using mutually incoherent beams which consist of two-component azimuthons [27]. In nonlinear optics, the XPM-mediated interaction between mutually incoherent or orthogonally polarized waves leads to the formation of bound states that are also vector solitons.
A field in which spatiotemporal solitons appear regularly is the Bose-Einstein condensation [28]. Recently, a theoretical analysis and numerical study of quantized vortices in a rotating Bose-Einstein condensate (BEC) with spatiotemporally modulated interaction in a harmonic potential, was reported in [29]. In [30], 2D vector matter waves in the form of soliton-vortex and vortex-vortex pairs have been investigated for the case of attractive intracomponent interaction in two-component BECs. Further, in two-component BECs, the superflow of atomic spinor BECs optically trapped in a ring-shaped geometry were considered in [31,32]. In our previous studies, 3D approximate but still analytical spatiotemporal vector LBs were built with the help of spherical harmonics, including multipole solutions and necklace rings [33]. In [34], we constructed exact self-similar soliton solutions of 3D coupled Gross-Pitaevskii equations for two-species BECs in a combined time-dependent harmoniclattice potential. We have investigated the control and manipulation of solitary waves for three kinds of BECs with changing diffraction and nonlinearity coefficients; the solutions include Ma breathers, and Peregrine and Akhmediev solitons [34].
However, thus far no exact 3D spatial vector LBs in two-component NLSE with spatially modulated coefficients and a trapping potential have been found, either theoretically or experimentally. Here, we extend the analysis performed in [33][34][35][36], to demonstrate such exact solutions of the coupled NLSE in (3 + 1)D with varying coefficients and a self-consistent trapping potential.
The paper is organized as follows. The coupled (3 + 1)D NLSE with the modulated coefficients and a composite trapping potential are analyzed in Sec. II, where the nonlinearity coefficients are allowed to vary in space. Also presented in Sec. II are the vortex-vortex LB pairs, obtained using the self-similarity method. The stability of these vortex LB pairs is discussed in Sec. III. Characteristic distributions of LBs are provided in Sec. IV. The concluding remarks with a simple summary are given in Sec. V.

The model and exact soliton solutions
In this paper we consider vector solitons consisting of N mutually incoherent optical components in a nonlinear Kerr medium. The propagation of the slowly-varying field components ( , , , ) j u r z ϕ θ ( 1, 2, 3 j N = ⋅⋅⋅ ) is described by a system of coupled dimensionless (3 + 1)D NLSEs with changing coefficients, under the influence of an optical trapping potential [37, 38]: is the spatiotemporal radius, andϕ is the azimuthal coordinate in the transverse ( , ) x y plane. The evolution variable is the propagation distance z, not time. The retarded timeτ is the third 'transverse' coordinate in this coordinate system.
where l is the separation constant, cast in a convenient form. As a consequence of the selfconsistency condition, the total intensity I = I(r,z) should be radially symmetric at any z; hence, in our search for vector solutions of Eq.
(1) we are restricted to the solutions that maintain this property. Equation (2)  We chose the simplest case of 3D vector solitons, with two components, 2 N = . The complex coefficients 1 2 , a a and 1 2 , b b are calculated from Eqs. (4) and (5) below, and satisfy the following relations [38,39]: where the parameter [0,1] q ∈ determines the modulation depth of the beam [29]. Hence, the coefficients a n and b n depend on q only.
Our aim in finding analytical solutions of Eq. (3) is to connect them with the stationary NLSE To find exact solutions of Eqs. (8)-(11), we introduce a self-similarity transformation [33,35,40] 3 2 1 Here, κ is the normalization constant, ( ) w z is the beam width, 1 ( , ) z r θ is the similarity variable to be determined, and ( ) a z and ( ) b z are the wave front curvature and the phase offset, respectively. These variables vary with the distance z. Substituting the presumed solutions for the amplitude and the phase into Eqs.
Here,η is constant and n is assumed to be a non-negative integer, which can be considered as the radial mode number. Finally, the trapping potential V is defined more closely; it is composed of a harmonic and an inverse harmonic term, in addition to a term depending on the transformed radial function R. Thus, V is determined self-consistently with the solution presumed. Differential Eq. (12) The third variable z in the Laplacian, in this case is τ. This eigenvalue problem is computed by the Fourier Collocation Method for discrete eigenvalues, which allows the whole spectrum to be calculated at once [43]. Once parameters m , n , l , q andη are specified, one gets the imaginary and real parts of δ from Eqs. (13) and (14). If imaginary parts of all eigenvalues δ are equal to zero or are positive, the soliton solutions can be stable; otherwise, the perturbed solution would grow exponentially with z, and thus, the corresponding matter waves would become linearly unstable [35,43].

Characteristic distributions of LBs
From Eq. (13), it is evident that the matter wave solitons are characterized by five parameters: the mode numbers n, m, l, the harmonic potential width η, and the modulation depth q. In this section, we present typical distributions of LBs, for some values of these parameters. When 0 w w = , the auxiliary function R, the nonlinearity coefficient χ , and the trapping potential V are the functions of the radial variable r only. χ is localized in the form of asymmetric double annuluses; when 1 m = , localized three amplitude peaks can be seen in Fig. 1(a). In form, the external potential V (r) is similar to the 2D pattern in [35]. Figures 1(c) and 1(d) illustrate the amplitude of the radial field distributions j ψ for different mode numbers 0,1, 2 m = (c) and n = 1, 2, 3 (d). The optical field of LBs approaches 0 when 0 r = and r = ∞ , forming different vortex structures. The ring size depends on the mode number n. The results show that the solutions in Eq. (13) indeed are localized. Figure 2 displays the intensity distributions of 3D LBs with 1 n m = = , 0.4 η = , and 0.1 q = , for different l ( 1, 2,3 l = ) from top to bottom, at 70 z = . A variety of intensity distributions for solitons with different l are shown. For 1 l = there are two layers (two pulses) along the vertical τ -axis, and the soliton is shaped as a pair of deformed ellipsoids. In general, by increasing l, the solution forms 1 l + layers in the vertical direction, forming multipole multipulse solitons. However, the total intensity distribution (the right column in Fig. 2) is radially symmetric. Obviously, when both parameters (l and m) are 0, the soliton is the fundamental spherical LB. The physical origin of this arrangement can qualitatively be understood from the nature of the self-focusing cubic nonlinearity as a third-order susceptibility (3) χ ; it yields the nonlinearity polarization of the medium. Owing to the assumption of strong dispersion along the verticalτ axis, resulting in spherical harmonics, the associated Legendre polynomial 0 1 P will appear in the solution; it possesses a maximum value (plus one) and a minimum value (minus one). Therefore, the spherical harmonics will force a change of the sign near the zero point.  , and for different 0, 2,3 m = from top to bottom, also at 70 z = . The first row just depicts the same three figures, because m = 0. For given l, the larger the m, the larger the LB radius in the horizontal plane. Looking at the field components, the number of solitary beads in the same layer is 2m; it comes from the azimuthal angular dependence in Eq. (13). Because 1 q ≠ and 0 m ≠ , there are 2m maxima and minima in the range 0 2 mϕ π ≤ ≤ . In this sense, the ringnecklace components are produced steadily as m increases. The larger the m, the shorter the pulse in the vertical direction. The total vector intensity is still radially symmetric, as it must be.   from left to right. It is seen that the profiles of LBs possess polycyclic structure with several amplitude peaks covering the same ring. LBs are self-similar and composed of four symmetrical petals in each ring (the bottom row in Fig. 4). The number of outer bright rings is 1 n + , the optical intensity in the center is zero, and the intensity of rings surrounding the center decreases with the increasing radial distance.
Increasing the modulation depth q from 0 to 1, with fixed n, m, and l, will lead to less azimuthally modulated vortex rings (Fig. 5). One can observe that the intensity and the phase of the beam is modulated by the modulation depth q. When 0 q = , for each component we find a necklace LB, with the number of beads equal to 2m. By increasing q, the distance between petals decreases, and the multi-TC vortices change into vortex rings. It is apparent that the soliton has just one layer in the vertical direction.  1 m = , 2 n = . Other parameters are the same as in Fig. 1. The stability of the obtained exact solutions is a very important aspect of the present discussion, which has to be addressed numerically. In order to check the stability of Lbs, we perform direct numerical simulations, using the split-step Fourier method [44] and solve Eq.
(1) by taking the analytical solution (13) at 0 z = as an initial condition. In Fig. 6, we present the comparison of analytical (the first column) and numerical (the second and the third columns) intensity contour plots. The linear-stability spectra are displayed in the fourth column of Fig. 6. It is found that only when the topological charge is 0,1 m = and 1 n l = = , the numerical solutions of LBs are stable against perturbation with an initial Gaussian noise level of up to 10%. But when the topological charge is 1 m > , and mode numbers n and l are greater than one, the soliton solution (13) is unstable in propagation and splits into spreading quasi square-shaped band structures. Further, a quadruple complex eigenvaluesδ are visible in Figs. 6(h) and 6(i), thus these solitons might be linearly unstable.

Conclusion
In summary, we have studied the propagation of LBs in coupled nonlinear Schrödinger equations with spatially modulated coefficients and a self-consistent trapping potential. Using the self-similarity method, exact vector solitary waves, or light bullets, are obtained. The stability of vortex 3D LB pairs is examined by direct numerical simulation. Our results show that the only stable low-order vortex LB pairs with 1 n ≤ , 1 l ≤ and 0 m = can be supported by the spatially modulated interaction in the composite trap. Higher-order LBs are found unstable over prolonged distances. Our approach may prove useful for other related problems, e.g., light propagation in Bose-Einstein condensates and plasmas.