Reflectionless propagation of Manakov solitons on a line:A model based on the concept of transparent boundary conditions
K.K. Sabirov, J.R. Yusupov, M.M. Aripov, M. Ehrhardt, D.U. Matrasulov
aa r X i v : . [ n li n . S I] F e b Reflectionless propagation of Manakov solitons on a line:A model based on the concept of transparent boundary conditions
K.K. Sabirov , , J.R. Yusupov , M.M. Aripov , M. Ehrhardt and D.U. Matrasulov Tashkent University of Information Technologies,108 Amir Temur Str., 100200, Tashkent Uzbekistan Yeoju Technical Institute in Tashkent, 156 Usman Nasyr Str., 100121, Tashkent, Uzbekistan National University of Uzbekistan, 2 Universitet Str., 100095, Tashkent, Uzbekistan Bergische Universit¨at Wuppertal, Gaußstrasse 20, D-42119 Wuppertal, Germany Turin Polytechnic University in Tashkent, 17 Niyazov Str., 100095, Tashkent, Uzbekistan Tashkent State Technical University, 2 Universitet Str., 100095, Tashkent, Uzbekistan
We consider the problem of absence of backscattering in the transport of Manakov solitons ona line. The concept of transparent boundary conditions is used for modeling the reflectionlesspropagation of Manakov vector solitons in a one-dimensional domain. Artificial boundary condi-tions that ensure the absence of backscattering are derived and their numerical implementation isdemonstrated.
I. INTRODUCTION
The Manakov system is an integrable system of cou-pled nonlinear Schr¨odinger equations (NLSE) which al-lows different soliton solutions. The main applicationof the Manakov system comes from nonlinear optics,where it describes optical vector solitons propagating inKerr media [1, 2]. Manakov-type vector solitons also ap-pear in optical fibers, in Bose-Einstein condensation, andin other areas of physics (see, e.g., Refs. [1, 2] for anoverview). So far, different aspects of the vector solitonsdescribed by the Manakov system have been studied [3–12]. In [13] an experimental realization of such solitonsin crystals was studied. The effect of small perturba-tions on the collision of vector solitons and its applicationto ultrafast soliton switching devices was investigated in[14]. The realization of logic gates and computationaloperations using Manakov vector solitons was discussedin [15]. The suppression of Manakov soliton interferencein optical fibers caused by the interaction of two vec-tor solitons in the Manakov equations that govern pulsetransmission in randomly birefringent fibers was studiedin [16]. Ref. [17] studies the rogue waves described by theManakov system with variable coefficients and externalpotential.In most cases, vector solitons are used as signal car-riers in optics and optoelectronic technologies. For ef-fective signal transmission in such devices and optimiza-tion of their functional properties, signal losses must beavoided or minimized by achieving a minimum of soli-ton backscattering, i.e., by propagating the solitons with-out reflections. The successful solution of such a prob-lem requires the construction of mathematical models de-scribing the reflectionless transport of solitons in a givenmedium. One of the effective mathematical tools for solv-ing the problem of reflectionless soliton propagation is theimposition of so-called transparent boundary conditions(TBCs) on a wave equation describing soliton transport.For special cases of the NLSE it is possible to formulatethe exact TBC in closed form, cf. [18]. For nonlinear wave equations, the concept of TBCs often relies on theso-called unified approach [19] that is based on a splittingprocedure of the linear and nonlinear part. However, ingeneral, for NLSE-type equations, it has turned out thatthe so-called potential approach [20] is the most tractableone. Here we extend this promising concept to the inte-grable Manakov system, which is a coupled system ofNLSE with vector soliton solutions.This paper is organized as follows. In the next sectionwe briefly recall soliton solutions and conserving quan-tities for the Manakov system on a line. Section IIIpresents the derivation of the transparent boundary con-ditions for the Manakov system. In Section IV we demon-strate our numerical implementation of such complicatedboundary conditions. Finally, Section V provides theconcluding remarks.
II. SOLITON SOLUTIONS OF THE MANAKOVSYSTEM
The
Manakov system can be written explicitly asi ∂ t Ψ + 12 ∂ x Ψ + ( | Ψ | + | Ψ | )Ψ = 0 , i ∂ t Ψ + 12 ∂ x Ψ + ( | Ψ | + | Ψ | )Ψ = 0 , (1)where (Ψ , Ψ ) = (cid:0) Ψ ( x, t ) , Ψ ( x, t ) (cid:1) , x ∈ R , t >
0. Itwas introduced first by Manakov [21] to describe statio-nary self-focusing electromagnetic waves in homogeneouswaveguide channels. The one-soliton solution of theManakov system can be written as(Ψ ∗ , Ψ ∗ ) = i (cid:16) c | c | (cid:17) η exp | η − ξ ) t − xξ | cosh[2 η ( x + x + 2 ξt )] , (2)where the initial position of the soliton is given by x =ln(2 η/ | c | ) / η and the unit vector (cid:0) c ≡ ( c , c )] (cid:1) deter-mines the polarization of the soliton. The parameters ξ and η denote the speed and amplitude of the soliton,respectively.Multi-soliton solutions of the Manakov system can beobtained using Hirota’s bilinearization method [22, 23].Eq. (1) approves two conserving quantities, such as thenorm determined as N = Z ∞−∞ (cid:0) | Ψ | + | Ψ | (cid:1) dx and the energy, which is given by, cf. [24] E = + ∞ Z −∞ X m =1 (cid:12)(cid:12)(cid:12)(cid:12) ∂ Ψ m ∂x (cid:12)(cid:12)(cid:12)(cid:12) − X m =1 | Ψ m | − | Ψ | | Ψ | ! dx. (3)In the following we will use these quantities for con-firming reflectionless transmission of Manakov solitonsthrough a given boundary. III. TRANSPARENT BOUNDARYCONDITIONS FOR THE MANAKOV SYSTEM
The reflection of nonlinear waves at the boundary of agiven domain is a practically important problem, the so-lution of which requires the use of an explicit solution of awave equation describing these waves. The mathematicaldescription of the absence of reflection at the boundary isa rather complicated task, since unlike quantum mechan-ics, no S-matrix theory exists for nonlinear waves. Oneof the effective approaches to solve such a problem canbe done in the framework of the concept of transparentboundary conditions. Such transparent boundary condi-tions for evolution equations can be constructed by cou-pling the solutions of the initial value boundary problems(IVBPs) in the interior and exterior domains [25–38].To construct transparent boundary conditions (TBCs)for a given PDE, one must first split the original waveequation into coupled equations determined in the inte-rior (Ω int ) and exterior (Ω ext ) domains. Then one appliesa Laplace transform in time to the exterior problems andobtains the solution of the ordinary differential equationsin the spatial variable x . Moreover, if one allows only”outgoing” waves by choosing the asymptotically decay-ing solution as x → ±∞ and matching the Dirichlet andNeumann values on the artificial boundaries of the inte-rior domain, one should apply (numerically) the inverseLaplace transform to complete the full derivation of theTBC [31].Here we will apply the above concept and procedurefor the derivation of TBCs for the Manakov systemEq. (1) and their numerical implementation at the ar-tificial boundary points x = 0, x = L . For this purpose,we use the so-called potential approach , which was previ-ously used to derive TBCs for the nonlinear Schr¨odingerequation [20, 39] and the sine-Gordon equation [40]. In[41–43] the TBC concept was used to develop transpar-ent quantum graphs model, which was later implementedto describe reflectionless transport of charge carriers in branched conducting polymers [44]. Within such an ap-proach, the Manakov system (1) is formally reduced to asystem of linear PDEs by introducing the following po-tential V ( x, t ) = | Ψ | + | Ψ | . This procedure allows usto rewrite Eq. (1) into the following linear form:i ∂ t Ψ + 12 ∂ x Ψ + V ( x, t )Ψ = 0 , i ∂ t Ψ + 12 ∂ x Ψ + V ( x, t )Ψ = 0 . (4)We also introduce the following new vector function as v ( x, t ) = e − i ν ( x,t ) Ψ( x, t ) , Ψ( x, t ) = (cid:18) Ψ ( x, t )Ψ ( x, t ) (cid:19) , v ( x, t ) = (cid:18) v ( x, t ) v ( x, t ) (cid:19) , (5)where ν ( x, t ) = Z t V ( x, τ ) dτ = Z t (cid:0) | Ψ ( x, τ ) | + | Ψ ( x, τ ) | (cid:1) dτ. (6)Taking here partial derivatives we obtain ∂ t Ψ = e i ν ( ∂ t + i V ) v,∂ x Ψ = e i ν (cid:0) ∂ x + 2i ∂ x ν · ∂ x + i ∂ x ν − ( ∂ x ν ) (cid:1) v. Inserting these latter equations in (4), we get L ( x, t, ∂ x , ∂ t ) v = i ∂ t v + 12 ∂ x v + A∂ x v + Bv = 0 , (7)where A = i ∂ x ν , B = (cid:0) i ∂ x ν − ( ∂ x ν ) (cid:1) . LinearizingEq. (7) using pseudo-differential operators we have L = (cid:16) √ ∂ x + iΛ − (cid:17)(cid:16) √ ∂ x + i Λ + (cid:17) = 12 ∂ x + i √ + + Λ − ) ∂ x + i √ ∂ x λ + ) − Λ − Λ + . (8)From Eqs. (7) and (8) we obtain the following system ofoperators i √ + + Λ − ) = A, i √ ∂ x λ + ) − Λ − Λ + = i ∂ t + B, (9)which yields the symbolic system of equationsi √ λ + + λ − ) = a, i √ ∂ x λ + − + ∞ X α =0 ( − i) α α ! ∂ ατ λ − ∂ αt λ + = − τ + b, (10)where the setting a = A , b = B is used. An asymptoticexpansion in the inhomogeneous symbols is defined as λ ± ∼ + ∞ X j =0 λ ± / − j/ . (11)Inserting the expansion (11) into Eq. (10) we can identifythe terms of order 1/2 in the first relation of the system(10): λ − / = − λ +1 / , λ +1 / = ±√− τ . (12)The Dirichlet-to-Neumann operator corresponds to thechoice λ +1 / = −√− τ . For the zero order terms we obtain λ − = − λ +0 − i √ a, i √ ∂ x λ +1 / − ( λ − λ +1 / + λ +0 λ − / ) = 0 . (13)From Eq. (13) we get λ +0 = − i √ a = √ ∂ x ν,λ − = − λ +0 − i √ a = √ ∂ x ν. (14)For the terms of order -1/2 we havei √ λ + − / + λ −− / ) = 0 , i √ ∂ x λ +0 − ( λ − / λ + − / + λ − λ +0 + λ −− / λ +1 / ) = b, (15)since ∂ αt λ ±− / = ∂ ατ λ ± = 0 , α ∈ N . From (15) we get λ ±− / = 0 . (16)Furthermore, one can obtain the next order terms as λ −− = − λ + − , λ + − = i √ τ ∂ x V. (17)Then the first order approximation reads1 √ ∂ x Ψ | x =0 − e − i π e i ν · ∂ / t ( e − i ν Ψ ) | x =0 = 0 , √ ∂ x Ψ | x =0 − e − i π e i ν · ∂ / t ( e − i ν Ψ ) | x =0 = 0 , (18)1 √ ∂ x Ψ | x = L + e − i π e i ν · ∂ / t ( e − i ν Ψ ) | x = L = 0 , √ ∂ x Ψ | x = L + e − i π e i ν · ∂ / t ( e − i ν Ψ ) | x = L = 0 . (19)The second order approximation is1 √ ∂ x Ψ | x =0 − e − i π e i ν · ∂ / t ( e − i ν Ψ ) | x =0 − i √ ∂ x V e i ν I t ( e − i ν Ψ ) | x =0 = 0 , √ ∂ x Ψ | x =0 − e − i π e i ν · ∂ / t ( e − i ν Ψ ) | x =0 − i √ ∂ x V e i ν I t ( e − i ν Ψ ) | x =0 = 0 , (20)1 √ ∂ x Ψ | x = L + e − i π e i ν · ∂ / t ( e − i ν Ψ ) | x = L + i √ ∂ x V e i ν I t ( e − i ν Ψ ) | x = L = 0 , √ ∂ x Ψ | x = L + e − i π e i ν · ∂ / t ( e − i ν Ψ ) | x = L + i √ ∂ x V e i ν I t ( e − i ν Ψ ) | x = L = 0 , (21)where I t ( f ) = R t f ( τ ) dτ . Unlike the standard the stan-dard Dirichlet, Neumann or Robin boundary conditions,the boundary conditions given by Eqs. (19) and (21), arequite complicated and can be implemented only numeri-cally. Therefore, we will provide in the next section theirnumerical implementation for the Manakov system (1). IV. DISCRETIZATION OF THE MANAKOVSYSTEM AND TRANSPARENT BOUNDARYCONDITIONS
For the numerical solution of the system (1) by impos-ing transparent boundary conditions given by Eqs. (19)and (21) one must use an effective discretization schemefor both, Eq. (1) and the boundary conditions. In thecase of transparent boundary conditions, the accuracyand stability of the numerical solution is very sensitiveto the choice of a discretization scheme. Here we presenta numerical scheme for Eq. (1) and a procedure for im-plementing the transparent boundary conditions.
A. Discretization of the equation
The numerical solution of coupled Schr¨odinger equa-tions is a well-studied problem and different high accu-racy numerical methods have been developed in the lit-erature so far (see, e.g. [45–49]).Here we use the explicit mid point rule [50], the so-called leap-frog finite difference method which is given asi Ψ n +11 ,j − Ψ n − ,j t + 12 D x Ψ n ,j + V nj Ψ n ,j = 0 , i Ψ n +12 ,j − Ψ n − ,j t + 12 D x Ψ n ,j + V nj Ψ n ,j = 0 , (22)with the standard second order difference quotient D x Ψ nj = 1∆ x (Ψ nj +1 − nj + Ψ nj − ) , and the discrete potential term V nj = | Ψ n ,j | + | Ψ n ,j | , where ∆ t and ∆ x are the time and space discretizationsteps, respectively. We note that there are other implicitor semi-implicit methods with higher accuracy [48–53].But the complexity of TBC forces to choose between ac-curacy and computational cost. And this method waschosen because of its simple implementation and low costper step.This leads to the following explicit finite differencescheme: U n +1 j = U n − j + i∆ tD x U nj + 2i∆ tV nj U nj , (23)where U nj = (cid:18) Ψ n ,j Ψ n ,j (cid:19) . Furthermore, we have to implement the transparentboundary conditions given by Eqs. (18)-(19) or Eqs. (20)-(21) in the above numerical methods.
B. Implementation of the TBC
The discretization of the TBC given by the Eqs. (18)-(19) and its subsequent implementation in the nume-rical scheme for the Eq. (1) and ensuring the stabilityof the whole numerical scheme is a rather complicatedtask. The presence of the fractional derivative also makesthe discretization scheme very complicated and unstable.Here we give an effective discretization scheme for theTBC, which can be implemented with high accuracy andstability in combination with the discretization schemefor Eq. (1). We state the scheme only for x = L by sayingthat for the left boundary (at x = 0) the implementationcan be done in the same way.The approximation of the fractional differential oper-ator is given by the numerical quadrature formula [20] ∂ / t f ( t n ) ≈ √ t n X k =0 β k f n − k , where { f n } n ∈ N is a sequence of complex values appro-ximating { f ( t n ) } n ∈ N and ( β k ) k ∈ N denotes the sequencesdefined by( β , β , β , β , β , β , . . . ) = (cid:18) , − , , − , · · , − · · , . . . (cid:19) . The function ν ( x, t ) given by (6) can be discretized usingthe trapezoidal rule as ν nj = ∆ t " n − X k =1 V kj + 12 V nj , for n ≥ , with ν j = 0 and ν j = ∆ t V j . Let us note, that the term V n was dropped, because the initial data is assumed tobe compactly supported and hence V n is zero. Then, wediscretize the function e i ν ( x,t ) from (5) as E nj = exp (i ν nj ) =exp (cid:16) i∆ t n − X k =1 V kj (cid:17) · exp (cid:16) i∆ t V nj (cid:17) = ^ E n − j · exp (cid:16) i∆ t V nj (cid:17) , where ^ E n − j = exp (cid:16) i∆ t n − P k =1 V kj (cid:17) .Thus, the TBC operator of the first order approxima-tion (19) at the right boundary j = J can be approxi-mated by the discrete convolutions(Λ n ) I = e − i π r t E nJ n X k =0 β k E n − kJ Ψ n − k ,J , (Λ n ) I = e − i π r t E nJ n X k =0 β k E n − kJ Ψ n − k ,J , where E nJ denotes the complex conjugate of E nJ .Then the values of the wave function at the rightboundary can be obtained by solving the system of non-linear equations with respect to (Ψ n ,J , Ψ n ,J ) ⊤ , given asΨ n ,J − Ψ n ,J − ∆ x + e − i π √ ∆ t h Ψ n ,J + ^ E n − J · exp (cid:16) i ∆ t (cid:0) | Ψ n ,J | + | Ψ n ,J | (cid:1)(cid:17) n X k =1 β k E n − kJ Ψ n − k ,J i = 0 , Ψ n ,J − Ψ n ,J − ∆ x + e − i π √ ∆ t h Ψ n ,J + ^ E n − J · exp (cid:16) i ∆ t (cid:0) | Ψ n ,J | + | Ψ n ,J | (cid:1)(cid:17) n X k =1 β k E n − kJ Ψ n − k ,J i = 0 . Using the same approach, we can proceed with the dis-cretization of the second order approximation. We recallthat V ( x, t ) = Ψ Ψ + Ψ Ψ and approximate ∂ x V ( x, t )at the right boundary x = L (i.e. j = J ) with dV nJ = 1∆ x (cid:16) | Ψ n ,J | − Ψ n ,J − Ψ n ,J − Ψ n ,J − Ψ n ,J +2 | Ψ n ,J | − Ψ n ,J − Ψ n ,J − Ψ n ,J − Ψ n ,J (cid:17) , where Ψ is the complex conjugate of Ψ. Then, again us-ing the trapezoidal method, we approximate the integralterm I t ( · ) in (21) with I nm,t = ∆ t n − X k =1 E kJ Ψ km,J + 12 E nJ Ψ nm,J ! , m = 1 , . Thus, the TBC operator of the second order approxima-tion (21) can be approximated as(Λ n ) II = (Λ n ) I + i √ dV nJ E nJ I n ,t , (Λ n ) II = (Λ n ) I + i √ dV nJ E nJ I n ,t . Again, the values of the wave function at the rightboundary can be obtained by solving the system of non-linear equations with respect to (Ψ n ,J , Ψ n ,J ) ⊤ , given asΨ n ,J − Ψ n ,J − ∆ x + e − i π √ ∆ t h Ψ n ,J + ^ E n − J · exp (cid:16) i ∆ t (cid:0) | Ψ n ,J | + | Ψ n ,J | (cid:1)(cid:17) n X k =1 β k E n − kJ Ψ n − k ,J i +i ∆ t x (cid:16) | Ψ n ,J | − Ψ n ,J − Ψ n ,J − Ψ n ,J − Ψ n ,J +2 | Ψ n ,J | − Ψ n ,J − Ψ n ,J − Ψ n ,J − Ψ n ,J (cid:17) · ^ E n − J · exp (cid:16) i ∆ t (cid:0) | Ψ n ,J | + | Ψ n ,J | (cid:1)(cid:17) · n − X k =1 E kJ Ψ k ,J + 12 E nJ Ψ n ,J ! = 0 , Ψ n ,J − Ψ n ,J − ∆ x + e − i π √ ∆ t h Ψ n ,J + ^ E n − J · exp (cid:16) i ∆ t (cid:0) | Ψ n ,J | + | Ψ n ,J | (cid:1)(cid:17) n X k =1 β k E n − kJ Ψ n − k ,J i +i ∆ t x (cid:16) | Ψ n ,J | − Ψ n ,J − Ψ n ,J − Ψ n ,J − Ψ n ,J +2 | Ψ n ,J | − Ψ n ,J − Ψ n ,J − Ψ n ,J − Ψ n ,J (cid:17) · ^ E n − J · exp (cid:16) i ∆ t (cid:0) | Ψ n ,J | + | Ψ n ,J | (cid:1)(cid:17) · n − X k =1 E kJ Ψ k ,J + 12 E nJ Ψ n ,J ! = 0 . V. NUMERICAL EXPERIMENT
We solve the system of coupled nonlinear Schr¨odingerequations, given by Eq. (1) on the finite interval [0 , L ]and impose the TBC at the right boundary (at x = L ).As initial condition we choose a single soliton from theexact solution given byΨ ( x,
0) = √ α sech (cid:2) √ α ( x − x ) (cid:3) exp (cid:2) i √ p ( x − x ) (cid:3) Ψ ( x,
0) = √ α sech (cid:2) √ α ( x − x ) (cid:3) exp (cid:2) i √ p ( x − x ) (cid:3) . In our experiments we selected the following systemparameters: L = 40, the parameters of the initial con-dition α = 1, p = 1, x = 20 and the discretization xt| Ψ |
020 4010 200 01
FIG. 1: Evolution of a right-travelling single soliton simu-lated with the finite difference scheme (23). The first orderapproximation of the TBC is imposed at the right boundary( x = 40). time E ne r g y FIG. 2: Time evolution of the soliton energy in the interior(computational) domain [0 , L ]. parameters ∆ x = 0 .
05, ∆ t = 0 . , L ], Eq. (3) can be rewritten as E = 12 L Z X m =1 (cid:12)(cid:12)(cid:12)(cid:12) ∂ Ψ m ∂x (cid:12)(cid:12)(cid:12)(cid:12) − (cid:16) | Ψ | + | Ψ | (cid:17) ! dx. (24)In our calculations we use the discrete analog of theenergy ( E ( t n ) = E n ) given by Eq. (24): E n = 14∆ x J − X j =1 h(cid:12)(cid:12) Ψ n ,j +1 − Ψ n ,j − (cid:12)(cid:12) + (cid:12)(cid:12) Ψ n ,j +1 − Ψ n ,j − (cid:12)(cid:12) − x (cid:16)(cid:12)(cid:12) Ψ n ,j (cid:12)(cid:12) + (cid:12)(cid:12) Ψ n ,j (cid:12)(cid:12) (cid:17) i . The time evolution of the soliton energy within thelimits of the calculation interval is shown in Fig. 2. It isclear that the energy disappears when the soliton crossesthe boundary, i.e. the absence of backscattering. time A b s o l u t e E rr o r × -3 FIG. 3: Plot of the absolute error ’Err’ in the L -norm. At the end, to show that the numerical solution is reli-able, we calculate the absolute error ’Err’ measured with the L -norm (discretized by the trapezoidal rule): || Err( t n ) || = ∆ x J − X j =1 (cid:16)(cid:12)(cid:12) ∆Ψ n ,j (cid:12)(cid:12) + (cid:12)(cid:12) ∆Ψ n ,j (cid:12)(cid:12) (cid:17) +∆ x (cid:16)(cid:12)(cid:12) ∆Ψ n , (cid:12)(cid:12) + (cid:12)(cid:12) ∆Ψ n , (cid:12)(cid:12) + (cid:12)(cid:12) ∆Ψ n ,J (cid:12)(cid:12) + (cid:12)(cid:12) ∆Ψ n ,J (cid:12)(cid:12) (cid:17) , where ∆Ψ nk,j = Ψ nk,j − Ψ exact k ( x j , t n ). The plot of thiserror over time is presented in Fig. 3. VI. CONCLUSIONS
In this work, the reflectionless transmission of Man-akov solitons described by the Manakov system on aline subject to so-called transparent boundary conditionshas been studied. These transparent boundary condi-tions (TBCs) are derived analytically and an effectivediscretization for the TBCs and its implementation inthe numerical method for the Manakov system are pre-sented. The absence of backscattering for Manakov soli-tons, when TBCs are imposed is demonstrated by a di-rect numerical experiment. The results of the work canbe used for modeling the tunable transport of Manakovsolitons in optical media and for optimization problemsof optical devices, where such solitons appear as signalcarriers. The above model can be extended for modelingthe reflectionless propagation of Manakov vector solitonsin higher dimensional and branched domains, which arenow the tasks under process. [1] Y. S. Kivshar and G. P. Agrawal,
Optical Solitons: FromFibers to Photonic Crystals (Academic Press, San Diego,2003).[2] G. P. Agrawal,
Applications of Nonlinear Fiber Optics (Academic Press, San Diego, 2001).[3] R. Radhakrishnan, M. Lakshmanan, and J. Hietarinta,Phys. Rev. E , 2213 (1997).[4] T. Kanna and M. Lakshmanan, Phys. Rev. Lett. , 5043(2001).[5] M. J. Ablowitz, B. Prinari, and A. D. Trubatch, Inv.Probl. , 1217 (2004).[6] R. Radhakrishnan and M. Lakshmanan, J. Phys. A:Math. Gen. , 2683 (1995).[7] A. P. Sheppard and Y. S. Kivshar, Phys. Rev. E , 4773(1997).[8] M. Vijayajayanthi, T. Kanna, and M. Lakshmanan,Phys. Rev. A , 013820 (2008).[9] B. F. Feng, J. Phys. A: Math. Theor. , 355203 (2014).[10] Y. Ohta, D. S. Wang, and J. Yang, Stud. Appl. Math. , 345 (2011).[11] P. G. Kevrekidis and D. J. Frantzeskakis, Rev. Phys. ,140 (2016).[12] D. J. Frantzeskakis, J. Phys. A: Math. Theor. , 213001(2010).[13] J. U. Kang, G. I. Stegeman, J. S. Aitchison, and N. Akhmediev, Phys. Rev. Lett., , 3699 (1996).[14] J. Yang, Phys. Rev. E, , 2393 (1999).[15] K. Steiglitz, Phys. Rev. E, , 016608 (2000).[16] J. Yang, Phys. Rev. E, , 036606 (2002).[17] W-P. Zhong, M. Beli, and B. Malomed, Phys. Rev. E, , 053201 (2015).[18] C. Zheng, J. Comput. Phys. , 036707(2011).[20] X. Antoine, Ch. Besse, and S. Descombes, SIAM J. Nu-mer. Anal., , 2272 (2006).[21] S. V. Manakov, Sov. Phys. JETP
248 (1974).[22] S. Gancsan and M. Lakshmanan, J. Phys. A: Math. andGen. L1143 (1987).[23] M. Lakshmanan, Int. J. Bifurcation Chaos, (4), 532 (2008).[25] A. Arnold and M. Ehrhardt, J. Comput. Phys., ,611 (1998).[26] M. Ehrhardt, VLSI Design, , 325 (1999).[27] M. Ehrhardt and A. Arnold, Riv. di Math. Univ. diParma, , 57 (2001).[28] M. Ehrhardt, Acta Acustica united with Acustica, ,711 (2002).[29] A. Arnold, M. Ehrhardt, and I. Sofronov, Commun.Math. Sci., , 501 (2003). [30] S. Jiang and L. Greengard, Comput. Math. Appl., ,955 (2004).[31] X. Antoine, A. Arnold, C. Besse, M. Ehrhardt, and A.Sch¨adle, Commun. Comput. Phys., , 729 (2008).[32] M. Ehrhardt, Appl. Numer. Math. , 660 (2008).[33] A. Zisowsky and M. Ehrhardt, Math. and Comput. Mod-ell., , 1264 (2008).[34] L. ˘Sumichrast and M. Ehrhardt, J. Electr. Engineering, , 301 (2009).[35] X. Antoine, C. Besse, and P. Klein, J. Comput. Phys., , 312 (2009).[36] M. Ehrhardt, Numer. Math.: Theor. Meth. Appl., ,295 (2010).[37] P. Klein, X. Antoine, C. Besse, and M. Ehrhardt, Com-mun. Comput. Phys., , 1280 (2011).[38] A. Arnold, M. Ehrhardt, M. Schulte, and I. Sofronov,Commun. Math. Sci., , 889 (2012).[39] J.R. Yusupov, K.K. Sabirov, M. Ehrhardt, and D.U. Ma-trasulov, Phys. Rev. E, , 032204 (2019).[40] K.K. Sabirov, J.R. Yusupov, M. Ehrhardt, and D.U. Ma-trasulov, arXiv:2011.13299[41] J.R. Yusupov, K.K. Sabirov, M. Ehrhardt, and D.U. Ma-trasulov, Phys. Lett. A, , 2382 (2019).[42] M.M. Aripov, K.K. Sabirov, and J.R. Yusupov, Nanosys-tems: physics, chemistry, mathematics , (5), 501 (2019).[43] J.R. Yusupov, K.K. Sabirov, Q.U. Asadov, M. Ehrhardt,and D.U. Matrasulov, Phys. Rev. E, , 110861 (2020).[45] M.S. Ismail and T.R. Taha, Math. Comp. Simul., (6),547 (2001).[46] M.S. Ismail and S.Z. Alamri, Int. J. Comp. Math., (3), 333 (2004).[47] M.S. Ismail and T.R. Taha, Math. Comp. Simul., (4-5), 302 (2007).[48] M.S. Ismail, Math. Comp. Simul., (1), 273 (2008).[49] Q.-B. Xu and Q.-S. Chang, Qs. Acta Math. Appl. Sin.Engl. Ser., , 205 (2010).[50] J. de Frutos and J.M. Sanz-Serna, J. Comp. Phys., (1), 160-168 (1992).[51] S. Zhou and X. Cheng, Math. Comp. Simul., (12),2362 (2010).[52] J. Chen and L.-M. Zhang, Acta Math. Appl. Sinica, Engl.Ser., (2), 435 (2017).[53] Y. Hu, H. Li, and Z. Jiang Appl. Numer. Math.,153