Breather stripes and radial breathers of the two-dimensional sine-Gordon equation
P.G. Kevrekidis, R. Carretero-González, J. Cuevas-Maraver, D.J. Frantzeskakis, J.-G. Caputo, B. A. Malomed
aa r X i v : . [ n li n . PS ] J u l Breather stripes and radial breathers of the two-dimensional sine-Gordon equation
P. G. Kevrekidis,
R. Carretero-Gonz´alez, J. Cuevas-Maraver,
D. J. Frantzeskakis, J.-G. Caputo, ∗ and B. A. Malomed
8, 9, † Department of Mathematics and Statistics, University of Massachusetts Amherst, Amherst, MA 01003-4515, USA Mathematical Institute, University of Oxford, Oxford, OX2 6GG, UK Nonlinear Dynamical Systems Group, ‡ Computational Sciences Research Center,and Department of Mathematics and Statistics, San Diego State University, San Diego, California 92182-7720, USA Grupo de F´ısica No Lineal, Universidad de Sevilla. Departamento de F´ısica Aplicada I,Escuela Polit´ecnica Superior, C/ Virgen de ´Africa, 7, 41011 Sevilla, Spain Instituto de Matem´aticas de la Universidad de Sevilla. Avda. Reina Mercedes, s/n. Edificio Celestino Mutis, 41012 Sevilla, Spain Department of Physics, National and Kapodistrian University of Athens, Panepistimiopolis, Zografos, Athens 15784, Greece Laboratoire de Math´ematiques, INSA de Rouen Normandie76801 Saint-Etienne du Rouvray, France. Department of Physical Electronics, 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
We revisit the problem of transverse instability of a 2D breather stripe of the sine-Gordon (sG) equation. Anumerically computed Floquet spectrum of the stripe is compared to analytical predictions developed by meansof multiple-scale perturbation theory showing good agreement in the long-wavelength limit. By means of directsimulations, it is found that the instability leads to a breakup of the quasi-1D breather in a chain of interacting 2Dradial breathers that appear to be fairly robust in the dynamics. The stability and dynamics of radial breathersin a finite domain are studied in detail by means of numerical methods. Different families of such solutionsare identified. They develop small-amplitude spatially oscillating tails (“nanoptera”) through a resonance ofhigher-order breather’s harmonics with linear modes (“phonons”) belonging to the continuous spectrum. Theseresults demonstrate the ability of the 2D sG model within our finite domain computations to localize energy inlong-lived, self-trapped breathing excitations.
PACS numbers:
I. INTRODUCTION
One of the most influential models in studies of solitarywaves is the sine-Gordon (sG) equation, which has beenextensively explored in numerous volumes [1–3] and re-views [4]. The one-dimensional (1D) version is the first non-linear equation whose integrability was found (in the form ofthe B¨acklund transform, ca. 140 years ago [5]). Later, its com-plete integrability was systematically investigated by means ofthe inverse scattering transform (IST) [6]. Commonly knownexact solutions to the sG equation are topological kink soli-tons and kink-antikink bound states in the form of breathers.Modes of the latter type are uncommon in continuum systems,but find broad realizations in their discrete counterparts [7, 8].While the 1D sG equation has become a textbook etalonof integrable models, far less is known about higher dimen-sional versions of the same equation, which are not integrableby the IST. In particular, some effort has been dedicated tothe kinematics and dynamics of kinks in two and three dimen-sions (2D and 3D) [9–12], including their ability to producebreathers as a result of collisions with boundaries [13], andproneness to be stably pinned by local defects [14].Nevertheless, despite the relevance of physical realizations ‡ URL : http://nlds.sdsu.edu ∗ Email:[email protected] † Email:[email protected] of the sG model in higher dimensions (see, e.g., Ref. [15]for a recent example), the study of breathers in such settingswas scarce. While effective equations of motion for breatherswere derived long ago [16], these equations do not apply tosteady isotropic breather profiles. Thus, it remains unknownat present whether such robust solutions may exist indefinitelylong, or they eventually decay.Moreover, there is a relevant related question concerningquasi-1D breather solutions in the higher-dimensional prob-lem, i.e., breather stripes and planes in the 2D and 3D settings,respectively. The seminal work of Ref. [17] had predicted thatsuch quasi-1D solutions are subject to transverse instabilities.This is a natural property due to the relation of the sG model,in the limit when it produces broad small-amplitude breathers,to the nonlinear Schr¨odinger (NLS) equation [1, 2], in whichboth bright and dark quasi-1D solitons are prone to transverseinstabilities [18] —see, e.g., theoretical work of Ref. [19] andRefs. [20–22] for experimental realizations in atomic and op-tical physics. On the other hand, quasi-1D sG kinks are not vulnerable to the transverse instability. In this work, we aimto further explore the dynamics of quasi-1D sG breathers inthe 2D geometry, including their transverse instability and itsdevelopment into radially shaped breathing modes.We start the analysis by revisiting the key result of Ref. [17]for quasi-1D breathers in the 2D sG equation. We test thisresult for breathers with different frequencies ω b and demon-strate that their instability in the long-wavelength limit is ac-curately captured by the analytical treatment. Our numeri-cal calculation of the respective Floquet multipliers (FMs) ad-ditionally permits the systematic identification of spectra oftransverse-instability modes as a function of the correspond-ing perturbation wavenumber, k y . These results indicate thatthe quasi-1D sG breathers do not only maintain a main bandof unstable wavenumbers, similarly to what is known for theNLS equation, where only a band of low wavenumbers is un-stable. Importantly, they also display, for appropriate frequen-cies ω b , instability bubbles for higher wavenumbers. Near theNLS limit of ω b → , we demonstrate how the NLS reduc-tion, as well as a novel, to the best of our knowledge, varia-tional approximation involving the transverse degree of free-dom are able to capture instability and dynamical features, in-cluding the basic necking phenomenon occurring in the courseof the breather stripe’s evolution.Studying the dynamics ensuing from the generic long-wavelength instabilities of the breather stripes, we identifyspontaneous formation of localized breather waveforms. Fur-ther, “naively” initializing a quasi-1D radial breather profile,with x replaced by the radial variable, r , we find that the so-lutions adjust themselves and their frequencies in the courseof a transient stage of the evolution, to form potentially verylong-lived localized modes trapping the energy in a localizedspatial region. This finding motivates us to study the radialtime-periodic breathers, which we are able to find, as genuinesteady states over a period, using fixed-point iteration tech-niques. Stability of these radial breathers is studied numeri-cally in the framework of the Floquet theory.The manuscript is organized as follows. Section II reportsthe results concerning the breather stripes. In particular weexpand the instability results for the stripes from Ref. [17].The results, which are valid from the long-wavelength limit tolarger wavenumbers, are obtained by connecting the sG dy-namics to the NLS limit, using multiple-scale expansions, andalso by dint of a variational approach, which is relevant forlarge ω b . We present numerical results demonstrating that thedevelopment of the breather-stripe instability nucleates long-lived radial breathers. The existence and stability of such ra-dial breathers is the subject of Sec. III. Finally, in Sec. IVwe present a summary of the results, concluding remarks, andpoint out avenues for further research. II. BREATHER STRIPESA. Modulational instability of breather stripes
The sG equation is a special case of the more general Klein-Gordon class of models which in 2D is written as u tt − ∇ u + V ′ ( u ) = 0 , (1)where u ( x, y, t ) is a real field depending on coordinates ( x, y ) and time t , subscripts with respect to the independent vari-ables denote partial derivatives, and V ( u ) is the potential thatdefines the particular model. The most studied among thesemodels pertain to V ( u ) = 1 − cos( u ) and V ( u ) = ( u − / ,which are referred to as the sG [3] and “ φ ” [23] equations,respectively. From now on, we focus on the former one: u tt − ∇ u + sin u = 0 . (2) The 1D version of the sG equation admits commonly knownkink/antikink and breather solutions. Kinks are (generallyspeaking, traveling) 1D profiles that asymptote to the uniformbackground states, u ±∞ = { , π } , of the form u k = 4 arctan h exp (cid:16) s ( x − ct ) / p − c (cid:17)i , (3)where s = ± is the topological charge, which distinguisheskinks and antikinks, and c is the velocity, which may take val-ues | c | < . On the other hand, breathers correspond to pro-files which are oscillatory in time and localized in space. Inparticular, the exact 1D sG breather solution is u b ( x, t ) = 4 tan − (cid:20) βω b sech ( β x ) cos ( ω b t ) (cid:21) , (4)where β ≡ p − ω and the band of the breather frequenciesis ≤ ω b < .In the spirit of the work done in Ref. [14] for a kink stripe,we consider the transverse instability of the breather stripe u ( x, y, t ) = u b ( x, t ) . (5)We then perturb the breather stripe according to u ( x, y, t ) = u b ( x, t ) + w ( x, y, t ) where the perturbation w is assumed small. Enforcing thissolution to satisfy Eq. (2) yields, to first order, the followingevolution equation for the perturbation: w tt − ( w xx + w yy ) + cos( u b ( x, t )) w = 0 . (6)This equation is a linear, inhomogeneous in x , wave equationso that we can study each transverse wave number k y sepa-rately and assume w ( x, y, t ) = ξ ( x, t ) exp( ik y y ) , (7)where ξ ( x, t ) describes the x -dependent shape (eigenmode)of the perturbation. Plugging (7) into (6), yields ξ tt − ξ xx + (cid:0) k y + cos( u b ( x, t ) (cid:1) ξ = 0 . (8)This is a periodically forced wave equation and can be solvedusing Floquet theory, see Appendix for details. The final resultis the linear operator M (cid:18) { ξ ( x, T ) }{ ξ t ( x, T ) } (cid:19) = M (cid:18) { ξ ( x, }{ ξ t ( x, } (cid:19) , (9)relating the solutions of Eq. (8) at t = 0 and t = T . If the max-imal Floquet multiplier Λ (obtained from the eigenvalues of M ) is such that | Λ | > , ξ will grow and the stripe will be un-stable. In the case of instability (where at least one eigenvaluehas modulus larger than 1) the eigenmode ψ corresponding tothe maximum | Λ | gives the solution ξ ( x, t ) and the modula-tion of the stripe through its perturbation w ( x, y, t ) .In the seminal work of Ref. [17], Ablowitz and Kodama(AK) predicted that, in the long-wavelength limit, breatherstripes are subject to modulational instability (MI). Their mainresult provides the instability growth rate λ for the stripe, as afunction of the perturbation wavenumber, k y : λ = σk y , σ = p − ω ω b arcsin (cid:18)q − ω (cid:19) . (10)Given that the AK result is valid in the k y → limit, wewish to investigate the instability for larger values of k y . Tothis end, one needs to solve Eq. (8) numerically, by means ofa suitable fixed-point method (see details in the Appendix),taking into regard that a finite size (and a finite discretizationparameter) of the solution domain induces gaps in the spec-trum of the linear modes and, also that breathers may fea-ture modified tails generated by “hybridization” (mixing) ofthe third and higher-order odd harmonics with linear modes(“phonons”), reminiscent of the formation of the so-called“nanoptera” in the φ equation [3]. The hybridization de-pends on the location of the linear modes, which, in turn, de-pends on the specifics of the domain size and boundary condi-tions. We have used periodic boundary conditions in domaininvolving x ∈ [ − L x , L x ) with L x = 100 and computed theperiodic steady states corresponding to zero initial velocity u t ( x, y, t = 0) = 0 . Finally, it is important to note that, ad-dressing the dynamics with the underlying time-periodic so-lution, one needs to resort to Floquet analysis and the compu-tation of FMs, in order to extract the (in)stability eigenvalues(see details in Appendix). The FMs, Λ , are translated into(in)stability eigenvalues λ (in particular, for the comparisonto the AK theory), via Λ = exp( λT ) , where T = 2 π/ω b isthe breather’s period. FIG. 1: (Color online) Modulational-instability eigenvalues for thebreather stripe vs. wavenumber k y of the modulational perturbations,for the indicated values of the breather’s frequency. Solid curvesdepict the eigenvalues numerically determined by unstable Floquetmultipliers (see the text). The dashed line represents the analyticalAblowitz-Kodama (AK) result for the long-wavelength limit, k y → [17], given by Eq. (10). Figure 1 depicts the dependence on k y of the MI growthrate for selected values of the breather stripe’s frequency, ω b .The eigenvalues have been rescaled by σ [see Eq. (10)] to fa-cilitate the comparison to the AK analytical prediction, given FIG. 2: (Color online) Snapshots of the breathers u b ( x, t = 0) andcorresponding eigenfunctions ψ for different MI windows/bubbles.The left and right sets of panels correspond to ω b = 0 . and ω b = 0 . , respectively. Snapshots of the breathers u b ( x, t = 0) aredisplayed in the top panels. The subsequent ones, from top to bottom,display the u (left) and u t (right) components of the eigenfunctions ψ corresponding to the highest-MI wavenumbers in the bubbles ofFig. 1. For ω b = 0 . , two MI bubbles exist (with maxima located at k y ≈ . and 1.28), while for ω b = 0 . there are three MI bubbles(with maxima located at k y ≈ . , . and 1.30). The computa-tion domain is x ∈ [ − , , half of which is shown here. by Eq. (10) and shown by the dashed line in the figure. TheMI spectrum indeed follows the AK prediction as k y → ,for all values of ω b , i.e., for long perturbation wavelengths,the MI eigenvalue is λ/σ ≈ k y , independent of ω b . Further,Fig. 1 demonstrates that, beyond the AK limit, λ grows with k y until reaching a maximum, and then decreases, falling tozero at some k y = ˜ k y . We refer to this first instability window( < k y < ˜ k y ) as the main MI “bubble”. The value ˜ k y , aswell as the maximum of λ/σ , increase with ω b for ω b ≤ . ,but this trend changes to a decrease for . < ω b < . Thisresult, which is in stark contrast with properties of MI for theNLS equation [24, 25], has its origin in increasing contribu-tions from the third and higher harmonics of the breather’soscillatory waveform when p − ω decreases and, conse-quently, truncating at the first harmonic is no longer a goodapproximation for expression (4) for ω b . . .Furthermore, beyond the main MI bubble, Fig. 1 suggestsanother noteworthy feature: the existence of secondary MIbubbles for k y > ˜ k y . These secondary bubbles get narrowerand more spaced with the increase of ω b ; in fact, they are notfound for ω b > . . To understand MI in the secondary bub-bles, and to compare it to the main MI window, k y < ˜ k y , weexplore the shape of the corresponding eigenfunctions of themodulational perturbations. In particular, Fig. 2 depicts thebreather profile and eigenfunctions corresponding to unstableeigenvalues in the main MI window and secondary bubbles for ω b = 0 . and ω b = 0 . . In the case of ω b = 0 . (the set ofleft panels in the figure), MI takes place in the main windowand in one secondary bubble, while in the case of ω b = 0 . (the set of right panels) two secondary MI bubbles are present.All perturbation eigenfunctions in these MI bubbles are real asthe corresponding eigenvalue Λ responsible for the instabilityis also real (we depict the respective components u and u t inthe left and right plots, respectively). Figure 2 exhibits eigen-functions that bear qualitatively similar shapes. Namely, theeigenfunctions in the main MI window are all localized with-out oscillating tails, while in the secondary bubbles, tails areattached to the core of the eigenfunctions, getting stronger inhigher-order bubbles. B. Connection to NLS
To better understand the MI spectrum for the breatherstripes, we may use the connection between the sG and theNLS equations, in the limit of ω b → . Such a connection canbe accurately established via the multiple-scales perturbationmethod [26]. Below we will briefly present the methodologyand provide the main results (see also Refs. [1, 2] for a de-tailed analysis). First we note that, in the limit of ω b → , thesG breather of Eq. (4) has a small amplitude and, hence, onecan expand the sG nonlinearity as sin u ≈ u − u / . Thus,we seek for solutions of the resulting φ -equation in the formof the asymptotic expansion: u ( x, t ) = N X n =1 ǫ n u n ( x , . . . , x N , t , . . . , t N ) , (11)where u n are functions of the independent variables x j = ǫ j x and t j = ǫ j t ( j = 0 , , , . . . N ), while < ǫ ≪ is a formalsmall parameter. Substituting ansatz (11) in the equation, andcollecting coefficients in front of different powers of ǫ yieldsa set of equations, the first three ones being O ( ǫ ) : L u = 0 , (12) O ( ǫ ) : L u + L u = 0 , (13) O ( ǫ ) : L u + L u + L u − u = 0 . (14)where L ≡ ∂ t − ∂ x + 1 , L ≡ ∂ t ∂ t − ∂ x ∂ x ) , L ≡ ∂ t − ∂ x + 2 ( ∂ t ∂ t − ∂ x ∂ x ) . To leading order, O ( ǫ ) , the solution to Eq. (12) is u = A ( x , x , . . . , t , t , . . . ) exp( i Φ) + c . c , (15)where A is a yet-to-be-determined function of the slow vari-ables, c . c stands for the complex conjugate, and Φ = kx − ωt , with wavenumber k and frequency ω obeying the lineardispersion relation, ω = k + 1 . (16) Substituting solution (15) in Eq. (13), it is evident that thesecond term is secular and has to be removed. This yieldsequation A t + v g A x = 0 , where the v g = k/ω is the sGgroup velocity, as per Eq. (16). This result suggests that theamplitude function A depends on variables x and t onlythrough a traveling coordinate, X = x − v g t , (17)i.e., A = A ( X, x , x , . . . , t , t , . . . ) . Furthermore, as con-cerns the solution of Eq. (13) [which is now reduced to L u = 0 ], we take the trivial solution u = 0 , because anontrivial one would be of the form as Eq. (15), hence its am-plitude function can be included in field A . Finally, at the nextorder, O ( ǫ ) , upon substituting u from Eq. (15) and u = 0 ,Eq. (14) becomes: L u = (cid:2)(cid:0) − v g (cid:1) A XX + 2 iω ( A t + v g A x ) (cid:3) exp( iθ )+ 16 A exp(3 iθ ) + 12 | A | A exp( iθ ) + c . c = 0 . (18)It is observed that secular terms ∝ exp( iθ ) arise on the right-hand side of Eq. (18) [term ∝ exp(3 iθ ) is not secular, be-cause it is out of resonance with the uniform solution]. Re-moving the secular terms leads to the following NLS equationfor complex amplitude A : iA t + (cid:0) ω (cid:1) − A XX + (4 ω ) − | A | A = 0 , (19)where we have removed the group-velocity term by redefining x x − v g t , and made use of identity − v g = 1 /ω ,resulting from Eq. (16).The fact that the NLS equation (19) is valid in the band ω ≥ , according to Eq. (16), allows us to consider the limitof ω → , i.e., k → . In this case, with the vanishing groupvelocity v g → , variable X , defined by Eq. (17), carries overinto x , Eq. (19) reduces to the form of iA t + (1 / A x x + (1 / | A | A = 0 . (20)The stationary soliton solution of the reduced NLS, expressedin terms of original variables x and t , is A = η sech[(1 / ǫηx )] exp[( i/ ǫ η t ] , (21)where η is an arbitrary O (1) parameter. Thus, an approximatesG breather solution, valid up to order O ( ǫ ) near the edge ofthe phonon band, is produced by Eqs. (15) and (21): u ( x, t ) ≈ ǫη sech (cid:18) ǫηx (cid:19) cos (cid:20)(cid:18) − ǫ η (cid:19) t (cid:21) . (22)Note that the frequency of soliton solution (22), ω s = 1 − ( ǫ η / , is smaller than the cutoff frequency, ω = 1 , whichnaturally means that this self-trapped excitation belongs to thephonon bandgap , as is known for the sG breather. This fact in-dicates strong connection between the stationary NLS solitonand the sG breather in the limit of < − ω b ≪ . (23) FIG. 3: (Color online) Comparison between the (scaled) instabilityspectrum of the sG breather stripe (shown by colored dots) and theinstability spectrum for a stationary NLS bright-soliton stripe (shownby the gray solid line, as per Refs. [24, 25]). The wavenumber andinstability growth rate are rescaled as: ˜ k y ≡ k y / p − ω and ˜ λ ≡ λ/ (1 − ω ) . Indeed, in this limit, observing that the argument of tan − in solution (4) is small, it may be approximated by u b ≈ (4 β/ω b ) sech( βx ) cos( ω b t ) . Then, letting ǫη ≡ β , whichmeans that ω b ≈ − ( ǫ η / ≡ ω s , we see that the approx-imate NLS soliton solution, given by Eq. (22), is identical tothe sG breather in the limit case defined by Eq. (23). Accord-ingly, it is expected that the breather stripe’s dynamics can beadequately approximated by the NLS equation (19), where ω is identified as the breather stripe’s frequency, ω b .In particular, the connection between the NLS equation andthe sG stripe allows us to approximate the stability spectrumof the sG breather stripe by that for the NLS bright-solitonstripes, which was studied in detail previously [24, 25, 27].Specifically, in Fig. 3 we compare the rescaled numeri-cally computed instability spectra of the sG breather stripe(dots) for different breather frequencies from region (23) withthose for a stationary NLS bright-soliton stripe, taken fromRefs. [24, 25]. It is evident that, in the limit of ω b → , the(scaled) instability spectrum of the sG stripe indeed smoothlyapproaches its NLS counterpart. C. Variational approach for large frequencies
In this Section we present a study of the evolution ofsG breather stripes, based on the variational approximation,which is valid for relatively large breather’s frequencies. Tothis end, we note that the 2D sG equation (2) can be derivedfrom the Lagrangian, L = RR L dx dy , with density L = 12 (cid:0) u t − u x − u y (cid:1) − (cid:16) u (cid:17) . (24)We rewrite the exact 1D sG breather solution (4) as u b = 4 tan − [(tan µ ) sin ( t cos µ ) sech ( x sin µ )] , (25) where parameter µ ≡ arcsin (cid:16)p − ω (cid:17) takes values ≤ µ ≤ π/ . We now focus on broad low-amplitude breathers inlimit of µ ≪ π/ , approximating the exact 1D breather (25)by the following ansatz : u b = χ ( y, t ) sech ( µx ) . (26)Here, χ is a variational parameter accounting for the pos-sibility of variations of the width of the stripe in the y di-rection, which may represent the so-called necking instabil-ity [27, 28]. According to the VA, evaluating the Lagrangianwith ansatz (26) substituted in Eq. (24), and deriving theEuler-Lagrange equation from it [29] should yield an effec-tive evolution equation for the width parameter χ ( y, t ) . In thepresent form, the Lagrangian density (24) cannot be integratedanalytically for ansatz (26). Nonetheless, focusing on small-amplitude breather stripes, one may expand the nonlinear termto simplify the Lagrangian density: L = 12 (cid:0) u t − u x − u y (cid:1) − u + 124 u . (27)Upon substituting ansatz (26) in Eq. (27), the integration canbe performed in the x direction, yielding the effective La-grangian: µL eff = Z + ∞−∞ (cid:18) χ t − χ y − χ + 118 χ (cid:19) dy. (28)Finally, the variational equation for χ ( y, t ) , which governsthe evolution of necking perturbations of the sG stripe un-der the condition (23), is immediately derived from the La-grangian (28): χ tt − χ yy + (cid:18) − µ (cid:19) χ − χ = 0 . (29)Before we compare the reduced dynamics, presented byEq. (29), to the evolution of the full sG breather stripe, weprovide a generalization of the VA methodology to incorpo-rate the full nonlinearity of the system. If instead of approx-imating sin( u ) up to third order as done above, one rewritesthe nonlinearity in the effective Lagrangian as sin ( u/
2) =(1 / (cid:2) − P ∞ m =0 ( − m u m / (2 m )! (cid:3) and uses the identity, Z + ∞−∞ sech m ( x ) dx = 2 m − [( m − (2 m − , (30)it is possible to evaluate the effective Lagrangian for the fullsinusoidal nonlinearity with ansatz (26). The respective ex-tended Euler-Lagrange equation is χ tt − χ yy − µ χ − ∞ X m =1 ( − m (cid:20) ( m − m − (cid:21) (2 χ ) m − = 0 . (31)Although Eq. (31) generalizes Eq. (29) in that it takes the fullnonlinearity into account, it does not provide an exact the-ory, because the underlying ansatz (26) is only valid for smallamplitudes. Thus, the extended equation may provide for a FIG. 4: (Color online) Stability spectra of the sG breather (thicklines) and as produced by reduced model (31) (thin lines) for dif-ferent breather’s frequencies. The top panel depicts the spectrumproduced by Eq. (31) with only N t = 2 [i.e., by Eq. (29)], instead ofthe infinite sum in the full sG equation. The second, third, and fourthpanels display results produced by Eq. (31), keeping N t = 3 , , and , respectively. a more accurate VA model of the necking dynamics, but thismodel retains the approximate nature of the above analysis(based on the proposed ansatz).To gain insight into the validity of the VA models (29) and(31), in Fig. 4 we compare stability spectra provided by theseapproximations with the numerically found main instabilitywindow of the full sG model. Different panels correspond toreplacing the infinite sum in Eq. (31) by the truncated sum P N t m =1 , that keeps the first N t terms in the expansion. In thisnotation, Eq. (31) reduces to Eq. (29) when N t = 2 . Theresults presented in Fig. 4 suggest that the MI is indeed ad- equately captured by the VA in the region (23). This is truequalitatively in that region and even quantitatively as ω b → .It is seen too that, as − ω b increases, the match between thefull sG model and its VA counterparts deteriorates. Nonethe-less, the VA reproduces the correct overall trend and the shapeof the instability spectra. As regards keeping more terms inthe expansion of the nonlinearity in the Lagrangian, Fig. 4shows that proceeding from N t = 2 to N t = 3 does not pro-duce a significant improvement in the accuracy, and there isno discernible difference between N t = 4 and N t > ei-ther. Therefore, from now on, we use the simplest reducedmodel (29) with N t = 2 for the comparison with the full sGmodel. D. Dynamical evolution of breather stripes
Having tested the reduced models in terms of the stabilityspectra, we now extend the comparison by numerically fol-lowing the necking dynamics in the framework of the full sGmodel and the VA for a number of scenarios. The evolutionwas initiated by inputs u ( x, y,
0) = u b ( x, (cid:20) ε cos (cid:18) n y πL y y (cid:19)(cid:21) ,u t ( x, y,
0) = 0 , (32)where n y is the perturbation’s integer mode number (associ-ated with the wavenumber n y π/L y ) and ε its relative strength,while u b ( x, is the exact 1D sG breather profile (4). Re-sults of the comparison are depicted in Fig. 5. The dynam-ics is depicted by means of surface plots of u ( x, y, t ) at theindicated times. The input was a slightly perturbed exactsG breather stripe, to which one [(a)–(f)] or several [(g) and(h)] perturbation modes with wavenumbers k y were added asfollows. The perturbations were constructed by seeding thenecking modulation as per Eq. (32), with perturbation strength ε = 0 . . Panels (a)–(d) depict the comparison for the pertur-bations with small wavenumbers n y = 1 , , , and , corre-sponding to the necking MI. Panels (e) and (f) correspond tolarge perturbation wavenumbers, with n y = 17 and n y = 19 ,which are taken, respectively, just below and above the in-stability threshold. For the stable case (f) with n y = 19 , astronger perturbation with ε = 0 . was used to stress the sta-bility of the sG breather against this perturbation wavenumber.Finally, panels (g) and (h) correspond to cases when a mix ofperturbation modes with random strengths (chosen in interval − . ≤ ε ≤ . ) is introduced in the initial conditions.Case (g) contains perturbations with modes ≤ n y ≤ ,while in case (h) we used ≤ n y ≤ . In all the cases,to match the underlying oscillation frequencies of the exactbreather stripe and its VA counterpart, the time for the VAresults is scaled by a fitting factor. This factor, for all thecases shown here, amounts to a reduction of the time be-tween and . As seen in the figure, and as predictedby the stability spectra, perturbations with small wavenum-bers always lead to necking MI. The perturbation eventuallybreaks the stripe into a chain of self-trapped “blobs” (local-ized modes). These time-oscillating “blobs” are identified as FIG. 5: (Color online) Comparison of the development of the modulational (necking) instability for the sG breather stripe, as produced by thefull sG model, and by VA, in the form of Eq. (29), for ω b = 0 . . The corresponding set of top and bottom subpanels respectively depict theevolution of u ( x, y, t ) for the full sG dynamics and the VA. The input was a slightly perturbed exact sG breather stripe, to which one [(a)–(f)]or several [(g) and (h)] perturbation modes with wavenumbers k y are added. The system was integrated in the domain | x | ≤ , | y | ≤ (only the segment | x | ≤ is shown), with periodic boundary conditions along y and homogeneous Dirichlet boundary conditions along x .Please see the main text for more details. highly perturbed radial breathers, this conclusion being thesecond key aspect of our study (see the following Section III,and also Ref. [13]). All the modes with perturbation indices ≤ n y ≤ ( k y = n y π/L y ) [see Eq. (32)] are unstable,and the VA is able to accurately capture the corresponding in-stability and ensuing dynamics, see panels (a)–(e). Naturally,as the necking perturbation grows, the VA becomes less ac-curate, as the underlying variational ansatz (26) is no longeran appropriate approximation for the profile of the deformedstripe. Nonetheless, for times up to those at which the stripebreaks into a chain of “blobs”, the VA evolves in tandem withthe full sG dynamics. This is also true when, instead of usinga single mode as the perturbation, one perturbs the stripe by acollection of modes, as shown in panels (g) and (h) where theperturbations contain, respectively, a combination of the low-est modes and even as many as ones. Finally, in panel(f) of Fig. 5 we present an example of the evolution with theperturbation index ( n y = 19 ) chosen beyond the instabilityboundary of the main MI window. In this stability case, thesimulations demonstrate that the full sG dynamics is closelyfollowed by the VA for longer times, in comparison to the un-stable cases. Figure 6 expands on the results presented in Fig. 5 by con-sidering random perturbation of the breather stripes and al-lowing the dynamics to evolve for longer times. Interest-ingly, Fig. 6 shows the formation of robust spatially localized“blobs”. The top row of panels depicts snapshots of the field u ( x, y, t ) at the indicated moments of time. Since the solu-tions that we are inspecting are oscillatory in time, in the cor-responding bottom row we plot snapshots of the correspond-ing 2D energy density, E ( x, y, t ) = 12 (cid:0) u t + u x + u y (cid:1) + V ( u ) , (33)which allows us to partially eliminate the oscillations, andthus to better capture the emergence and dynamics of the“blobs” (localized breathers). The figure suggests that thebreathers generated by the MI from the unstable stripes aregenerically robust, although some of them disappear, whileothers spontaneously emerge from condensation of the back-ground random field. Note that the localized breathers traveland, as a result, collide or merge with neighboring ones. InFig. 7, we display 3D isocontour plots of the energy densityin the ( x, y, t ) spatiotemporal continuum. This figure makes FIG. 6: (Color online) Long term evolution of randomly perturbed breather stripes showing the generation of “blobs” (perturbed radial sGbreathers) as a result of necking MI. Depicted is the corresponding field u ( x, y, t ) and energy density E ( x, y, t ) (top and bottom plots,respectively) at times indicated in the panels for (a) ω b = 0 . , (b) ω b = 0 . , and (c) ω b = 0 . . The input is the exact sG breather stripeperturbed by a uniformly distributed random perturbation of amplitude . . The domain size corresponds to | x | , | y | ≤ (only the segment − ≤ x ≤ is shown). it possible to better follow the emergence, evolution, and in-teractions of the “blobs”. III. RADIAL BREATHERS
We now turn to the consideration of breathers with axialsymmetry, i.e., u ( x, y, t ) = u ( r, t ) , with r = p x + y ina circular domain of radius R . Their existence is suggestedby the findings reported in the previous Section, where neck-ing MI splits the sG breather stripe into a chain of persistent localized modes (“blobs”). A. Constructing radial breathers
With the aim to construct initial conditions that directlylead to radial breathers, we take, as an input, the exact 1Dbreather (4), with the x coordinate replaced by the radial one: u radial ( r,
0) = u b ( r, . The results are summarized in Figs. 8and 9, which depict, respectively, isocontours of the energydensity, and dynamics and power spectra of the field at the FIG. 7: (Color online) Isocontour plots of the energy density E ( x, y, t ) corresponding to the evolution of slightly perturbed sGbreather stripes shown in Fig. 6. In all cases, random perturbationsadded to the initial stripe trigger the development if its transverseMI. As a result, several robust “blobs” emerge, whose motion andinteractions dominate subsequent dynamics.FIG. 8: (Color online) Isocontour plots depicting the evolution ofthe energy density E ( x, y, t ) for radial waveforms initialized by the1D breather profile, with x replaced by r , for (a) ω b = 0 . , (b) ω b = 0 . , and (c)-(d) ω b = 0 . . The waveforms reshape themselvesinto radial breathers, which persist over an indefinitely large numberof oscillations. In (a) and (b), the isocontour levels correspond to E ( x, y, t ) = E , with E = 2 . ; in (c), three cuts are shown for E = 0 . , . , and . (more transparent isocontours correspond tosmaller E ). Apparently, the radial breather for ω b = 0 . seems todisperse. However, an isocontour drawn for an even smaller value, E = 0 . , in (d) suggests the presence of a much shallower radialbreather. central point, u (0 , , t ) . The results suggest that persistent ra-dial breathers are generated by appropriately crafted inputs.However, the resulting radial breathers do not keep the oscil-lation frequency corresponding to the the one seeded by theinitial condition. In fact, for ω b = 0 . and ω b = 0 . , theensuing radial breathers are very similar (in amplitude andfrequency), but they feature an oscillation frequency higherthan the one introduced by the input, see green vertical dashedlines in Figs. 9(d)–(f). Note that, as observed in panel (f) ofFig. 9, the input with larger ω b leads to a much shallower ra-dial breather with a dominant frequency close to , i.e., theedge of the phonon band. Similar results were obtained for ω b = 0 . and ω b = 0 . (not shown here). FIG. 9: (Color online) Time series depicting the evolution of fieldat the central point u (0 , , t ) (top panels) and the corresponding nor-malized power spectra ˆ u (0 , , t ) (bottom panels), for the cases shownin Fig. 8. For the cases with ω b = 0 . (a) and ω b = 0 . (b) theresulting radial breathers are very similar in terms of the field ampli-tude [see plots for u (0 , , t ) ] and dominant frequency [see plots for ˆ u (0 , , ω ) ]. In panel (c), corresponding to the case with ω b = 0 . , amuch shallower radial breather emerges with a dominant frequencytending to ω = 1 , which is the edge of the phonon band. The corre-sponding dominant frequencies are ω ≈ . in (a), ω ≈ . in (b), and ω ≈ . in (c). Red vertical dashed lines denotefrequencies ω b of the 1D sG breather used as the input. The numerical experiments shown above might suggest thatradial breathers exist as exact solutions. However, on thecontrary to the 1D sG equation, its 2D counterpart is notintegrable. Therefore, genuine radial breathers cannot existin the infinite spatial domain, as multiple breathing frequen-cies, nω with ω < but n > /ω , resonate with the fre-quencies belonging to the phonon spectrum, ω > , whichshould give rise to the radiative decay of the breather (onlyintegrable equations make an exception, nullifying the de-cay rate) [1, 2, 4]. However, in some cases, accurate simu-0 FIG. 10: (Color online) Numerically computed profiles of radialbreathers for ω b = 0 . , . , . , . , . (from top to bottom). Thesolution domain is ≤ r ≤ , only the central portion r ≤ being shown. All these profiles have oscillatory tails with a smallamplitude (see insets for ω b ∈ . , . , . ). lations show that the rate of the radiative decay is so smallthat the 2D sG breathers survive thousands of oscillation pe-riods, with a negligible loss of amplitude [10]. Neverthe-less, if the domain is finite , there will be gaps in the spec-trum, which admit the existence of intraband breathers simi-lar to the phantom breathers in discrete lattices introduced inRef. [30] and feature the aforementioned nanopteron wingsoscillating in space. For example, we plot in Fig. 10 profilesof numerically exact radial breathers with with frequencies ω b = 0 . , . , . , . , . , computed in the radial domain ≤ r ≤ and with zero initial velocity u t ( r, t = 0) = 0 (see below for more details of these numerics). As it is madeevident by the figure, the radial breathers indeed have distinc-tive oscillatory tails, which are more prominent for lower val-ues of ω b , and their amplitude is decaying to zero, in line withthe expectation predicted above by the NLS reduction, in thelimit of ω b → . Indeed, the 2D NLS equation with the cu-bic self-focusing term creates radial 2D Townes solitons ; see,e.g., Refs. [31, 32] for a review and Refs. [33, 34] for theirexperimental realizations, an earlier one in optics and a re-cently created Townes soliton in atomic gases. Although 1Dbreathers are also intraband ones and possess a wing in thecase of the finite-domain computation, the continuation of 2Dradial breathers within the band is a much more subtle prob-lem. A full analysis of the continuation (bifurcation) scenarioin the whole range of available frequencies, for which the non-
FIG. 11: (Color online) Existence branches of radial breathers for ω b close to . (a–c): Total energy E of radial breathers versus thebreather’s frequency, ω b . Panel (a) presents an overview of the differ-ent solution branches (see the text), while panels (b) and (c) displayzoomed in versions of the transition zones between A-B and B-Cbranches, respectively. Representative profiles for the solutions be-longing to each branch, corresponding to the large red dots, are pre-sented in Figs. 12 and 13. Vertical dashed lines denote the most rel-evant frequency values: (i) dashed red lines, corresponding to ω r,A and ω r,B , denote the interval outside of which the breather acquires aminimum at r = 0 ; (ii) dashed black lines, corresponding to ω A and ω B , denote the turning points for branches A and B, respectively;(iii) the dashed green line at ω b = ˜ ω denotes the location of the en-ergy minimum for branch B, above which the breathers belonging tothis branch are prone to exponential instabilities. The bottom paneldepicts the maximum amplitude of the main Fourier coefficients ofthe breathers belonging to branch B. All the computations were per-formed in region r ≤ R = 200 . integrability for the 2D equation must play a crucially impor-tant role, is outside the scope of the present work. Therefore,in the summary of the stability results that we present below,we focus on branches of radial breathers for large ω b , viz ., ω b > . , in line with our earlier consideration of the limitof ω b → . An advantage of this option is also that, whendealing with larger ω b , the numerics (performed over one pe-riod through the Floquet analysis, see below) are more man-ageable (faster) than for smaller frequencies.1 B. Radial breathers solution branches and stability
We now consider the stability of a “radial breather” u b ( r, t ) .For this purpose, we proceed as in the breather stripes section.First, we write the PDE (2) in polar coordinates and obtain theequation for the perturbation w corresponding to Eq. (6): w tt − r ( rw r ) r − r w θθ + cos( u b ( r, t )) w = 0 . (34)We then assume a solution w ( r, θ, t ) = ζ ( r, t ) exp( ik θ θ ) , (35)with integer angular wavenumber k θ , yielding the final time-periodic PDE for the perturbation ζ : ζ tt − r ( rζ r ) r + (cid:18) k θ r + cos( u b ( r, t )) (cid:19) ζ = 0 . (36)To deal with the formal singularity of the radial Laplacian at r = 0 we do not include this point in the domain so that theclosest point to r = 0 is r = h/ [that is, r ∈ ( h/ , R ) ], andassume u r ( r = 0) = 0 . This way of handling the singularityis standard, see for example Ref. [35], and it is equivalent toconsidering u ( r = − h/
2) = u ( r = h/ . (37)Under this scheme of the spatial discretization, and using zeroDirichlet boundary conditions at r = R (i.e., fixed edges), weemploy the standard shooting method in time —implementedby dint of a fourth-order explicit and symplectic Runge-Kutta-Nystrom method developed in Ref. [36], with time step equalto T / for the shooting method and
T / for the Floquetanalysis (see Appendix)— and seek for periodic solutions, foreach given temporal period. Then, for each periodic solutionwe compute the corresponding stability spectrum via Floquetanalysis (see Appendix for details).To gain insight into the stability (and the ensuing dynam-ics) of radial (intraband) breathers and their hybridization withphonon modes, we have considered a small interval around ω b = 0 . . We have verified that for the cases under consider-ation, where ω b is close to , the Fourier coefficients associ-ated with the breathers decay relatively fast with the increaseof their order. This allows us to resolve the solutions withsatisfactory resolution. Namely, we have corroborated thatincreasing the spatial resolution by including more Fouriermodes does not lead to any visible change in the obtained sta-bility spectra.The main stability results for radial breathers with frequen-cies close to ω b = 0 . are summarized in Fig. 11. Specifi-cally, panels (a)–(c) depict the total energy [ E = RR E dxdy ,with energy density defined as per Eq. (33)] versus ω b . Here,three distinct branches (A, B, and C) of the radial branches arefound. In analogy to the 1D phantom breathers of Ref. [30],there is branch A to the left, which terminates at a turningpoint, ω b = ω A = 0 . , and central branch B thatalso finishes at a turning point, ω b = ω B = 0 . .Herein, we mostly focus on branch B. Panels (b) and (c) depict FIG. 12: (Color online). Radial breather’s profiles at the points P ,P , P , and P , depicted by large red dots in panel (b) of Fig. 11.FIG. 13: (Color online). Radial breather’s profiles at the points Q ,Q , Q , and Q , depicted by large red dots in panel (c) of Fig. 11. zoomed areas from (a) around the frequencies where, respec-tively, branch A turns and branch B starts, or branch B turnsand branch C starts. These two bifurcations occur, due to theresonance of the seventh harmonic of the ω b with phonons,respectively, at points ω phonon / ω A b = 0 . and ω phonon / ω B b = 0 . . These phonon frequencies2 FIG. 14: (Color online) Largest Floquet multipliers for k θ = 0 (thesolid blue line) and k θ = 1 (the red dashed line) for breathers be-longing to branch B. The strongest instability, which takes place at ω b > ˜ ω = 0 . , is exponential (monotonously growing), otherinstabilities being oscillatory. produce, via the dispersion relation (16), the wavenumber ob-served in the tails of the breather, confirming that the reso-nance between the breather and phonons takes place. The leftand central branches join, via turning points, through addi-tional solution branches ( A ′ and B ′ , respectively, see Fig. 11).Note that breathers belonging to branch B- B ′ have the am-plitude maximum displaced from r = 0 (that is, the solutionsare shaped as ring breathers ) at ω b < ω r,A = 0 . and ω b > ω r,B = 0 . . Indeed, Figs. 12 and 13 dis-play the breather’s profiles —again for zero initial velocity u t ( r, t = 0) = 0 — close to the transitions between branches, viz ., respectively, A ⇄ B and B ⇄ C , and, specifically, inpanel (c) of Fig. 12 and panel (b) of Fig. 13, the profiles donot have a maximum at r = 0 . This qualitative change of thebreathers’ shape can be understood by noting that, as men-tioned above, the seventh harmonic of these breathers is res-onant with phonon waves, its amplitude being ≃ . timesthe amplitude of the fundamental (first) harmonic. As a con-sequence, the breather progressively resembles a delocalizedphonon wave when the frequency is decreased. Thus, the ringsin the breather’s shape at ω b < ω r,A are a consequence of thehybridization with phonons.Finally, we address the stability of the radial breather solu-tions. Our numerical Floquet stability analysis suggests thatthe breathers belonging to branch B are mainly affected by ra-dial perturbations, i.e., those with k θ = 0 in Eq. (35), and thatthe only exponential instabilities that occur are related to theexistence of minima in the energy-vs.-frequency dependence,i.e., ones of the kind predicted in Ref. [37]. Because of this,breathers belonging to branch B in Fig. 11 are subject to ex-ponential instability at ω b > ˜ ω = 0 . , where ˜ ω is theenergy-minimum point on B branch, denoted by the verticalgreen dashed line in the figure. The full stability spectrum forthe radial-breather solutions in branch B is shown in Fig. 14.As expected, the breather becomes exponentially unstable toradial perturbations (see the solid blue line representing the FM corresponding to radial perturbations) at ω b > ˜ ω . In ad-dition, the breathers belonging to this branch are also sub-ject to an oscillatory instability at ω b > . , althoughits magnitude is small ( | Λ | . . ). We have checked thatthis weak oscillatory instability is not a numerical artifact, aschanging k m [controlling the number of Fourier modes, perEq. (41)] or the integration step in the Floquet analysis donot alter the spectrum. Changing the discretization parame-ter would also lead to the appearance of new phonon modesand, at the same time, new resonances that would change theexistence interval of the radial breathers. In any case, theseweak oscillatory instabilities should not have a noticeable ef-fect on the breathers’ dynamics. The oscillatory instabilitieswere only found for k θ = 1 in Eq. (35). The breathers belong-ing to branch B ′ are exponentially unstable against the pertur-bations with k θ = 0 at all values of ω b , and feature weakoscillatory instabilities for k θ ≥ . A brief exploration of thestability for breathers belonging to branch A indicates that theresults are very similar to those reported here for branch B. IV. CONCLUSIONS
We have studied the extension of sG breathers to 2D, inthe form of quasi-1D breather stripes and localized radialbreathers. Starting from the long-wave MI (modulational in-stability) of the stripe, discovered in the seminal work byAblowitz and Kodama [17] and confirmed herein, we have ex-panded the MI analysis to arbitrary wavelengths of the trans-verse perturbations. In the limit of small-amplitude broadbreathers, with frequency ω b → , we have employed theasymptotic multiscale expansion method to match the dynam-ics of the breather stripes to that of bright soliton stripes inthe framework of the 2D NLS equation. Perhaps more im-portantly, we have developed a novel version of the VA (vari-ational approximation), whereby the 2D dynamics of the sGbreather stripe is reduced to a filament-type evolution equa-tion for the amplitude of the breather along the stripe. The VAallowed us to formulate an approximate reduced model thatnot only captures the necking MI (i.e., the dependence of theMI gain on the perturbation wavenumber), but is also able topredict the dynamics beyond the linear instability setting, andup to close to the breakup of the stripe by the growing neckingperturbations.It is precisely the necking MI of sG breather stripes thatresults in the nucleation of temporally oscillatory, spatiallylocalized “blobs” in the form of radial breathers. Numeri-cal results demonstrate that the blobs can collide or mergewith nearby ones. The apparent robustness of these modeshas prompted us to study their existence and stability in moredetail. It is important to note that, contrary to the 1D sG equa-tion, its 2D counterpart is not integrable, hence genuine radialbreathers cannot exist in infinite domains for infinitely longtimes due to resonances of multiple harmonics of the breatherwith the phonon (continuous) spectrum. Nonetheless, in thefinite domain it is possible to find intraband breathers, simi-larly to how it was done for 1D phantom breathers that werediscovered in Ref. [30]. In fact, our results show that such in-3traband radial breathers, possessing a “nanopteron” spatiallyoscillatory tail, do exist and may be stable in the finite do-main. The full bifurcation structure of such radial breatherswith arbitrary values of the oscillation frequency, ω b , is quitecomplex. Nonetheless, we have performed a detailed numer-ical study for the existence of broad small-amplitude radialbreathers for ω b close to , in the vicinity of some of theunderlying resonances. We have identified several solutionbranches and their stability, by means of numerical Floquetanalysis, finding that such breathers may indeed be stable inthe finite domain.There are numerous avenues deserving further study. Forinstance, a more complete characterization of the bifurcationscenarios for arbitrary values of ω b may be able to eluci-date further nontrivial bifurcations and destabilization scenar-ios. Also, the dynamics and interactions of radial breathersclearly merit further study: for instance, a challenging ob-jective is to investigate the outcome of interactions betweenradial breathers (including collisions between moving ones),as a function of their relative phase, velocity, and frequen-cies. Lastly, it is natural to generalize the consideration for 3Dsettings, and examine whether spherically symmetric radialbreathers may survive in a finite domain, and/or whether theymay be produced spontaneously via MI of quasi-2D breatherplanes embedded in 3D. These topics are presently under con-sideration, and results will be reported elsewhere. Acknowledgments
This material is based upon work supported by the US Na-tional Science Foundation under Grants PHY-1602994 andDMS-1809074 (P.G.K.). P.G.K. also acknowledges supportfrom the Leverhulme Trust via a Visiting Fellowship andthanks the Mathematical Institute of the University of Ox-ford for its hospitality during this work. R.C.G. gratefullyacknowledges support from the US National Science Foun-dation under Grant PHY-1603058. J.C-M. was supported byproject P18-RT-3480 (Regional Government of Andalusia).B.A.M. appreciates support provided by Israel Science Foun-dation through grant Np. 1286/17.
Appendix: Numerical implementation for periodic states andtheir Floquet stability analysis
In this Appendix, we outline the Floquet analysis thatwas deployed for quasi-1D breathers embedded in rectangu-lar domains and radial breathers in circular domains. Thesebreathers are T -periodic solutions corresponding to steadystates of the dynamics after time evolves from t = 0 to t = T .Furthermore, as these solutions only depend (non-trivially) ona single spatial coordinate —the x - or r -direction for stripeand radial breathers respectively—, they can be numericallycomputed as 1D solutions along this longitudinal direction.Nonetheless, it is crucial to note that the stability for thesesolutions has to be computed in the full
2D domain wherethese breather structures are embedded in. Namely, perturba- tions to these solutions need to be followed, not only alongthe longitudinal direction, but more importantly, along thetransverse direction. To this end, as described below, we fol-low the growth of plane-wave transverse perturbations withwavenumbers k y and k θ , along the y - and θ -directions for,respectively, the breather stripes in the Cartesian plane andradial breathers in polar coordinates. Thus, each wavenum-ber, k y or k θ , yields a 1D perturbation equation determiningthe stability of the corresponding modulations over one period T that is treated, as detailed below, using Floquet analysis toobtain its corresponding eigenvalues (eigenfrequencies) andeigenvectors (eigenmodes).The first step in our numerical computations is to discretizein space the solution and its corresponding derivatives. To thisend, we consider a finite-difference scheme with a (spatial)discretization parameter h = 0 . , which transforms the PDEfor u ( x, t ) into a set of N coupled ODEs for u n ≡ u ( x n , t ) onthe discrete grid { x n } . The discretized version of the quasi-1D sG equation (also known as the Frenkel-Kontorova model)reads: ¨ u n + sin( u n ) + 1 h ( u n +1 − u n + u n − ) = 0 , (38)where n = − N/ . . . N/ . Similarly, for the radial breatherthe corresponding discretized radial sG equation may be writ-ten as ¨ u n + sin( u n ) + 1 h ( u n +1 − u n + u n − )+ 12 nh ( u n − − u n +1 ) , n = 1 . . . N (39)To produce breathers in the numerical form, we have madeuse of two different techniques, based on the fact that the solu-tions are T -periodic, with T = 2 π/ω b : (i) a shooting method,based on the consideration of the map, Y (0) → Y ( T ) , Y ( t ) = (cid:20) { u n ( t ) }{ ˙ u n ( t ) } (cid:21) , (40)and (ii) a Fourier-transform implementation, based on ex-pressing the solution of the discretized dynamical equationsin the form of a truncated Fourier series: u n ( t ) = z ,n + 2 k m X k =1 z k,n cos( kωt ) , (41)with k m being the maximum of the absolute value of the k in our truncation of the full Fourier series. In the numerics, k m = 11 was chosen. Note that the spatial evenness of the sGpotential [ V ( u ) = 1 − cos( u ) ], the Fourier coefficients witheven k are zero and, consequently, only odd harmonics of ω b can resonate with phonons.After the introduction of Eq. (41) in the dynamical equa-tions, one gets a set of N × ( k m + 1) nonlinear, coupled al-gebraic equations. For a detailed explanation of these meth-ods, the reader is referred to Refs. [38, 39]. In both meth-ods, the continuation in the frequency is implemented viathe path-following (Newton-Raphson) method. The Fourier-transformed methods have the advantage, among others, of4providing an explicit analytical form of the Jacobian, whichmakes the calculations faster. However, for frequenciessmaller than ∼ . , the Fourier coefficients decay slowly anda large number of coefficients should be kept to produce anaccurate solution, which makes the Jacobian numerically ex-pensive. Because of this, we developed the analysis based onmethod (i) to quasi-1D breathers in the range . ≤ ω b < ,whereas the Fourier-transform methods (ii) were applied tothe radial breathers, when we focused on frequencies around ω b = 0 . .To study the spectral stability of the breathers, we introducea small perturbation ξ n ( t ) to a given solution u n, ( t ) of thediscretized dynamical equations as u n ( t ) = u n, ( t ) + ξ n ( t ) [8]. Then, defining the perturbation as ξ n = ξ x,n exp( ik y y ) in the quasi-1D breather stripe case yields ¨ ξ x,n + (cid:2) cos( u n ) + k y (cid:3) ξ x,n + 1 h ( ξ x,n +1 − ξ x,n + ξ x,n − ) = 0 . On the other hand, for the radial breather the correspondingdynamics for the perturbation may be written as ¨ ζ r,n + (cid:20) cos( u n ) + k θ n h (cid:21) ζ r,n + 1 h ( ζ r,n +1 − ζ r,n + ζ r,n − )+ 12 nh ( ζ r,n − − ζ r,n +1 ) = 0 , as ζ n = ζ r,n exp( ik θ θ ) . Due to the (temporal) periodicity ofthe solutions, Floquet analysis must be employed. In such acase, the stability properties are determined by the spectrumof the Floquet operator M (whose matrix representation is themonodromy), defined as: (cid:18) { ξ n ( T ) }{ ˙ ξ n ( T ) } (cid:19) = M (cid:18) { ξ n (0) }{ ˙ ξ n (0) } (cid:19) , (42)for the perturbation ξ n in the quasi-1D breather stripe caseand an identical equation for ζ n for the radial breather case.The N × N monodromy eigenvalues Λ ≡ exp( i Θ) aredubbed the FMs (Floquet multipliers), with Floquet exponents (FEs) Θ . A consequence of the fact that the Floquet opera-tor is real is that, if Λ is an FM, Λ ∗ is an FM too. Further,because of the symplecticity of the Floquet operator, / Λ isalso an FM. In other words, FMs always come in quadru-plets (Λ , Λ ∗ , / Λ , / Λ ∗ ) if the monodromy eigenvalues arecomplex, and in pairs (Λ , / Λ) if the eigenvalues are real.As a consequence, a necessary and sufficient condition for abreather to be linearly stable is that Θ must be real (i.e., thatthe corresponding FMs lie on the unit circle in the complexplane). Finally, note that FMs Λ are related to eigenvalues λ through relation Λ = exp( λT ) . [1] R. K. Dodd, J. C. Eilbeck, J. D. Gibbon, and H. C. Morris, Soli-tons and Nonlinear Wave Equations , Academic Press (London,1982).[2] Th. Dauxois and M. Peyrard,
Physics of Solitons , CambridgeUniversity Press (Cambridge, 2006).[3] J. Cuevas, P. G. Kevrekidis, and F. L. Williams (Eds.),
The sine-Gordon Model and its Applications: From Pendula and Joseph-son Junctions to Gravity and High Energy Physics , Springer-Verlag (Heidelberg, 2014).[4] Yu. S. Kivshar and B. A. Malomed, Rev. Mod. Phys. , 763(1989).[5] C. Rogers and W. K. Schief, B¨acklund and Darboux transfor-mations (Cambridge University Press: Cambridge, 2002).[6] M. J. Ablowitz, D. J. Kaup, A. C. Newell, and H. Segur, Phys.Rev. Lett. , 1262 (1973).[7] S. Flach and A. Gorbach, Discrete breathers – Advances in the-ory and applications , Phys. Rep. , 1 (2008).[8] S. Aubry, Phys. D , 201 (1997).[9] P. L. Christiansen and O. H. Olsen, Phys. Scr. , 531 (1979);Phys. Lett. A , 185 (1978); P. S. Lomdahl, O. H. Olsen, and P.L. Christiansen, Phys. Lett. A , 125 (1980); P. L. Christiansenand P. S. Lomdahl, Phys. D , 482 (1981).[10] J. Geicke, Phys. D , 197 (1982); Phys. Scr. , 431 (1984);Phys. Lett. A , 147 (1983).[11] I. L. Bogolubsky and V. G. Makhankov, JETP Lett. , 12(1976).[12] M. R. Samuelsen, Phys. Lett. A , 21 (1979).[13] J.-G. Caputo and M. P. Soerensen, Phys. Rev. E , 022915(2013).[14] P. G. Kevrekidis, I. Danaila, J.-G. Caputo, and R. Carretero-Gonz´alez, Phys. Rev. E , 052217 (2018).[15] E. G. Fedorov, A. V. Zhukov, R. Bouffanais, B. A. Malomed,H. Leblond, D. Mihalache, N. N. Rosanov, M. B. Belonenko, and T. F. George, Opt. Express , 27592 (2019).[16] E. M. Maslov, Physica D , 433 (1985).[17] Y. Kodama and M. J. Ablowitz, J. Math. Phys. , 928 (1980).[18] Yu. S. Kivshar and D. E. Pelinovsky, Phys. Rep. , 117(2000).[19] E. A. Kuznetsov and S. K. Turitsyn, Zh. Eksp. Teor. Fiz. ,119 (1988) [Sov. Phys. JETP , 1583 (1988)].[20] B. P. Anderson, P. C. Haljan, C. A. Regal, D. L. Feder, L. A.Collins, C. W. Clark, and E. A. Cornell, Phys. Rev. Lett. ,2926 (2001).[21] V. Tikhonenko, J. Christou, B. Luther-Davies, and Yu. S.Kivshar, Opt. Lett. , 1129 (1996).[22] S.-P. Gorza, B. Deconinck, Ph. Emplit, T. Trogdon, and M.Haelterman, Phys. Rev. Lett. , 094101 (2011).[23] J. Cuevas-Maraver and P. G. Kevrekidis, in A dynamical per-spective on the φ model. Past, present and future , P. G.Kevrekidis and J. Cuevas Maraver (eds.) (Springer, Cham,2019).[24] L. A. Cisneros-Ake, R. Carretero-Gonz´alez, and P. G.Kevrekidis, Phys. Rev. Research , 033043 (2019).[25] J. Yang, Nonlinear waves in integrable and nonintegrable sys-tems, Society for Industrial and Applied Mathematics (2010).[26] A. Jeffrey and T. Kawahara, Asymptotic Methods in NonlinearWave Theory , Pitman, London (1982).[27] E. A. Kuznetsov, A. M. Rubenchik, and V. E. Zakharov, Phys.Rep. , 103 (1986).[28] S.-P. Gorza, P. Kockaert, Ph. Emplit, and M. Haelterman, Phys.Rev. Lett. , 134101 (2009).[29] B. A. Malomed, Progr. Optics , 71 (2002).[30] A. M. Morgante, M. Johansson, S. Aubry, and G. Kopidakis, J.Phys. A: Math. Gen. , 4999 (2002).[31] G. Fibich, The Nonlinear Schr¨odinger Equation: Singular So-lutions and Optical Collapse (Springer: Heidelberg, 2015). [32] B. A. Malomed, Eur. Phys. J. Special Topics s225 , 2507 (2016).[33] K. D. Moll and A. L. Gaeta, in Nonlinear Optics: Materials,Fundamentals and Applications , A. Sawchuk (ed.) (Optical So-ciety of America, 2002).[34] Cheng-An Chen and Chen-Lung Hung, arXiv:1907.12550(2019).[35] H. Sakaguchi and B. A. Malomed, Phys. Rev. E , 012211(2020). [36] M.P. Calvo and J.M. Sanz-Serna, SIAM J. Sci. Comp. , 094101 (2016).[38] J. F. R. Archilla, R. S. MacKay, and J.L. Mar´ın, Phys. D ,406 (1999).[39] J. L. Mar´ın and S. Aubry, Nonlinearity9