arXiv:2006.00958v1 [nlin.PS] 1 Jun 2020 Kuznetsov-Ma breather-like solutions in the Salerno model J. Sullivan,1, ∗E. G. Charalampidis,2, †J. Cuevas-Maraver,3, ‡P. G. Kevrekidis,1, 4, §and N. I. Karachalios5, ¶ 1Department of Mathematics and Statistics, University of Massachusetts Amherst, Amherst, MA 01003-4515, USA 2Mathematics Department, California Polytechnic State University, San Luis Obispo, CA 93407-0403, USA 3Grupo de F´ısica No Lineal, Departamento de F´ısica Aplicada I, Universidad de Sevilla. Escuela Polit´ecnica Superior, C/ Virgen de ´ Africa, 7, 41011-Sevilla, Spain Instituto de Matem´aticas de la Universidad de Sevilla (IMUS). Edificio Celestino Mutis. Avda. Reina Mercedes s/n, 41012-Sevilla, Spain 4Mathematical Institute, University of Oxford, Oxford, UK 5Department of Mathematics, University of the Aegean, Karlovassi, 83200 Samos, Greece (Dated: June 2, 2020) The Salerno model is a discrete variant of the celebrated nonlinear Schr¨odinger (NLS) equation interpolating between the discrete NLS (DNLS) equation and completely integrable Ablowitz-Ladik (AL) model by appropriately tuning the relevant homotopy parameter. Although the AL model possesses an explicit time-periodic solution known as the Kuznetsov-Ma (KM) breather, the existence of time-periodic solutions away from the integrable limit has not been studied as of yet. It is thus the purpose of this work to shed light on the existence and stability of time-periodic solutions of the Salerno model. In particular, we vary the homotopy parameter of the model by employing a pseudo-arclength continuation algorithm where time-periodic solutions are identified via fixed-point iterations. We show that the solutions transform into time-periodic patterns featuring small, yet non-decaying far-field oscillations. Remarkably, our numerical results support the existence of previously unknown time-periodic solutions even at the integrable case whose stability is explored by using Floquet theory. A continuation of these patterns towards the DNLS limit is also discussed. PACS numbers: I. INTRODUCTION The study of rogue wave patterns has been a focal point of recent research in dispersive nonlinear wave systems [1–7]. Whether they may arise “spontaneously out of nowhere and disappear without a trace” [8], or generated gradually through energy transfer in multiple soliton collision, they are studied intensely in a diverse range of settings [9–12]. Relevant studies have appeared from superfluid helium [13], to plasmas [14] and from nonlinear optics [15–21] to the, arguably, most natural venue of water waves [22–25]. In the discrete realm, there is far fewer studies. This is in a sense relatively natural to expect. Much of what we know about rogue waves is intimately connected to techniques stemming from integrable systems and the discrete setting of differentialdifference equations is no exception. The most prototypical model associated with integrability here is the so-called AblowitzLadik (AL) model [26, 27]. In this setting, rogue waves in the form of the prototypical Peregrine soliton [28], but also in the form of the Kuznetsov [29], Ma [30], and Akhmediev [31] breathers have been identified in the work of [32]. The Peregrine soliton is a solution that is localized in both space and time and it is a special member of a family extending between the spatially periodic Akhmediev solutions and the temporally-periodic Kuznetsov-Ma (KM) ones. Integrability has also provided higher order solutions [33]. However, beyond these findings there is very little. A question that is natural to ask is what of this structure persists in the more physically relevant (for settings such as waveguide arrays or Bose-Einstein condensates in optical lattices) discrete nonlinear Schr¨odinger (DNLS) [34] setting? A partial result has stemmed from the statistical analysis of [35], which concluded that there is a large propensity towards freak waves near the integrable limit. Over the past few years, we have attempted to address some of the relevant extensions of the understanding of rogue waves past the strict realm of integrable systems. On the one hand, some of the present authors have attempted to develop rogue-wave identifying methods that go beyond the integrable realm [36, 37] such as tracing these solutions as fixed points in space-time. On the other hand, a different subset of the present authors has developed techniques for understanding the stability of these states, considering the Floquet analysis of the KM waves and examining the Peregrine states as natural limiting states thereof [38]. ∗Email: jksulliva[email protected] †Email:
[email protected] ‡Email:
[email protected] §Email:
[email protected] ¶Email:
[email protected]
2 Dynamical studies of the evolution of generic initial data (including experimental ones) [39–41] are ongoing and such examples exist also in the case of discrete systems [42]. The purpose of the present work is multi-fold. We aim to use the so-called Salerno model [43] as a vehicle for going controllably beyond the integrable AL limit and towards the physical realm of the DNLS limit. In principle, this is possible since the model itself involves a homotopy parameter, identified as ghereafter (see, also Eq. (1)) which allows such a (continuous) deformation. The hope is that one can take this path “all the way”, arriving at physically relevant solutions of the latter limit. Indeed, we will show that in special cases an unprecedented possibility exists to carry through this program. However, generically our study will show that solutions in the form of KM solitons will exist; we seek specifically the latter because not only can we obtain the Peregrines as a special limit thereof, but also the stability is computationally tractable. We carry out the continuation of such states and find surprisingly convoluted bifurcation diagrams. We also observe that these bifurcation diagrams (over g) can be quite different for different breather frequencies and in special cases may even continue to the DNLS limit. Importantly also, we find that the looping bifurcation diagrams identified come back to “intersect” at a different point the integrable limit, presumably giving rise to previously unknown, to the best of our understanding, integrable system solutions. Hence, our work may be of interest both to the applications of the system (bearing the first example of such states “reaching” the DNLS limit), but also to the mathematical analysis of the integrable system (providing previously unknown solutions with oscillatory tails in the AL limit). Our presentation is structured as follows. In the next section, we offer the general framework of the model and the analytical results from the stability analysis of the plane wave solutions. Then in Section III, we examine our numerical findings for different frequencies of the periodic state ωand over the Salerno-model homotopic continuation strength g. Finally, in Section IV we summarize our findings and offer some possible directions for future study. II. MODEL AND THEORETICAL SETUP A. Model and breather solutions The model that this work focuses on is the so-called Salerno model [43] given by i˙ Ψn+C(Ψn+1 −2Ψn+ Ψn−1) + g(Ψn+1 + Ψn−1)|Ψn|2+ 2 (1 −g)|Ψn|2Ψn= 0,(1) where Ψn:= Ψn(t) : Z×R→Cis the complex wavefunction of the nth lattice site (n∈Z) and overdot denotes differentiation with respect to time. The parameter g∈[0,1] is the homotopy parameter whence the DNLS [34] and AL [26, 27] models are obtained from Eq. (1) for values of g= 0 and g= 1, respectively. If hstands for the distance between adjacent nodes, then C= 1/h2in Eq. (1). Recall that the model is used due to its versatility. On the one hand, we can benefit from the knowledge of analytical solutions for the Peregrine soliton or the KM waveform at the AL limit, yet we can also attempt to connect (through continuations in g) these findings to the physically applicable (non-integrable) DNLS limit. Upon inserting the separation of variables ansatz Ψn=ψne2iq2t(2) into Eq. (1), we arrive at i˙ ψn+C(ψn+1 −2ψn+ψn−1) + g(ψn+1 +ψn−1)|ψn|2+ 2 (1 −g)|ψn|2−q2ψn= 0,(3) where the parameter qfixes the background amplitude. At the AL limit, i.e., when g= 1 and C= 1, Eq. (3) possesses a time-periodic solution, the Kuznetsov-Ma breather [32] given explicitly by ψn(t) = qcos (ωt +iθ) + Gcosh (rn) cos (ωt) + Gcosh (rn),(4) with frequency ω(related to the period of the solution via T= 2π/ω), θ=−arcsinh (ω),r= arccosh ([2 + cosh (θ)] /3) and G=−ω/ √3 sinh (r). Away from the integrable AL limit, explicit analytical solutions are no longer available. Thus, we identify time-periodic solutions of period Tby considering a temporal discretization in terms of the Fourier series expansion: ψn(t) = ∞ X m=−∞ φn,meimωt,(5)
3 where φn,m are the Fourier coefficients. Upon plugging Eq. (5) into Eq. (3), we arrive at a set of algebraic equations: −(mω + 2q2)φn,m +C(φn+1,m −2φn,m +φn−1,m)−2q2 +gX m′X m′′ {φn,m′φ∗ n,m′′ (φn−1,m−m′+m′′ +φn+1,m−m′+m′′ ) +2(1 −g)φn,m′φ∗ n,m′′ φn−1,m−m′+m′′ }= 0. (6) This system is solved by a fixed-point iteration method, e.g., the Newton-Raphson method. To this aim, one must firstly truncate the infinite spatial lattice into a finite one to render the system tractable. Here we choose n∈[−N/2, N/2],N∈2Z+ supplemented with periodic boundary conditions ψ−N/2(t) = ψN/2(t),∀t. This way, the total number of nodes is N+ 1. On the other hand, one must truncate the infinite (temporal) Fourier series of Eq. (5). We use |m| ≤ 41 in our computations presented below. For given ω,Cand g, we identify time-periodic solutions with high accuracy by imposing a strict tolerance criterion on the norm of successive iterates which is within 10−11. It should be noted in passing that upon convergence of the Newton-Raphson method, we construct the time-periodic solution by means of Eq. (5) [38]. The resulting, by construction, timeperiodic solution can be used as an initial condition (at t= 0) for a time-stepping scheme to examine its dynamical evolution. Additionally, it can also be used in a Floquet analysis-based stability computation to assess the spectral stability of the solution as we now discuss. B. Modulational instability and Floquet analysis The modulational instability (MI) of the asymptotic state (constant background) Ψn=qas |n| → ∞ of Eq. (1) (see, the seminal works of [44, 45] as well as the recent work of [46] on the subject) is investigated by considering the plane wave solution Ψn=qei(kn−ωt), k = 2πM √C N!, M ∈Z,(7) with frequency ωand wavenumber k. If we insert Eq. (7) to Eq. (1) we obtain the following dispersion relation: ω= 4 C+gq2sin2k 2−2q2.(8) We can now explore the linear stability of the plane wave solution of Eq. (7) by introducing the ansatz e Ψn=hq+ε(a+ib)ei(Qn−ωt)iei(kn−ωt), ε ≪1, Q = 2πM′ √C N!, M′∈Z,(9) where both aand bare time-dependent, real-valued functions, and Qand Ωcorrespond to the wavenumber and frequency of the perturbation, respectively. Upon plugging Eq. (9) into Eq. (1), we obtain at order O(ε)the MI dispersion relation given by Ω2= 2 gq2+C[cos (Q+k)−cos (k)] 2gq2+Ccos (Q+k) + 2 gq2−Ccos (k) + 4 (1 −g)q2.(10) If this condition is satisfied, yielding real frequencies Ωfor a given perturbation wavenumber Q, then the relevant wavenumber is stable. On the other hand, the existence of Q’s associated with complex Ωleads to dynamical instability of the background. The stability of time-periodic solutions with period Tidentified via fixed-point iterations (see, Sec. III) and denoted as ψ0 nis examined by considering the perturbation ansatz e ψn=ψ0 n+ε ξn(t), ε ≪1,(11) where ξn∈Cis the perturbation imposed at the nth site of the lattice. Then, we insert Eq. (11) into Eq. (3) and obtain at order O(ε)the governing equation for the perturbation ξn: i˙ ξn=−C(ξn+1 −2ξn+ξn−1)−2gψ0 n+1 +ψ0 n−1Re ξn(ψ0 n)∗−g(ξn+1 +ξn−1)|ψ0 n|2 −2 (1 −g)2|ψ0 n|2ξn+ξ∗ n(ψ0 n)2+ 2q2ξn.(12) Then, the eigenvalues λof the so-called monodromy matrix Mstemming from: Re(ξn(T)) Im(ξn(T)) =MRe(ξn(0)) Im(ξn(0)) .(13)
4 -1 0 1 2 -1 -0.5 0 0.5 1 5 10 15 20 0 0.2 0.4 0.6 0.8 1 FIG. 1: (Color online) The left panel presents the spectrum of a(n unstable) time-periodic solution with ω= 8 (filled black circles). Note that the theoretically predicted unstable modes via Eq. (10) are shown with red open circles therein. The right panel demonstrates a comparison between the Floquet exponents Re(iΩ) of time-periodic solutions (solid black lines) and the ones of the MI analysis (dashed-dotted red lines) of Eq. (10) over ω. Note that both panels correspond to the case with g= 1 (as well as C= 1 and q= 1/√2). determine the stability trait of of a time-periodic solution ψ0 n. Those eigenvalues are the so-called Floquet multipliers. In particular, as the system is symplectic and Hamiltonian, a solution is deemed neutrally stable if all the Floquet multipliers λ=λr+iλiof Mlie on the unit circle. If |λ|>1, then two types of instabilities can arise. If the multipliers arise in real pairs that are away from the unit circle, the instability is considered exponential due to the exponential growth of the associated perturbations. On the contrary, if the multipliers arise in complex quartets (inside and outside the unit circle), then the instability is deemed oscillatory. We conclude this section by mentioning in passing that the eigenfrequency of the perturbation Ωconnects with the Floquet multipliers λvia λ=eiΩT.(14) III. NUMERICAL RESULTS The availability of the exact solution (4) for the case with g= 1, i.e., the AL limit allows us to not only benchmark our numerical methods but also to compute the Floquet multipliers directly from Eq. (13). Hereafter, we consider a lattice with N= 100 sites and set C= 1 and q= 1/√2in Eqs. (6) and (12) (all of our numerical results discussed in this section were obtained for these choices). At first, we compute the Floquet multipliers from Eq. (13) (for g= 1) by using the initial-valueproblem (IVP) solver DOP853 [47] with (relative and absolute) tolerances 10−11. Although a non-stiff and explicit IVP solver, DOP853 can perform a time step-size adaptation (to satisfy the user-specified tolerance criteria) for stiff regions by reducing the time step-size. Indicatively, Fig. 1 summarizes our findings on the Floquet multipliers for the exact solution (4). The left panel of the figure presents the Floquet multipliers of the exact time-periodic solution with ω= 8 (filled black circles) together with the theoretically predicted unstable modes from Eq. (10) via Eq. (14) (open red circles). As far as the KM breather is concerned, the presence of Floquet multipliers with |λ|>1render the solution unstable. In addition, it is evident from this panel that a subset of the unstable modes of the KM breather coincides with those of the modulationally unstable background. However, there exist other unstable eigenmodes that deviate from the latter due to the presence of the localized solution perturbing the modes of the background. This is natural to expect given the modulational instability of the background (analyzed in the previous section) on top of which the periodic solution “lives”. The right panel of Fig. 1 complements our analysis of the exact KM breather and demonstrates the dependence of its unstable Floquet exponents Re(iΩ) on ω(solid black line) together with the respective results of the MI analysis (dashed-dotted red line). It can be discerned from this panel that the Floquet exponents of the KM breather approach the asymptotic values theoretically predicted by Eq. (10). We now wish to explore the existence and stability of time-periodic solutions over the relevant parameter space by going beyond the well-defined and analytically tractable case of the integrable limit of g= 1. We thus vary the parameter g(which is used as our bifurcation parameter) by employing a pseudo-arclength continuation algorithm [48, 49], and consider different fixed values of the period T. It should be noted that the pseudo-arclength continuation algorithm is capable of passing through turning points (where the Jacobian of the system of equations is singular and thus non-invertible) and tracing (connected) branches of solutions. The numerical results reported in this work were obtained by using a fixed and relatively small arclength step-size of 10−2in order to prevent the predictor-corrector step from converging to a solution of a nearby branch. As per the direct numerical simulations reported below, we use again the DOP853 [47] method for advancing Eq. (1) forward in time (and with the same tolerance criteria as before). Figure 2 summarizes our results on the existence of of time-periodic solutions to the Salerno model [cf. Eq. (1)]. In particular, we demonstrate the dependence of max (|ψ|)on gfor ω= 5,ω= 25, and ω= 10 in the left, middle, and right panels,
5 0.98 0.99 1 1.01 6 8 10 12 0.99 1 1.01 15 20 25 30 35 0 0.2 0.4 0.6 0.8 1 5 10 15 0.9 0.95 1 1.05 10 12 14 16 FIG. 2: (Color online) The l∞-norm of |ψn|as a function of gcorresponding to time-periodic solutions for values of ω= 5,ω= 25 and ω= 10, respectively. The labels (a)-(e) in the left and right panels are associated with Figs. 3 and 4 that follow. respectively. Note that the labels (a)-(e) in the left and right panels are connected to Figs. 3 and 4, respectively. Based on the panels of Fig. 2, a cascade of turning points is clearly evident although all of them suggest an intriguing finding that we discuss now. We consider first the left panel of the figure corresponding to ω= 5. Starting from g= 1, the time-periodic solution departs from the integrable limit (i.e., the AL limit), heading to smaller values of guntil it reaches a turning point at g≈0.981 upon which it comes back to g= 1, i.e., the AL limit, and then follows a “snake” pattern with several crossings of the AL limit. It should be noted that we stopped our continuation algorithm at g= 1 (the terminal solution is labeled by (e) therein). Based on our additional numerical investigations, this is en route to multiple additional crossings of g= 1 (results are not shown). The labels of the left panel are connected with the prototypical configurations for this case in Fig. 3. In particular, its left column presents the spatial distribution of the densities of the time-periodic solutions, i.e., |ψn|2for values of g= 1 (panels (a), (c)-(e)), and g= 0.982 (panel (b)), respectively. The insets therein correspond to the associated Floquet multipliers which themselves suggest that all computed solutions are highly unstable with a dominant unstable mode of the order of λr∼ O(1). Again, this can be naturally expected on the basis of the instability of the background. The striking feature of the solutions presented in panels (b)-(e) is the formation of an oscillatory background (of small amplitude) in contrast with the KM breather of panel (a) where the localized wave sits atop a constant background. Such profiles featuring small in-amplitude wave trains are strongly reminiscent of nanoptera, and to the best of our knowledge are first reported for the Salerno model in the present work (see, also the recent works of [50, 51] for Toda lattices and the DNLS with saturation). The other striking feature of our findings is that the KM breather (4) of the AL model (g= 1) seems not to be the only solution at that limit, as this has already been evident in panels (b)-(e) of Fig. 3. Although the investigation of those extra solutions by using integrable systems techniques is of fundamental importance, it is beyond the scope of our present work. We summarize our results in this case by presenting the respective spatio-temporal dynamics in the right panel of Fig. 3. In particular, contour plots of the spatio-temporal evolution of the density of time-periodic solutions are shown for one period (T= 2π/5in this case). Though the solutions remain robust for one period, thus validating our numerical approach for identifying them (via fixed-point iterations), we observed the emergence of the instability which happens at later times (results not shown) leading to a breakdown of the localized waveform. Similarly, the middle panel of Fig. 2 corresponds to the case of ω= 25. It presents a type of spiral structure around the KM solution that eventually leads to a progressively more dense loop structure that keeps repeating itself in the course of the pseudoarclength continuation. See also the inset therein which suggests that complex structures may arise also at a fine scale within the bifurcation diagram. The identified periodic states bear similar features as before, including progressively more pronounced (along the relevant branches) tails in the relevant waveforms. Furthermore, in line with earlier (MI-based and computational) findings and general expectations, the solutions are generically found to be unstable, due to the instability of their respective background. It is perhaps important to highlight two pervading features regarding the nature of our results (both the ones above, as well as the ones that will follow). The first one is that to the best of our knowledge the results above (and below) constitute the first definitive identification of rogue-like patterns in discrete (nonlinear dynamical lattice) systems beyond the extremely important, yet practically limited realm of integrable systems. Secondly, our findings are not solely of interest to non-integrable dispersive system practitioners, but they are also a motivation for further integrable system investigations. The panel (b) of Fig. 3, for example, suggests the existence of KM-type solutions on top of a stationary nanopteronic background which is time-independent in the spirit of recent works of [52, 53] in corresponding continuum limit problems. These are states that we believe are eminently relevant to explore in an analytical form within the framework of integrable systems (although this is outside the scope of the present study). We conclude this section of our findings by going through the right panel of Fig. 2 and the associated results shown in Fig. 4 corresponding to the case of ω= 25. Based on the former, we observe again a cascade of crossings happening at the AL limit (see, the inset therein), however here the bifurcation curve starts heading towards smaller values of guntil it reaches g= 0,
6 -20 -10 0 10 20 0 20 40 60 (a) 0 -1 0 1 -1 0 1 0123 -1 0 1 -20 -10 0 10 20 0 0.5 1 0 10 20 30 40 50 -20 -10 0 10 20 0 10 20 30 40 50 (b) -1 0 1 -1 0 1 -10 -5 0 -1 0 1 0 0.3 0.6 0.9 -20 -10 0 10 20 0 0.5 1 0 10 20 30 40 -20 -10 0 10 20 0 10 20 30 40 50 (c) -1 0 1 -1 0 1 -1 0 1 2 3 -1 0 1 0 0.3 0.6 0.9 -20 -10 0 10 20 0 0.5 1 0 10 20 30 40 -20 -10 0 10 20 0 30 60 90 (d) -1 0 1 -1 0 1 -1 0 1 2 3 -1 0 1 0 1 2 -20 -10 0 10 20 0 0.5 1 0 20 40 60 80 -20 -10 0 10 20 0 40 80 120 (e) -1 0 1 -1 0 1 -1 0 1 2 3 -1 0 1 -30 -25 -20 -15 0 0.3 0.6 0.9 -20 -10 0 10 20 0 0.5 1 0 40 80 120 FIG. 3: (Color online) Numerical results on time-periodic solutions of the Salerno model (3) for ω= 5 associated with the bifurcation curve of the left panel of Fig. 2 (see, the labels therein). The left column corresponds to the spatial distribution of |ψn|2for values of (a), (c)-(e) g= 1 and (b) g= 0.982, respectively. The insets shown on the left and right of each panel provide close-ups of the profiles as well as Floquet spectra. The right column presents the spatio-temporal evolution of |ψn|2of the profiles shown in the left column for one period (T= 2π/5).
7 -20 -10 0 10 20 0 50 100 150 (a) -1 0 1 2 -1 0 1 -20 -10 0 10 20 0 0.2 0.4 0.6 0 50 100 150 -20 -10 0 10 20 0 20 40 60 80 (b) -40-30-20-10 0 0.3 0.6 0.9 -1 0 1 -1 0 1 -20 -10 0 10 20 0 0.2 0.4 0.6 0 20 40 60 80 -20 -10 0 10 20 0 90 180 270 (c) -1 0 1 -1 0 1 -50 -25 0 0 0.4 0.8 1.2 -20 -10 0 10 20 0 0.2 0.4 0.6 0 90 180 270 -20 -10 0 10 20 0 20 40 60 80 (d) -1 0 1 -1 0 1 -50-35-20 -5 0 0.4 0.8 1.2 -20 -10 0 10 20 0 0.2 0.4 0.6 0 35 70 -20 -10 0 10 20 0 2 4 6 8(e) -1 0 1 2 -1 0 1 -20 -10 0 10 20 0 2 4 60 2 4 6 FIG. 4: (Color online) Same as Fig. 3 but for the case of ω= 10 associated with the right panel of Fig. 2. Again, the left and right columns correspond to the spatial distribution of |ψn|2and its spatio-temporal evolution, respectively. The integration was performed over one period for (a)-(d) and over 10 periods for panel (e), respectively. The insets shown on the left and right of each panel provide close-ups of the profiles as well as Floquet spectra. Panels (a)-(c) correspond to g= 1 whereas the ones of (d) and (e) to g= 0.9and g= 0 (DNLS case), respectively.
8 i.e., the DNLS limit. Notice how this appears to be an intermediate case between the less ordered and more expanded diagrams of small ωand the more ordered and confined diagram of ω= 25. The profiles shown in panels (a)-(c) in Fig. 4 correspond to values of g= 1 whereas the ones of (d) and (e) to g= 0.9and g= 0, respectively. Again, a common finding is that all time-periodic solutions identified in this work feature an oscillatory background as soon as we depart from the integrable limit (see, the left insets in the left column of Fig. 4) as well as they are unstable (see, the right insets therein). This effect is even more pronounced for the case with g= 0 as is shown in panel (e) of the figure. To the best of our knowledge, such time-periodic solutions of DNLS on a background have not appeared in the literature so far; indeed we are only aware of such effectively quasi-periodic solutions of the model on top of a vanishing background per the work of [54]. This solution is also unstable as is evident in the Floquet spectrum shown in the inset of the panel, although it remains robust over many periods of time integration. Indeed, the right panel therein demonstrates the spatio-temporal evolution of the density again over 10 periods where |ψ0(t)|2 breathes over time. It should also be noted that the oscillatory background is rather visible in this case (as well as the one of panels (b) and (d)) but again it remains steady over the time evolution, i.e., corresponding to a definitive nanopteronic state in this setting. IV. CONCLUSIONS AND FUTURE CHALLENGES In this work, we made an attempt to explore the existence, stability and dynamics of time-periodic solutions of the Salerno model. This was with a three-fold scope in mind: firstly, to establish that relevant solutions such as the KM are not unique or particular to the integrable limit, but can be continued to generic non-integrable values of the homotopic parameter g. Secondly, we wished to explore whether additional intriguing solutions could arise in the integrable model, a feature that was brought forth from our pseudo-arclength continuation results for periodic orbits. Lastly, we intended to examine whether some of the relevant solutions could be continued to the DNLS limit of g= 0; here, we found that for suitable choice of the breather frequency indeed that was also possible. Upon employing fixed-point methods, we identified the pertinent waveforms using Newton’s method (for periodic orbits) and their stability was inferred by performing a Floquet analysis. The use of pseudo-arclength continuation allowed us to perform a parametric continuation over the homotopy parameter gand this proved to be crucial in unraveling the complexity of the possible solutions in the Salerno model. Additionally for the integrable limit of g= 1 we identified multiple time-periodic solutions that sit atop of an oscillatory background being strongly reminiscent of nanoptera. Another striking finding of our work was the time-periodic solution which was identified at the DNLS limit, i.e., g= 0. To the best of our knowledge, such a waveform for the DNLS (on top of a non-vanishing background) was not previously reported. Based on the above findings and computational techniques that we have developed in this work, there are clearly many directions for future studies. At the AL limit, a potential analysis of the Zhakarov-Shabat problem for the nanopteronic solutions reported in Section III will be of paramount importance in order to derive them in possibly closed form (analogously to their identification in continuum problems [52, 53]). On the other hand, it is worth investigating the g= 0 solution reported in the present work, and in particular to study the configuration space of solutions as a function of ωin this case. This will pave the way towards potentially identifying Peregrine-like entities for the DNLS as time-periodic solutions at the limit T→ ∞(ω≫1). Another important path to explore is the continuation of the present time-periodic solutions as the distance between adjacent sites decreases, thus approaching the continuum limit when C= 1/h2increases. This way, we will be able to connect our findings with the results for the (focusing) NLS. A subset of the above computational studies can be carried out quite efficiently with the use of state-of-the-art bifurcation packages such as AUTO [55] and COCO [56] which additionally allow branch switching, among other features. Such directions are presently under consideration and will be reported in future publications. Acknowledgments J.C.-M. was supported by MAT201679866-R project (AEI/FEDER,UE). PGK acknowledges support from the U.S. National Science Foundation under Grants no. PHY-1602994 and DMS-1809074 (PGK). He also also acknowledges support from the Leverhulme Trust via a Visiting Fellowship and the Mathematical Institute of the University of Oxford for its hospitality during part of this work. [1] M. Onorato, S. Residori, U. Bortolozzo, A. Montina, and F. T. Arecchi, Phys. Rep. 528, 47 (2013). [2] P. T. S. DeVore, D. R. Solli, D. Borlaug, C. Ropers, and B. Jalali, J. Opt. 15, 064001 (2013). [3] Z. Yan, J. Phys. Conf. Ser. 400, 012084 (2012). [4] N. Akhmediev et al., J. Opt. 18, 063001 (2016). [5] S. Chen, F. Baronio, J. M. Soto-Crespo, P. Grelu, and D. Mihalache, J. Phys. A: Math. Theor. 50, 463001 (2017). [6] D. Mihalache, Rom. Rep. Phys. 69, 403 (2017).
9 [7] B. A. Malomed and D. Mihalache, Rom. J. Phys. 64, 106 (2019). [8] N. Akhmediev, A. Ankiewicz, and M. Taki, Phys. Lett. A 373, 675 (2009). [9] E. Pelinovsky and C. Kharif (eds.), Extreme Ocean Waves, Springer-Verlag (New York, 2008). [10] C. Kharif, E. Pelinovsky, and A. Slunyaev, Rogue Waves in the Ocean, Springer-Verlag (New York, 2009). [11] A. R. Osborne, Nonlinear Ocean Waves and the Inverse Scattering Transform, Academic Press (Amsterdam, 2010). [12] M. Onorato, S. Residori, and F. Baronio, Rogue and Shock Waves in Nonlinear Dispersive Media, Springer-Verlag (Heidelberg, 2016). [13] A. N. Ganshin, V. B. Efimov, G. V. Kolmakov, L. P. Mezhov-Deglin, and P. V. E. McClintock, Phys. Rev. Lett. 101, 065303 (2008). [14] H. Bailung, S. K. Sharma, and Y. Nakamura, Phys. Rev. Lett. 107, 255005 (2011). [15] D. R. Solli, C. Ropers, P. Koonath, and B. Jalali, Nature 450, 1054 (2007). [16] B. Kibler, J. Fatome, C. Finot, G. Millot, F. Dias, G. Genty, N. Akhmediev, and J. M. Dudley, Nature Phys. 6, 790 (2010) . [17] B. Kibler, J. Fatome, C. Finot, G. Millot, G. Genty, B. Wetzel, N. Akhmediev, F. Dias, and J. M. Dudley, Sci. Rep. 2, 463 (2012). [18] J. M. Dudley, F. Dias, M. Erkintalo, and G. Genty, Nat. Photon. 8, 755 (2014). [19] B. Frisquet, B. Kibler, P. Morin, F. Baronio, M. Conforti, G. Millon, and S. Wabnitz, Sci. Rep. 6, 20785 (2016). [20] C. Lecaplain, Ph. Grelu, J. M. Soto-Crespo, and N. Akhmediev, Phys. Rev. Lett. 108, 233901 (2012). [21] G. Genty, C. M. de Sterke, O. Bang, F. Dias, N. Akhmediev, and J. M. Dudley, Phys. Lett. A 374, 989 (2010). [22] A. Chabchoub, N. P. Hoffmann, and N. Akhmediev, Phys. Rev. Lett. 106, 204502 (2011). [23] A. Chabchoub, N. Hoffmann, M. Onorato, and N. Akhmediev, Phys. Rev. X 2, 011015 (2012). [24] A. Chabchoub and M. Fink, Phys. Rev. Lett. 112, 124101 (2014). [25] J.N. Steer, A.G.L. Borthwick, M. Onorato, A. Chabchoub, and T.S. van den Bremer, Phys. Rev. Lett. 123, 184501 (2019). [26] M.J. Ablowitz and J.F. Ladik, J. Math. Phys. 16, 598 (1975). [27] M.J. Ablowitz and J.F. Ladik, J. Math. Phys. 17, 1011 (1976). [28] D. H. Peregrine, J. Austral. Math. Soc. B 25, 16 (1983). [29] E. A. Kuznetsov, Sov. Phys.-Dokl. 22, 507 (1977). [30] Ya. C. Ma, Stud. Appl. Math. 60, 43 (1979). [31] N. N. Akhmediev, V. M. Eleonskii, and N. E. Kulagin, Theor. Math. Phys. 72, 809 (1987). [32] A. Ankiewicz, N. Akhmediev, and J. M. Soto-Crespo, Phys. Rev. E 82, 026602 (2010). [33] Y. Ohta and J. Yang, J. Phys. A 47, 255201 (2014). [34] P. G. Kevrekidis, The discrete nonlinear Schr¨odinger equation: Mathematical analysis, Numerical Computations and Physical Perspectives, (Springer-Verlag, Heidelberg, 2009). [35] A. Maluckov, Lj. Hadzievski, N. Lazarides, and G. P. Tsironis, Phys. Rev. E 79, 025601(R) (2009). [36] C.B. Ward, P.G. Kevrekidis, and N. Whitaker, Phys. Lett. A 383, 2584 (2019). [37] C.B. Ward, P.G. Kevrekidis, T.P. Horikis, and D.J. Frantzeskakis, Phys. Rev. Research 2, 013351 (2020). [38] J. Cuevas-Maraver, P. G. Kevrekidis, D. J. Frantzeskakis, N. I. Karachalios, M. Haragus, and G. James, Phys. Rev. E 96, 012202 (2017). [39] M. Bertola and A. Tovbis, Comm. Pure Appl. Math. 66, 678 (2013). [40] A. Tikan, C. Billet, G. El, A. Tovbis, M. Bertola, T. Sylvestre, F. Gustave, S. Randoux, G. Genty, P. Suret, and J.M. Dudley Phys. Rev. Lett. 119, 033901 (2017). [41] E.G. Charalampidis, J. Cuevas-Maraver, P.G. Kevrekidis, and D.J. Frantzeskakis, Rom. Rep. Phys. 70, 504 (2018). [42] C. Hoffmann, E.G. Charalampidis, D.J. Frantzeskakis, and P.G. Kevrekidis, Phys. Lett. A 382, 3064 (2018). [43] M. Salerno, Phys. Rev. A 46, 6856 (1992). [44] Y.S. Kivshar and M. Peyrard, Phys. Rev. A 46, 3198 (1992). [45] F.Kh. Abdullaev, A. Bouketir, A. Messikh, and B.A. Umarov, Physica D 232, 54 (2007). [46] F.II Ndzana and A. Mohamadou, Chaos 27, 073118 (2017). [47] E. Hairer, S.P. Nørsett and G. Wanner, Solving ordinary differential equations I (Springer-Verlag, Berlin, 1993). [48] A.H. Nayfeh and B. Balachandran, Applied Nonlinear Dynamics: Analytical, Computational and Experimental Methods (Wiley Series in Nonlinear Science, 1995). [49] E. Doedel, H.B. Keller and J.P. Kern´evez, Internat. J. Bifur. and Chaos, 01, 493-520 (1991). [50] C.J. Lustri and M.A. Porter, SIAM J. Appl. Dyn. Sys. 17, 1182 (2018). [51] G.L. Afimov, A.S. Korobeinikov, C.J. Lustri, and D.E. Pelinovsky, Nonlinearity 32, 3445 (2019). [52] J. Chen, D.E. Pelinovsky, R.E. White, arXiv:1905.11638. [53] J. Chen, D.E. Pelinovsky, J. Nonlinear Sci. 29, 2797 (2019). [54] M. Johansson and S. Aubry, Nonlinearity 10, 1151 (1997). [55] E. Doedel. AUTO. http://indy.cs.concordia.ca/auto/ [56] H. Dankowicz and F. Schidler. COCO. https://sourceforge.net/projects/cocotools