Physica D 175 (2003) 127–138
Cylindrical solitary pulses in a two-dimensional stabilized Kuramoto–Sivashinsky system Bao-Feng Feng a,∗ , Boris A. Malomed b , Takuji Kawahara c a
Department of Mathematics, The University of Kansas, Lawrence, KS 66045, USA Department of Interdisciplinary Studies, Faculty of Engineering, Tel Aviv University, Tel Aviv 69978, Israel Department of Aeronautics and Astronautics, Graduate School of Engineering, Kyoto University, Sakyo-ku, Kyoto 606-8501, Japan b
c
Received 27 December 2001; received in revised form 27 September 2002; accepted 4 October 2002 Communicated by Y. Kuramoto
Abstract By linearly coupling a generalized Zakharov–Kuznetsov equation (alias the two-dimensional (2D) Benney equation) to an extra linear dissipative equation, a 2D extension of a recently proposed stabilized Kuramoto–Sivashinsky system is developed. The model applies to the description of surface waves on 2D liquid layers in various physical settings. The extra equation provides for the stability of the zero state in the system, thus paving a way to the existence of stable 2D localized solitary pulses (SPs). A perturbation theory, based on a family of cylindrical solitons existing in the conservative counterpart of the system, is developed by treating dissipation and gain in the model as small perturbations. It is shown that the system may select two steady-state solitons from the continuous family provided by the conservative counterpart, of which the one with larger amplitude is expected to be stable. Numerical simulations support the analytical predictions quite well. Additionally, it is found that a shallow quasi-one-dimensional (1D) trough is attached to the stable SP if the integration domain is not very large, and an explanation to this feature is proposed. Stable double-humped bound states of two pulses are found too. © 2002 Elsevier Science B.V. All rights reserved. PACS: 05.45.Yv; 46.15.Ff Keywords: 2D solitary pulse; Stabilized Kuramoto–Sivashinsky equation; Zakharov–Kuznetsov equation
1. Introduction The Korteweg–de Vries (KdV) equation and one of its two-dimensional (2D) generalizations, the Kadomtsev– Petviashvili (KP) equation, are classical paradigms of integrable nonlinear evolution equations, which played an important role in the development of the soliton theory [1]. They arise in many physical contexts, such as surface water waves and ion-acoustic waves in plasmas [2,3], where dispersion and nonlinearity dominate, while dissipative effects are small enough to be neglected in the lowest-order approximation. However, in many cases dissipation must be taken into account, at least, as a perturbation. In this case, an active (amplifying) ingredient is necessary in the system, in order to provide for a loss-compensating gain. In plasma physics and hydrodynamics, the gain ∗
Corresponding author. E-mail address:
[email protected] (B.-F. Feng). 0167-2789/02/$ – see front matter © 2002 Elsevier Science B.V. All rights reserved. PII: S 0 1 6 7 - 2 7 8 9 ( 0 2 ) 0 0 7 2 1 - 2
128
B.-F. Feng et al. / Physica D 175 (2003) 127–138
is usually provided by intrinsic instability of the system. Search for SP solutions and analysis of their stability in media where the dispersion and nonlinearity and, separately, loss and gain may be in balance are subjects of great interest. In particular, it has been shown that stable SPs exist and are stable in a number of one-dimensional (1D) dissipative models built as perturbed versions on the nonlinear Schrödinger equations, such as the cubic–quintic Ginzburg–Landau (CQGL) equation [4–6], and two linearly coupled GL equations [7]. Stable 2D pulses of various types have been found in isotropic [8] and strongly anisotropic [9] 2D versions of the CQGL model. Recently, a 1D model consisting of a mixed Kuramoto–Sivashinsky–KdV equation, linearly coupled to an extra linear dissipative equation, was proposed. The model had the form ut + uux + uxxx + αuxx + γ uxxxx = vx ,
(1)
vt + cvx − Γ vxx = ux ,
(2)
where u and v are the two real wave fields, α the gain, γ and Γ the loss constants, and c the group-velocity mismatch between the two components. The parameters α, γ and Γ are positive, while c may have any sign or be zero. The extra equation (2) makes it possible to stabilize the zero solution u = v = 0, which opens a way to the existence of stable SPs. The existence of stable SPs was predicted by treating the dissipation and gain as small perturbations, and was confirmed by direct simulations [10], thus providing for the first example of stable SP within the framework of a system of the Kuramoto–Sivashinsky (KS) type. These results suggest a natural question, whether stable 2D (completely localized) pulses may be possible in a 2D version of a model of this type. The aim of the present work is to introduce such a physically meaningful model and find stable pulses in it. To develop a 2D generalization of Eqs. (1) and (2), it is natural to begin with a 2D generalization of the KdV equation proper. As is well known, there are two different natural 2D extensions of the KdV equation, the KP equation: (ut + uux + uxxx )x = σ 2 uyy
(3)
and the Zakharov–Kuznetsov (ZK) equation: ut + uux + ux = 0,
(4)
where ≡ ∂x2 + ∂y2 is the isotropic 2D Laplacian, y being the transverse coordinate, and σ 2 a real constant (Eq. (3) with σ 2 = +1 and −1 are called, respectively, the KP-I and KP-II equations). The KP equation is a weakly 2D generalization in the sense that it accounts for slowly varying transverse perturbations of unidirectional KdV solitons moving along the x-direction [11]. The ZK equation, which is a “more isotropic” 2D model, was originally derived to describe weakly nonlinear ion-acoustic waves in a magnetized lossless plasma in two dimensions [12]. It is quite interesting to study 2D modifications of the stabilized KS systems (1) and (2) of both KP and ZK types. The corresponding system of the KP type is considered elsewhere [13], while here we aim to focus on the 2D extension of the systems (1) and (2) based on the ZK equation. Before coupling the ZK equation to the extra one, it must be supplemented by dissipation and gain terms. In fact, such a generalization of the ZK equation was obtained earlier and is known as a 2D version of the Benney equation [14]: ut + uux + ux + αuxx + γ 2 u = 0,
(5)
where α > 0 and γ > 0 have the same meaning as in Eq. (1), and 2 ≡ (∂x2 + ∂y2 )2 is the fourth-order linear isotropic operator. Eq. (5) describes a variety of physical phenomena in two dimensions (chiefly, of hydrodynamic origin), for example, long waves on a thin liquid film [15–17] (for the 1D case, see a review in Ref. [18]), the Rossby waves in rotating atmosphere [19] and the drift waves in plasma [20].
B.-F. Feng et al. / Physica D 175 (2003) 127–138
129
Much attention has been paid to the investigation of SP solutions to the ZK equation (4) and the corresponding Benney equation (5) [21–23]. It was found that, although it is not integrable by means of the inverse scattering transform, the ZK equation has a family of steady-shape 2D SP solutions, moving at an arbitrary velocity in the positive x-direction [22,23]. 2D pulses in the Benney equation were numerically identified in the limiting case of zero dispersion [21]. However, it is obvious that SPs in Eq. (5) with α = 0 and γ = 0 can never be stable in an infinitely long system, since the zero solution, which is the background of the SP, is always unstable against long-wave perturbations in this equation. Indeed, numerical simulations have shown an instability of such SPs even when the dissipation and gain are very weak (α, γ ∼ 10−2 ) [24]. Similar to the 1D system of Eqs. (1) and (2), we develop our 2D model by linearly coupling the generalized ZK (Benney) equation (5) to an extra linear dissipative equation, arriving at a system ut + uux + ux + αuxx + γ 2 u = vx ,
(6)
vt + cvx − Γ vxx = ux ,
(7)
where the parameters have the same meaning as in Eqs. (1) and (2), and no new independent parameter is necessary. It is conjectured that, in the additional equation, the dissipation in the x-direction dominates over the dispersion and nonlinearity, therefore terms of these types are not included in Eq. (7), but they can be added, if necessary, without affecting the results in an essential way. Note here that the system conserves two “masses”: +∞ +∞ +∞ +∞ u(x, y) dx dy, N= v(x, y) dx dy. (8) M= −∞
−∞
−∞
−∞
As concerns physical applications of the system (6) and (7), they constitute a natural model to describe coupled surface and interface waves in two-layered 2D fluid layers, in the case when the surface waves on the top layer are, in isolation, governed by the 2D Benney equation (5) (recall these may be long waves on the surface of a weakly viscous liquid film [15–17], or Rossby waves [19]), while the bottom layer is strongly viscous, cf. a similar interpretation of the 1D system (1) and (2) presented in Ref. [10]. In particular, the linear coupling via the first derivatives between Eqs. (6) and (7) is a known feature of coupled internal waves propagating in multi-layered fluids [25]. Another physical realization of Eqs. (6) and (7), which can be presented here in full detail, assumes that the same fluid layer which is described by the known 2D equation (5), see above, traps a surfactant on its surface (provided that the lateral size of the layer is such that the capillary-wave scale is not negligible). We notice that Eq. (5) was originally derived from, essentially, the x-component of the Navier–Stokes equation for the flow on the surface of the corresponding liquid layer [15–17,19]. If the surfactant is distributed on the surface with a density c + v(x, t), where c and v are, respectively, its constant and small variable parts (|v| c), the gradient of v creates, through the variation of the surface tension, an additional force with the x-component ∼ vx which drives the flow, hence the Navier–Stokes equation adds the term vx to the right-hand side of Eq. (5) (the coefficient in front of this term may be scaled to be 1), so that the equation takes the form of Eq. (6). Further, the evolution of the surfactant density is governed, in a quasi-1D situation, by an obvious advection–diffusion equation: vt + [u(c + v)]x = Γ vxx , where Γ is the surface diffusion constant. With regard to the condition |v| c, the latter equation takes the form (2). We stress that, on top of the above-mentioned particular physical realizations, Eqs. (6) and (7) furnish one of the generic systems that may give rise to stable 2D localized pulses. Stable solitary pulses (SPs) in two dimensions are objects of considerable interest in their own right, and the aim of the present work is to find them in the new system. The paper is organized as follows. In Section 2, we analyze the stability of the zero solution, which is a necessary condition for the stability of SPs in the infinitely long system. In Section 3, an analytical perturbation theory, based on a family of SP solutions of the unperturbed (conservative) counterpart of the system, is developed by treating the gain and dissipation constants α, γ , and Γ in Eqs. (6) and (7) as small parameters. In the conservative system with
130
B.-F. Feng et al. / Physica D 175 (2003) 127–138
α = γ = Γ = 0, a family of 2D SP solutions can be derived from the known cylindrical soliton solutions of the ZK equation proper. Making use of the known approach based on the balance equation for the wave momentum [26], we demonstrate that the combination of the perturbation terms in Eqs. (6) and (7) may select two stationary solitons out of the continuous family existing in the conservative system. Following general principles of the bifurcation theory (see a detailed description and references in Ref. [10]), one may then expect that the soliton with a larger amplitude (and smaller width) is stable, while its broader counterpart with a smaller amplitude is unstable, playing the role of a separatrix between attraction domains of the stable soliton and stable zero state. In Section 4, we present results of direct numerical simulations of the full systems (6) and (7), which confirm the analytical predictions, i.e., the existence of stable 2D SPs. In fact, simulations in a finite domain with periodic boundary conditions produce stable 2D SPs, to which shallow quasi-1D troughs are attached along the x-direction. The existence of this shallow trough is explained by the finite size of the simulation domain. Stable bound states (BSs) of two identical pulses are also found numerically in Section 4.
2. Stability of the zero solution As it was mentioned above, to guarantee the existence of stable SPs, it is first of all necessary to investigate the stability of the trivial solution of Eqs. (6) and (7), u = v = 0. To this end, inserting a 2D perturbation in the form u ∼ exp(ikx + iqy + λt), v ∼ exp(ikx + iqy + λt), into the linearized equations (6) and (7), where k and q are arbitrary real wave numbers of the perturbation, and λ is the corresponding instability growth rate, we arrive at a dispersion equation: [λ − ik(k 2 + q 2 ) − αk 2 + γ (k 2 + q 2 )2 ](λ + ick + Γ k 2 ) + k 2 = 0.
(9)
The stability condition states that both solutions of the quadratic equation (9) for λ must satisfy the inequality Re[λ(k, q)] ≤ 0
(10)
at all real values of k and q. Direct verification of the condition (10) for arbitrary values of the parameters α, γ , Γ , and c is a very cumbersome problem. However, we have found a simple result, linking the stability condition (9) to a similar one studied in the 1D counterpart of the model in Ref. [10]. Namely, if the condition Re[λ(k)] ≤ 0 holds at all real values of k in the case considered in Ref. [10], then the condition (10) holds at all real values of k and q. In other words, the stability of the zero solution in the 1D case guarantees the stability in 2D case. This property for linear stability can be demonstrated by means of a symbolic mathematical software, e.g., Maple. Because the proofs is straightforward but fairly tedious, we only outline it without going into details. The quadratic equation (9) has two roots: λ1,2 = 21 ((Γ − α)k 2 − γ (k 2 + q 2 )2 ) + 21 ik(k 2 + q 2 − c) ± 21 {[γ (k 2 + q 2 )2 − (α + Γ )k 2 ]2 − k 2 [4 + (k 2 + q 2 + c)2 ] − 2ik(k 2 + q 2 + c)[γ (k 2 + q 2 )2 − (α + Γ )k 2 ]}1/2 . By splitting λ1,2 into real parts and imaginary parts, and differentiating the real parts with respect to q, it is found that the derivative of the real part vanishes only at q = 0, hence Re[λ(k, q)] is a monotonically varying function of q in the region 0 < q < ∞. Further, it is easy to see that Re[λ(k, q)] → −∞ as q → ∞, therefore Re[λ(k, q)] must attain its maximum at q = 0, which exactly corresponds to the above-mentioned 1D problem. The proof is complete. This property can also be illustrated by Fig. 1, which is a plot of the instability growth rate Re λ vs. k and q for a set of typical values of parameters. Obviously, when k or q is large, Eq. (9) is satisfied due to the fourth-order
B.-F. Feng et al. / Physica D 175 (2003) 127–138
131
Fig. 1. A three-dimensional (3D) plot of the growth rate Re λ vs. the longitudinal and transverse wavenumbers k and q for α = 0.2, γ = 0.05, c = −1.0 and Γ = 0.55.
dissipative term in Eq. (6). The most dangerous case is that when k and q are relatively small. For any cross-section k = const. in Fig. 1, it can be seen that the growth rate attains its maximum value at q = 0, in compliance with the above proof. Thus, in the present 2D system, the zero solution is stable in the same parameter region which had been identified in Ref. [10] as the stability region of the zero solution in the 1D counterpart of the model. 3. The perturbation theory for SPs By setting α = γ = Γ = 0 but keeping an arbitrary value of c, we arrive at a conservative (unperturbed) version of the systems (6) and (7): ut + uux + ux = vx ,
vt + cvx = ux .
(11)
Our first aim here is to find a soliton solution, moving at a constant velocity s in the x-direction, within the framework of Eq. (11). To this end, we set u(x, y, t) = u(ξ, y),
v(x, y, t) = v(ξ, y)
with ξ ≡ x − st.
(12)
It is easy to show that, as well as in the 1D case [10], the v-component of soliton can be eliminated: v(ξ, y) = (c − s)−1 u(ξ, y).
(13)
132
B.-F. Feng et al. / Physica D 175 (2003) 127–138
Substitution of Eqs. (12) and (13) into the first equation (11) leads to −suξ + uuξ + uξ = (c − s)−1 uξ .
(14)
Integrating (14) with respect to ξ and using the fact that the soliton wave fields vanish at infinity, we arrive at the equation
u + 21 u2 − [s + (c − s)−1 ]u = 0.
(15)
A simple scaling analysis shows that Eq. (15) admits a family of 2D axially symmetric SP solutions of the form A (16) , r 2 ≡ ξ 2 + y 2 , A ≡ 2[s + (c − s)−1 ], u = Aw 2r where w(r) is exactly the known cylindrical soliton solution of the ZK equation, moving at the velocity +1. A numerical expression for it in terms of mapped-Chebyshev polynomials can be found in Ref. [27]. It is necessary to stress that, in contrast with the 2D soliton solutions of the ZK equation proper, which may have only positive velocities, the velocity of the present solution may be negative, provided that the coefficient A defined in Eq. (16) remains positive. For example, in the case c = 0, the condition A > 0 implies that either s > 1, or −1 < s < 0. Note that, if the velocity is negative, belonging to the latter interval, both fields v and u are positive, and, furthermore, v(x, y) > u(x, y), as it follows from Eq. (13). We have checked by direct simulations of Eq. (11) that all the cylindrical soliton solutions given by Eqs. (13) and (16) are stable within the framework of the unperturbed equation (11). We have also performed simulations of both catch-up and head-on collisions between the 2D solitons moving at different velocities. The results clearly show that the collisions are strongly inelastic, hence the conservative system (11) is not an integrable one, as well as the ZK equation proper is not integrable. A typical example of the head-on collision between two solitons with s1 = 2.0 and s2 = −0.8 is shown in Fig. 2. It should be noted here that the panels (c) and (d) in Fig. 2 indicate the possibility of a resonance, i.e., generation of additional pairs of pulses after the head-on collision of positive and negative pulses. Such phenomena were observed in the 1D regularized-long-wave equation [28], as well as in a 2D version of it [29]. In the case of the regularized-long-wave equation, the resonance took place under the condition of vanishing of the lowest-order conserved quantity. In the present case, M and N given by Eq. (8) are corresponding quantities. M is always positive since u is positive for either s > 1 or −1 < s < 0, but N may vanish since v is negative for s > 1 and positive for −1 < s < 0. It is expected that the resonance occurs for vanishing N , but detailed investigation of this problem is left for future study. Next, we go back to the full system of Eqs. (6) and (7), treating the gain and loss parameters as small perturbations. As well as in the 1D model, the relative velocity, δ ≡ c − s, will be regarded as an intrinsic parameter of the unperturbed soliton family. Because the amplitude A in the solution (16) must be positive, the range of meaningful values of δ is restricted by the condition c − δ + δ −1 > 0.
(17)
To select solitons which remain stationary ones in the presence of the gain and losses, we demand that these solutions provide for an equilibrium between the gain-induced pump of the wave momentum into the soliton and the dissipation. The momentum, which is a conserved quantity in the absence of the gain and losses, is 1 +∞ +∞ 2 P = [u (x, y) + v 2 (x, y)] dx dy. (18) 2 −∞ −∞
B.-F. Feng et al. / Physica D 175 (2003) 127–138
133
Fig. 2. A typical inelastic collision between two stable cylindrical solitons with different velocities, simulated within the framework of the conservative system (11) with c = 0. The panels (a) and (b) show the u- and v-fields before the collision, at t = 0, and (c) and (d) display the fields after the collision, at t = 18. The initial velocities of the two solitons are s1 = 2.0 and s2 = −0.8.
An expression for the gain- and dissipation-induced rate of the change of the momentum, i.e., the balance equation for it, can be easily derived from Eqs. (6) and (7): +∞ +∞ dP [αu2x − γ (uxx + uyy )2 − Γ vx2 ] dx dy. (19) = dt −∞ −∞ In the first approximation, we assume that the SP solution still has the form given by Eqs. (13) and (16). Then, substituting the unperturbed solution into Eq. (19) and imposing the balance condition, dP /dt = 0, we obtain an equation to select the value of the parameter δ corresponding to stationary solitons: Γ 1 1 2 C1 α − 2 c−δ+ = 0, (20) − C2 c − δ + δ δ δ where C1 and C2 are two integrals: +∞ dw 2 r dr ≈ 2.47π, C1 = π dr 0
C2 = 2π 0
+∞ d2 w
dr 2
1 dw + r dr
2 r dr ≈ 6.01π,
(21)
which were calculated numerically, using the above-mentioned unperturbed cylindrical soliton solutions w(r) for the ZK equation. Similar to what has been done in the 1D case [10], the balance condition (20) can be eventually cast in the form of a cubic equation for δ: δ 3 + (0.41α˜ − c)δ 2 − δ − 0.41Γ˜ = 0, where α˜ ≡ α/γ , Γ˜ ≡ Γ /γ .
(22)
134
B.-F. Feng et al. / Physica D 175 (2003) 127–138
The general balance equation (19) can be used not only to predict what steady-state pulses are selected by the competition of the losses and gain, but also to describe evolution of nonstationary pulses. However, we do not consider the latter issue in detail, as it is quite difficult to compare the predictions of this sort for the nonstationary states with results of direct numerical simulations. The physical roots of Eq. (22), which select particular cylindrical solitons that remain steady-shape SP solutions to the full system, are those which not only are real, but also satisfy the condition (17). Generally speaking, there may exist up to three physical roots of Eq. (22). As well as in the 1D case, one may argue that the existence of at least two physical roots is a necessary condition for the stability of an SP, as stable and unstable ones may only appear in pairs. Then, as argued in Ref. [10] for the 1D case, the stationary soliton with a larger value of the amplitude is stable, while the other soliton, corresponding to the smaller amplitude, is unstable and plays the role of a separatrix in the present bistable system, the second stable solution being the trivial one, u = v = 0. Moreover, it can be shown that the soliton’s velocities s corresponding to the physical roots of Eq. (22) are always positive, i.e., the existence of a stable soliton moving in the negative x-direction seems impossible. Thus, we can numerically determine a parametric region in which there are exactly two physical roots to Eq. (22), and, simultaneously, the zero solution is stable, which is another condition necessary for the stability of the SP (see above). In particular, such a region of the assumed complete stability of the 2D SP in the parameter plane (α, Γ ) does not exist when c and γ take values 0 and 0.05, respectively. There is a very narrow stability region if c = −1.0 and γ = 0.05, as can be seen in Fig. 3. The zero solution is stable above the dashed line, and there are exactly two physical roots of Eq. (22) inside the region bounded by two solid lines. Under the lower solid line, Eq. (22) has three real roots, but, at most, one of them may be physical, and above the upper solid line, two physical roots bifurcate into a pair of complex roots. Not until α becomes larger than a certain value αcr , do two physical roots appear. For example, it is found that αcr = 0.13 for c = 0, γ = 0.05, and αcr = 0.15 for c = −1, γ = 0.05.
Fig. 3. The expected stability region for the SPs in the parametric plane (α, Γ ) for γ = 0.05 and c = −1.0. The zero solution is stable above the dashed line, and inside the region bounded by two solid lines, Eq. (22) produces two physical solutions.
B.-F. Feng et al. / Physica D 175 (2003) 127–138
135
When c = −1, a sharp wedge bounded by two solid curves is found above the dashed line (see Fig. 3). It is expected that, inside this narrow overlapping region, stable localized cylindrical SPs exist, which is confirmed by direct simulations reported in the next section. In the above-mentioned case c = 0 and γ = 0.05, the region bounded by two solid curves is located, as a whole, beneath the dashed line (not shown here), hence there is no overlapping between the existence of two physical roots of Eq. (22) and the stability of the zero background. The same analysis was performed for other values of the dissipation parameter γ . The results show that variation of γ produces little change in terms of the expected SP stability region. However, the group-velocity mismatch c affects the stability region significantly. We have found that there is no stability region for c ≥ ccr ≈ −0.85, and it appears below this critical value of c.
4. Numerical simulations of the SPs As stated above, it is necessary to directly check whether stable 2D SPs indeed exist in the stability region predicted by the analytical consideration. To this end, we extend the implicit pseudo-spectral scheme applied to the 1D systems (1) and (2) (see the appendix in Ref. [10]). It was found that the numerical integration, even with a coarse grid having x = y = 0.5 and relatively large time steps, t = 0.1 or 0.2, still leads to quite accurate results in most cases. These grid and time steps were used in a majority of simulations; in selected case, the numerical scheme was refined, using essentially smaller x, y and t, in order to check that improvement of the numerical accuracy does not produce any appreciable difference in the results. The initial conditions were taken as the cylindrical solitons of the conservative version of the model, with various values of their amplitude, in order to check whether strongly perturbed pulses relax to stable ones, i.e., whether the stable SPs are attractors in the full system. Results have shown that the cylindrical SPs are indeed stable everywhere inside the expected stability region. Moreover, all the stable pulses were found to be strong attractors. A typical 2D localized pulse with α = 0.2, γ = 0.05, c = −1.0, Γ = 0.53 is shown in Fig. 4. For these values of the parameters, it was found that the initial pulses definitely relax to a single stationary cylindrical SP if their initial amplitude A0 exceeds 0.53, or, in other words, the initial velocity s0 exceeds 0.7. For instance, starting with the initial amplitudes A0 = 2.36 and A0 = 7.87 (to which there correspond s0 = 1.0 and s0 = 2.0) at t = 0, a cylindrical SP develops with A0 = 4.75 and A0 = 5.25 (s = 1.42 and s = 1.51), respectively, at t = 400. For the same values of the parameters, the analytical prediction for the amplitude and velocity of the steady-state pulse (16) and (13) (the one with the larger value of the amplitude) yields Aanal ≈ 3.5 and sanal ≈ 1.2, which is in some disagreement with the numerical results. Actually, this inconsistency comes from Eq. (22). Indeed, as the cylindrical soliton solutions of the unperturbed system were obtained numerically, errors are unavoidable in the computation of the integrals C1 and C2 in Eq. (21). Therefore, the coefficient 0.41 in Eq. (22), which is the ratio of C1 /C2 , may be not very close to an exact value. It is this value which affects the magnitude of the physical roots significantly. For example, a change of the coefficient from 0.41 to 0.45 would correct the analytical prediction for the soliton’s amplitude from 3.5 to 4.54, bringing it much closer to the numerically found value. On the other hand, if the initial amplitude is too small, e.g., A0 = 0.14 (s0 = 0.64), the pulse decays to zero, which is expected too, as the stable zero solution has its own attraction basin. Note that for the second (smaller) steady-state pulse, which is expected to play the role of the separatrix between the attraction basins of the stable pulse and zero solution, the perturbation theory predicts, in the same case, the amplitude A˜ anal ≈ 0.21, so it seems quite natural that the initial pulses with A0 = 0.53 and A0 = 0.14 relax, respectively, to the stable pulse and to zero. It was also found that the relaxation to the stable cylindrical SP might sometimes be very slow. For instance, starting with the cylindrical soliton of the unperturbed system with the amplitude A0 = 0.53, we observed that it
136
B.-F. Feng et al. / Physica D 175 (2003) 127–138
Fig. 4. A stable cylindrical SP produced by direct simulations in the case α = 0.2, γ = 0.05, c = −1.0, Γ = 0.53. The panels (a) and (b) show established shapes of the u- and v-fields, respectively.
evolves into a cylindrical SP with amplitude 1.84 at t = 2000, and to the one with amplitude 4.49 at t = 3000 in the perturbed system with the same parameters as above. The simulations also show that, as a matter of fact, the cylindrical SPs may be stable even when the zero background is unstable, provided that the integration domain is not very large, and periodic boundary conditions are imposed. The explanation for a similar finding in the 1D model, given in Ref. [10], applies to the 2D case as well: the SP, periodically passing across the finite domain, suppresses perturbations faster than they grow. In Fig. 4, a shallow quasi-1D trough, attached to the soliton, is observed. This feature is not generated by numerical errors, and it occurs whenever a stable cylindrical SP develops. A cause for the appearance of the trough can be understood. As mentioned previously, an implicit pseudo-spectral method is used to solve the 2D model, so the Fourier transform converts Eqs. (6) and (7) into the following equations for the Fourier transforms of the fields u and v: duˆ k,q 1 + ikF (uk,q )2 − ik(k 2 + q 2 )uˆ k,q − ik vˆk,q = αk 2 uˆ k,q − γ (k 2 + q 2 )2 uˆ k,q , (23) dt 2 dvˆk,q + ickvˆk,q − ik uˆ k,q = −Γ k 2 vˆk,q , dt
(24)
F{·} standing for the Fourier transform. For the Fourier component u0,q with k = 0, it follows from Eq. (23) that duˆ 0,q = −γ q 4 uˆ 0,q . dt
(25)
Therefore, all the Fourier components uˆ 0,q with nonzero q eventually decay to zero, but it is obvious that these components in a finite domain are not zero. If they are deliberately eliminated, they reappear as a result of the
B.-F. Feng et al. / Physica D 175 (2003) 127–138
137
Fig. 5. A stable BS of two pulses found in the case α = 0.1, γ = 0.05, c = −1.0, Γ = 0.55. The panels (a) and (b) have the same meaning as in Fig. 4.
simulations. As the shallow quasi-1D trough occurs due to the finite size of the integration domain, it is expected to become more shallow when the integration domain is enlarged, and, effectively, to disappear if the integration domain is sufficiently large. This expectation is perfectly confirmed by the simulations. For the same values of the parameters as above, the depth of the quasi-1D trough in the u-component along the central line is about −0.13 when the domain of the size [64, 64] is used, and it drops to −0.06 in the domain of the size [128, 128]. Similar to the 1D case, stable BSs of two cylindrical pulses can be found in the present model, in addition to the single SP. To this end, we ran simulations, with two or more solitons borrowed from the unperturbed system (11) that were initially positioned, with some spacing, along the line y = 0. The solitons had different amplitudes and speeds, so they collided. The simulations have shown formation of stable BSs of two pulses after the collision. The BSs feature two peaks of an equal height (see a typical example in Fig. 5). It should be noted that the line connecting the centers of the peaks is not parallel to the direction y = 0, along which the initial solitons were set. Instead, the line is oblique. We have also checked that, if the amplitudes of the initial pulses and separations between them are changed, the same BS with the equal-amplitude peaks, and oriented in the same direction, finally develops, i.e., the BSs are fairly robust dynamical attractors. Unlike the 1D case, we have not been able to find BSs with three identical peaks. It is also worthy to mention that, in another generic 2D model of the present type, which is based on the KP-I equation, rather than the ZK one, in the lowest approximation, stable SPs were also found, but they cannot form any BS [13].
5. Conclusion In this work, we have extended the recently proposed stabilized KS system in 1D case to the 2D case based on a generalized ZK equation (5). The model includes this equation linearly coupled to an extra linear diffusion
138
B.-F. Feng et al. / Physica D 175 (2003) 127–138
equation. The model applies to the description of surface and interface waves on 2D liquid layers in various physical situations. The additional linear equation stabilizes the zero background in the system, paving a way to the existence of stable 2D localized SPs. In the conservative version of the model, a family of cylindrical solitons can be found, which are stable but collide inelastically. Treating the dissipation and gain in the model as small perturbations, and taking advantage of the balance equation for the wave momentum, we have found that the condition of the balance between the gain and dissipation may select two steady-state cylindrical solitons from their continuous family existing in the conservative system. When the zero solution is stable and, simultaneously, two pulses are picked up by the balance equation for the momentum, the pulse with the larger value of the amplitude is expected to be stable in the infinitely long system, while the other pulse must be unstable, playing the role of a separatrix between attraction domains of the stable pulse and zero solution. These predictions have been completely confirmed by direct simulations. If the integration domain is not very large, and periodic boundary conditions are imposed, the cylindrical pulse develops with a shallow quasi-1D trough attached to it along the x-direction. An explanation for the latter feature was proposed. Finally, stable BSs of pulses with two identical peaks were found too. References [1] M.J. Ablowitz, P.A. Clarkson, Solitons, Nonlinear Evolution Equations and Inverse Scattering, Cambridge University Press, Cambridge, 1991. [2] A. Jeffery, T. Kakutani, SIAM Rev. 14 (1972) 582. [3] R.M. Miura, SIAM Rev. 18 (1976) 412. [4] A.M. Sergeev, V.I. Petviashvili, Sov. Phys. Dokl. 29 (1984) 493. [5] W. van Saarloos, P. Hohenberg, Phys. Rev. Lett. 64 (1990) 749. [6] H.R. Brandt, R.J. Deissler, Phys. Rev. Lett. 63 (1989) 2801; R.J. Deissler, H.R. Brandt, Phys. Rev. Lett. 72 (1994) 478; R.J. Deissler, H.R. Brandt, Phys. Rev. Lett. 74 (1995) 4847; R.J. Deissler, H.R. Brandt, Phys. Rev. Lett. 81 (1998) 3856. [7] B.A. Malomed, H.G. Winful, Phys. Rev. E. 53 (1996) 5365; J. Atai, B.A. Malomed, Phys. Rev. E. 54 (1996) 4371; J. Atai, B.A. Malomed, Phys. Lett. A 246 (1998) 412. [8] L.-C. Crasovan, B.A. Malomed, D. Mihalache, Phys. Rev. E 63 (2001) 016605. [9] H. Sakaguchi, B.A. Malomed, Physica D 159 (2001) 91; H. Sakaguchi, B.A. Malomed, Physica D 123 (2002) 123. [10] B.A. Malomed, B.-F. Feng, T. Kawahara, Phys. Rev. E 64 (2001) 6304. [11] B.B. Kadomtsev, V.I. Petviashvili, Sov. Phys. Dokl. 15 (1970) 539. [12] V.E. Zakharov, E.A. Kuznetsov, Sov. Phys. JETP 39 (1974) 285. [13] B.-F. Feng, B.A. Malomed, T. Kawahara, Stable two-dimensional solitary pulses in linearly coupled dissipative Kadomtsev–Petviashvili equations, Phys. Rev. E 66 (2002), to appear. [14] D.J. Benney, J. Math. Phys. 45 (1966) 150. [15] J. Topper, T. Kawahara, J. Phys. Soc. Jpn. 44 (1978) 663. [16] S. Toh, H. Iwasaki, T. Kawahara, Phys. Rev. A 40 (1989) 5472. [17] S. Melkonian, S.A. Maslowe, Physica D 34 (1989) 255. [18] A. Oron, S.H. Davis, S.G. Bankoff, Rev. Mod. Phys. 69 (1997) 931. [19] V.I. Petviashvili, JETP Lett. 32 (1980) 619. [20] K. Nozaki, Phys. Rev. Lett. 46 (1981) 184. [21] V.I. Petviashvili, Physica D 3 (1981) 329. [22] V.I. Petviashvili, V.V. Yan’kov, Dokl. Akad. Nauk SSSR 267 (1982) 825. [23] H. Iwasaki, S. Toh, T. Kawahara, Physica D 43 (1990) 293. [24] K. Indireshkumar, A.L. Frenkel, Phys. Rev. E 55 (1997) 1174. [25] J.A. Gear, R. Grimshaw, Stud. Appl. Math. 70 (1984) 235. [26] T. Kawahara, S. Toh, Phys. Fluids 28 (1985) 1636; T. Kawahara, S. Toh, in: S. Takeno (Ed.), Dynamical Problems in Soliton Systems, Springer, Berlin, 1985, p. 153; T. Kawahara, S. Toh, Pure Appl. Math. 43 (1989) 95. [27] B.-F. Feng, T. Mitsui, T. Kawahara, J. Comput. Phys. 153 (1999) 467. [28] J. Courtenay Lewis, J.A. Tjon, Phys. Lett. A 73 (1979) 275. [29] T. Kawahara, K. Araki, S. Toh, Physica D 59 (1992) 79.
