# Localized Pulsating Solutions of the Generalized Complex Cubic-Quintic Ginzburg-Landau Equation

Authors: Ivan M. Uzunov; Zhivko D. Georgiev
Journal: Journal of Computational Methods in Physics
 (2014)

---

## Abstract

We study the dynamics of the localized pulsating solutions of generalized complex cubic-quintic Ginzburg-Landau equation (CCQGLE) in the presence of intrapulse Raman scattering (IRS). We present an approach for identification of periodic attractors of the generalized CCQGLE. Using ansatz of the travelling wave and fixing some relations between the material parameters, we derive the strongly nonlinear Lienard-Van der Pol equation for the amplitude of the nonlinear wave. Next, we apply the Melnikov method to this equation to analyze the possibility of existence of limit cycles. For a set of fixed parameters we show the existence of limit cycle that arises around a closed phase trajectory of the unperturbed system and prove its stability. We apply the Melnikov method also to the equation of Duffing-Van der Pol oscillator used for the investigation of the influence of the IRS on the bandwidth limited amplification. We prove the existence and stability of a limit cycle that arises in a neighborhood of a homoclinic trajectory of the corresponding unperturbed system. The condition of existence of the limit cycle derived here coincides with the relation between the critical value of velocity and the amplitude of the solitary wave solution (Uzunov, 2011).


---

## Body

## 1. Introduction

As is well known the complex cubic Ginzburg-Landau equation has been used to describe a variety of phenomena including second-order phase transitions, superconductivity, superfluidity, Bose-Einstein condensation, and strings in field theory [1, 2].Perturbative methods are used for the description of spatiotemporal pattern formations in systems driven away from equilibrium near the threshold where the nonlinearities are weak and the spatial and temporal modulations of the unstable modes are slow [1]. One of them is the method of “amplitude (model) equations” for the envelope function of unstable mode [1]. Model equations depend on the type of the linear instability. Universal amplitude equations are the real and complex cubic Ginzburg-Landau equation (CCGLE) as well as their generalizations, like the coupled complex cubic Ginzburg-Landau equation and the complex cubic-quintic Ginzburg-Landau equation (CCQGLE). CCGLE describes the evolution of the envelope function of unstable mode for any process exhibiting a Hopf bifurcation. As a model equation CCGLE is applied to the study of oscillatory uniform instability in lasers, oscillatory periodic instability in Rayleigh- Bénard convection in binary mixtures as well as electrohydrodynamic instabilities in nematic liquid crystals (see [1, 2] and references therein).CGLE also appears as a continuous limit of the amplitude equations that describe the chain of equal weakly nonlinear oscillators with the nearest neighbor interaction. The weakly nonlinear oscillators exhibit a limit cycle.It has been recently established that in optics the one-dimensional complex cubic-quintic Ginzburg-Landau equation (CCQGLE) can model soliton transmission lines [3, 4] as well as passively mode-locked laser systems [5, 6].CCQGLE has exact chirped solitary wave solutions [7–10]. Their mathematical nature and the consistent way of derivation have been deeply discussed in [10]. Numerical solutions of CCQGLE could be divided into two groups: localized fixed-shape solutions and localized pulsating solutions. Novel numerical pulsating solutions of the CCQGLE, namely, pulsating, creeping, snaking, and erupting solutions have been reported in [11]. Chaotic pulsating solutions and period doubling were reported in [12]. Detailed review of the analytical and numerical solutions of CCQGLE can be found in [13–15].A theoretical approach for analysis of the observed solutions of CCQGLE using the variational method has been reported [16, 17]. The resulting Euler-Lagrange equations have been analyzed for the existence of periodic, quasiperiodic, and chaotic attractors [16, 17]. It has been shown that the different numerically observed solutions of CCQGLE (dissipative solitons) may be related to the stable periodic attractors of the Euler-Lagrange equations [16, 17].The influence of the higher order effects, namely, the third-order of dispersion (TOD), intrapulse Raman scattering (IRS), and self-steeping effect (see for detailed description these effects [18–20]) on the description of fiber laser operation has been studied in [21]. In order to perform this analysis we consider a generalized CCQGLE that includes the higher-order effects [21]. The existence of the exact chirped solitary solution of this generalized CCQGLE has been reported for the first time in [22]. Very recently, the influence of these high order effects on the dynamics of pulsating, erupting, and creeping solutions using the generalized CCQGLE has been studied in [23, 24]. Generally, it has been shown that these higher-order effects can have strong impact on these solutions. It was established that under the influence of IRS and TOD, the plain pulsating and the creeping solutions can lose their pulsating behavior [23]. A further observation was made that in the presence of all higher-order effects the explosions of an erupting soliton can be reduced and even eliminated [23, 24]. Numerical findings of [23] suggest that in the presence of IRS the pulsating regime is stabilized (see Figures 2(a) and 2(d) in [16]), or, in other words, a periodic attractor in (1) appears.The main aim of this paper is to examine analytically the influence of the intrapulse Raman scattering on the localized pulsating solutions of the generalized CCQGLE in an attempt to explain the very recent numerical results of [23]. In fact we propose a theoretical approach for identification of periodic attractors of (1). We introduce a dynamical system with finite degrees of freedom or a system of ordinary differential equations (SODE) related to generalized CCQGLE. We next identify the periodic attractors of this SODE. The last step of this approach is to solve numerically the generalized CCQGLE in order to compare periodic attractors of SODE with those of CCQGLE (this step will be the subject of further investigation).In order to introduce SODE, we use ansatz of the travelling wave. Fixing some relations between the material parameters of (1), we have succeeded to derive the equation of strongly nonlinear Lienard-Van der Pol oscillator (see (11)) for the amplitude of the nonlinear wave. After identifying the possible equilibrium points of (11) in the general case, we had to fix some values of his coefficients in order to demonstrate our approach. The simpler version of the equation of strongly nonlinear Lienard-Van der Pol oscillator, namely, the equation of strongly nonlinear Duffing-Van der Pol oscillator, has been the object of intensive study by means of different perturbation methods [25]. Interestingly enough, even for the large values (larger than 1) of the small parameter, appearance of stable periodic attractors-limit cycles was observed [25]. Here we apply the Melnikov method [26–29] to analyze the possible existence of limit cycles in the equation of strongly nonlinear Lienard-Van der Pol oscillator and the equation of strongly nonlinear Duffing-Van der Pol oscillator in two cases, around closed oval trajectories having finite temporal period and around homoclinic trajectory with infinite temporal period. For the fixed values of the coefficients of the equation of strongly nonlinear Lienard-Van der Pol oscillator we prove the existence of a single limit cycle that arises around a closed phase trajectory of the unperturbed system. We could then expect that the observed limit cycle will be related to the corresponding periodical attractor of (1).The second aim of this paper is to explore analytically the influence of the intrapulse Raman scattering on the localized pulsating solutions of the particular case of generalized CCQGLE, namely, (1) with ε
=
μ
=
ν
=
0. Physically this equation describes the bandwidth limited amplification and has recently been investigated in [30, 31]. Equation of the strongly nonlinear Duffing-Van der Pol oscillator has been introduced [30] in order to describe the influence of IRS on the soliton solutions. Here we apply the Melnikov method to the equation of strongly nonlinear Duffing-Van der Pol oscillator and prove the existence of a limit cycle that arises in a neighborhood of a homoclinic trajectory of the unperturbed system. Due to the large period of periodic movement, the limit cycle transforms into perturbed soliton solution. The condition of existence of the limit cycle derived here coincides with the condition of existence of perturbed soliton solution derived by means of the soliton perturbation theory in [30] and the critical value of the homoclinic bifurcation parameter calculated in [30] by means of the hyperbolic perturbation method [32] and hyperbolic Lindstedt-Poincare method [33]. As a result of [30] which has been successfully confirmed by the numerical solution of the corresponding generalized CCQGLE in [31], we can expect that the obtained results here for (1) will also well correlate to the numerical solution of (1).The results presented here have been first announced (but not published) in [34] and in a brief and different form reported in [35].The paper is organized as follows. First, the physical meaning and applications of the generalized CCQGLE are presented in Section2. Second, the derivation of the equation of strongly nonlinear Lienard-Van der Pol oscillator (see (10) below) is given in Section 3. Next, in Section 4 the results of the analysis of (10) by means of the Melnikov method are given. In Section 5 we apply the Melnikov method to the investigation of the influence of the IRS on the bandwidth limited amplification [30]. In Section 6 we discuss the obtained results. Finally, we make our conclusions in Section 7. A classification of the equilibrium points of (10) is given in Appendix A. Appendix B contains a description of the calculation of Abelian integrals necessary for the application of the Melnikov method. Appendix C gives the condition of existence of perturbed solitary wave solution [30] derived by the hyperbolic perturbation method and the hyperbolic Lindstedt-Poincare method as well as the relationship between the equilibrium amplitude and the velocity of the perturbed soliton solution, obtained by means of soliton perturbation theory [30].


## 2. Basic Equation

The propagation of ultrashort pulses in the presence of spectral filtering, linear and nonlinear gain/loss, and intrapulse Raman scattering is described by the following generalized CCQGLE [23, 24]:
(1)
i
∂
U
∂
x
+
1
2
∂
2
U
∂
t
2
+
|
U
|
2
U
=
i
δ
U
+
i
β
∂
2
U
∂
t
2
+
i
ε
|
U
|
2
U
-
ν
|
U
|
4
U
+
i
μ
|
U
|
4
U
+
γ
U
∂
∂
t
(
|
U
|
2
)
,
where U is the normalized envelope of the electric field, x is the normalized propagation distance, t is the retarded time, δ is the linear gain or loss coefficient, β describes spectral filtering (gain dispersion), and ε is related to nonlinear gain-absorption process, μ, if negative accounts for the saturation of the nonlinear gain, ν, if negative corresponds to the saturation of the nonlinear refractive index. The last term in (1) describes the intrapulse Raman scattering (IRS)—the important nonlinear physical effect which has to be taken into account when femtosecond optical pulses propagate in optical fibers. It is related to the first moment of the nonlinear response function (the slope of the Raman gain spectrum) and leads to new physical phenomena as soliton self-frequency shift and breakup of the N-soliton bound states (for the review see [18–20]).Equation (1) has been used to model solving of the problem of linear-wave growth in bandwidth-limited-amplified soliton transmission systems [3, 4]. In the context of solid state lasers, (1) (ν
=
μ
=
γ
=
0) has been proposed for description of a fast saturable absorber and additive pulse mode-locking [5] and later has been used (γ
=
0) as a model for self-limited additive pulse mode-locking [36]. Next it turned out that (1) (ν
=
μ
=
γ
=
0) was also applicable as a main equation for the mode-locked fiber lasers [37]. Recently the all-normal-dispersion passively mode-locked fiber lasers have been successfully described by means of (1) (ν
=
γ
=
0) [38, 39]. Generally, the CCQGLE (γ
=
0) has proved to be a good model for the real mode-locked lasers (for the review see [40, 41]).As a result of intensive numerical investigation of CCQGLE, some areas in the space of physical parameters(
δ
,
β
,
ε
,
μ
,
ν
) have been established, in which there exist stable localized solutions of CCQGLE [14, 42].


## 3. Derivation of Lienard-Van der Pol Equation

Now we will look for the stationary pulse solutions of (1) in the form:
(2)
U
(
x
,
t
)
=
u
(
ξ
)
exp
⁡
⁡
(
i
(
f
(
ξ
)
+
K
x
)
)
,
where ξ
=
t
-
M
x and M and K are real numbers. M has a meaning of the unknown inverse equilibrium velocity. Inserting (2) into (1) the following nonlinear system of ordinary differential equations for the amplitude and phase functions u
(
ξ
) and f
(
ξ
) is obtained:(3a)
-
β
u
′′
+
f
′
u
′
-
M
u
′
-
δ
u
+
β
u
(
f
′
)
2
+
(
1
2
)
u
f
′′
-
ε
u
3
-
μ
u
5
=
0
,
(3b)
(
1
2
)
u
′′
+
2
β
u
′
f
′
-
2
γ
u
2
u
′
-
K
u
+
M
u
f
′
-
(
1
2
)
u
(
f
′
)
2
+
β
u
f
′′
+
u
3
+
ν
u
5
=
0
.Summing up (3a) multiplied with 1
/
(
2
β
) and (3b) we obtain
(4)
-
(
2
δ
+
4
β
K
)
u
+
4
β
M
u
f
′
+
(
1
+
4
β
2
)
u
f
′′
+
2
(
1
+
4
β
2
)
u
′
f
′
-
(
2
M
+
8
γ
β
u
2
)
u
′
+
(
4
β
-
2
ε
)
u
3
+
(
4
ν
β
-
2
μ
)
u
5
=
0
.Multiplying (4) with u
(
ξ
) we get(5a)
d
d
ξ
[
(
1
+
4
β
2
)
u
2
f
′
-
M
u
2
-
2
β
γ
u
4
]
+
[
-
(
2
δ
+
4
β
K
)
u
2
+
4
β
M
u
2
f
′
+
(
4
β
-
2
ε
)
u
4
1
+
(
4
ν
β
-
2
μ
)
u
6
]
=
0
.Assume the following conditions are satisfied:(5b)
1
+
4
β
2
=
4
M
β
;
M
=
2
(
δ
+
2
β
K
)
;
ε
=
β
(
γ
+
2
)
;
μ
=
2
β
ν
,and introduce the function(6)
φ
(
ξ
)
=
u
(
ξ
)
2
[
-
M
+
(
1
+
4
β
2
)
f
′
(
ξ
)
-
2
β
γ
u
(
ξ
)
2
]
.Equation (5a) can be written in the form, φ
′
(
ξ
)
+
φ
(
ξ
)
=
0, with corresponding solution: φ
(
ξ
)
=
φ
0
e
-
ξ
,
φ
0
=
const.The phase functionf
(
ξ
) then is completely determined by the amplitude function u
(
ξ
):
(7)
f
(
ξ
)
=
M
(
1
+
4
β
2
)
ξ
+
φ
0
(
1
+
4
β
2
)
∫
e
-
ξ
u
2
d
ξ
+
2
β
γ
(
1
+
4
β
2
)
∫
u
2
d
ξ
.If we replacef
(
ξ
) with (7) in (3a), the following equation for the amplitude function u
(
ξ
) is obtained:
(8)
u
′′
+
(
δ
β
-
M
2
(
1
+
4
β
2
)
2
-
4
φ
0
e
-
ξ
γ
β
(
1
+
4
β
2
)
2
)
u
+
(
ε
β
-
4
M
γ
β
(
1
+
4
β
2
)
2
)
u
3
+
(
μ
β
-
4
γ
2
β
2
(
1
+
4
β
2
)
2
)
u
5
-
φ
0
2
e
-
2
ξ
(
1
+
4
β
2
)
2
1
u
3
=
4
(
1
+
4
β
2
)
(
-
β
M
+
γ
u
2
)
u
′
.Let us introduce the following quantities:(9)
c
1
=
-
1
+
16
δ
β
16
β
2
,
c
2
=
4
β
γ
(
1
+
4
β
2
)
,
c
3
=
2
+
4
γ
β
2
(
1
+
4
β
2
)
,
c
4
=
1
(
1
+
4
β
2
)
2
,
c
5
=
μ
β
-
4
γ
2
β
2
(
1
+
4
β
2
)
2
,
σ
=
4
γ
(
1
+
4
β
2
)
.The coefficientsc
1, c
2, c
3, c
4, c
5, σ depend on δ, β, μ, γ. Using the quantities (9) and replacing K and M as well as ε by (5b), (8) can be written in the form:
(10)
u
′′
+
(
c
1
-
c
2
φ
0
e
-
ξ
)
u
+
c
3
u
3
+
c
5
u
5
-
c
4
φ
0
2
e
-
2
ξ
u
3
=
(
-
1
+
σ
u
2
)
u
′
.As we can see from (10), by increasing ξ, the role of the terms proportional to φ
0
e
-
ξ and φ
0
2
e
-
2
ξ will decrease. The constant φ
0 can be freely chosen to be smaller than unity. Then we could expect that the terms proportional to c
2 and c
4 can be safely neglected, in case we consider the dynamics of (10) at sufficiently long distances ξ. Neglecting the terms proportional to c
2 and c
4 we obtain
(11)
u
′′
+
c
1
u
+
c
3
u
3
+
c
5
u
5
=
(
-
1
+
σ
u
2
)
u
′
.Equation (11) is the equation of a strongly nonlinear Lienard-Van der Pol oscillator.In order to check the correctness of the approximation of (10) with (11), we have numerically solved the systems corresponding to these equations for some fixed values of parameters and compared the obtained results. As the aim of this work is to find the limit cycles, we should compare the abilities of (11) to describe the limit cycles predicted by (10).We have found that for the values of quantitiesc
1, c
2, c
3, c
4, c
5, σ of order of unity, and φ
0
≤
0.1, for properly long times both systems give comparable results. The difference appears at the initial stage and it is proportional to the value of φ
0. For very small values of φ
0 (φ
0
≪
1), the difference is so small that it can be safely neglected. For values of c
2 and c
4 smaller than unity, the acceptable value of φ
0 increases.In order to demonstrate these results we fix the following values of the coefficientsc
1
=
-
1, c
2
=
1.005, c
3
=
5
/
2, c
4
=
0.5024, c
5
=
-
1, and φ
0
=
0.1. The values of c
1, c
2, c
3, c
4, c
5 are determined by the values of δ
=
-
0.371895, β
=
0.497519, ε
=
1.49504, μ
=
-
0.371895, ν
=
-
0.37375, γ
=
1.005 given in Section 6. The value of σ
=
2.15 is chosen to satisfy the inequality given by (26). In this case according to our analysis performed in Section 4, (11) possesses a stable limit cycle. In order to demonstrate the existence and stability predicted by (11) limit cycle, we compare the results of the numerical integration of systems that correspond to (10) and to (11) for two different initial conditions: (a) u
(
0
)
=
0.2, d
u
/
d
ξ
(
0
)
=
0, and (b) u
(
0
)
=
1.0, d
u
/
d
ξ
(
0
)
=
0. The obtained results for functions u
(
ξ
) and d
u
/
d
ξ
(
ξ
) are shown in Figures 1(a) and 1(b) and Figures 2(a) and 2(b), respectively.Evolution ofu
(
ξ
) (a) and d
u
/
d
ξ (b), according to (9) (blue) and (10) (red), for c
1
=
-
1, c
2
=
1.005, c
3
=
5
/
2, c
4
=
0.5024, c
5
=
-
1, σ
=
2.15, φ
0
=
0.1. Initial conditions: u
(
0
)
=
0.2, d
u
/
d
ξ
(
0
)
=
0.
(a)
(b)Evolution ofu
(
ξ
) (a) and d
u
/
d
ξ (b), according to (9) (blue) and (10) (red), for c
1
=
-
1, c
2
=
1.005, c
3
=
5
/
2, c
4
=
0.5024, c
5
=
-
1, σ
=
2.15, φ
0
=
0.1. Initial conditions: u
(
0
)
=
1.0, d
u
/
d
ξ
(
0
)
=
0.
(a)
(b)As can be seen from Figures1(a) and 1(b), the discrepancy between the two equations appears at the initial stage (till ξ
≈
15). For larger ξ the behavior of both functions is very similar. The periodic behavior of the functions u
(
ξ
) and d
u
/
d
ξ
(
ξ
) shows the existence of the stable limit cycle.As can be seen from Figures2(a) and 2(b), the behavior of both functions with ξ is quite similar. Again, the periodic behavior of the functions u
(
ξ
) and d
u
/
d
ξ shows the existence of the stable limit cycle. Comparing the results shown in Figures 1(a) and 1(b) and Figures 2(a) and 2(b), we see that the stable limit cycle predicted by (11) also exists for (10). So, we believe that the approximation of (10) with (11) is acceptable in a proper region of parameters.From the definition ofγ
>
0 it follows that c
3
>
0 and σ
>
0. In our further investigation we will assume that these inequalities are satisfied. Each of the coefficients c
1 and c
5 can be either positive or negative. Moreover, we will assume that the right hand side of (11), which can be regarded as a dissipative term, is a small quantity, or small perturbation. Then the unperturbed equation corresponding to (11) is exactly the so-called cubic-quintic Duffing equation. Many methods for the analysis of (11) require in advance that we obtain the phase portrait of the cubic-quintic Duffing equation. For this reason a classification of the equilibrium points of this equation is given in Appendix A.


## 4. The Melnikov Method for Analysis of Lienard-Van der Pol Equation

Consider (11) under the following conditions: c
1
<
0, c
3
>
0 and c
5
<
0. After applying the scaling
(12)
u
=
c
1
c
5
4
X
,
ξ
=
-
1
c
1
τ
,(11) takes the form
(13)
d
2
X
d
τ
2
-
X
+
c
~
3
X
3
-
X
5
=
ε
~
(
-
1
+
b
~
X
2
)
d
X
d
τ
,
where
(14)
c
~
3
=
-
(
c
3
c
1
)
c
1
c
5
>
0
,
ε
~
=
(
1
-
c
1
)
>
0
,
b
~
=
σ
c
1
c
5
>
0
.Further we assume thatε
~ is a small parameter. Then the right hand side of (13) is a small perturbation. Taking into account the expression for c
~
1 given in (9), we can conclude that c
~
1 is a small parameter on condition that the spectral filtering β has a small value, or the module of linear gain or loss coefficient |
δ
| has a large value, or when both conditions are fulfilled.Our goal is to establish that (13) admits limit cycles and to analyze them. For this purpose we will use two main mechanisms for the arising of limit cycles: (i) arising of limit cycles around a closed trajectory belonging to a family of such closed trajectories of the unperturbed equation; (ii) arising of limit cycles around a homoclinic trajectory of the unperturbed equation. The limit cycles in the first case correspond to a periodic solution in the initial wave equation. The limit cycle around a homoclinic trajectory in the second case (with infinitely long period) corresponds to a soliton solution in the initial wave equation.The analysis of limit cycles will be carried out using the Melnikov theory [27], which also allows finding the conditions for the arising of limit cycles. A classification of the equilibrium points of the unperturbed cubic-quintic Duffing equation depending on the values of its coefficients is given in Appendix A for convenience. In this situation, one can easily obtain the phase portrait of the unperturbed equation and determine whether the preliminary requirements of the Melnikov theory (presence of a family of closed trajectories, or presence of homoclinic trajectories) are fulfilled. In order to be able to perform the computations and to get a better representation of the results we choose c
~
3
=
5
/
2. If we choose another value for c
~
3, the phase portrait of the unperturbed cubic-quintic Dufing equation is changed, but when the preliminary requirements of the Melnikov theory are fulfilled, the method of study of the limit cycles remains the same. Only the calculations leading to other limit cycles with other parameters are changed.Taking into account the value ofc
~
3, (13) can be written as the following system:
(15)
X
˙
=
Y
,
Y
˙
=
X
-
(
5
2
)
X
3
+
X
5
+
ε
~
(
-
1
+
b
~
X
2
)
Y
,
where the dot means differentiation with respect to τ. The perturbation functions in this system are
(16)
p
(
X
,
Y
)
=
0
,
q
(
X
,
Y
)
=
(
-
1
+
b
~
X
2
)
Y
.The unperturbed system (i.e., system (15) with ε
~
=
0
) is
(17)
X
˙
=
Y
,
Y
˙
=
X
-
(
5
2
)
X
3
+
X
5
.The last system is Hamiltonian with Hamiltonian function(18)
H
(
X
,
Y
)
=
(
1
2
)
Y
2
-
(
1
2
)
X
2
+
(
5
8
)
X
4
-
(
1
6
)
X
6
=
h
,
where the parameter h corresponds to a constant Hamiltonian level.System (17) has the following equilibrium points, (
0
,
0
), which is a “saddle point,” (
±
X
2
,
0
)
=
(
±
1
/
2
,
0
), which are “centers,” and (
±
X
4
,
0
)
=
(
±
2
,
0
), which are “saddle points.” The phase portrait of system (17) is symmetrical concerning X- and Y-axis and is shown in Figure 3. Note that this phase portrait corresponds to Subcase 6.1 given in Appendix A.Figure 3
Phase portrait of the unperturbed Hamiltonian system (17). The five equilibrium points, two centers (green points) and three saddle points (red points), can be clearly observed. The homoclinic trajectory forms figure “eight loop.”According to the Melnikov theory, limit cycles in a given perturbed Hamiltonian system can arise around a closed phase trajectory of the unperturbed Hamiltonian system. Because of this, we will further be interested only in the family of closed trajectories localized within the right half of the figure “eight loop,” that is, the closed trajectories localized between the center(
X
2
,
0
)
=
(
1
/
2
,
0
) and the homoclinic loop. Each phase trajectory corresponds to a constant value of the Hamiltonian function; that is, H
(
X
,
Y
)
=
h, h
=
const. The value of the Hamiltonian function on the homoclinic trajectory is h
hom
=
0 and in the “center” the value is h
c
=
-
11
/
96. It can be easily calculated that the homoclinic trajectory intersects the positive part of the X-axis at the point X
hom
=
(
15
-
33
)
/
8
=
1.0759. A given closed trajectory localized within the right half of the figure “eight loop” has the following equation:
(19)
Γ
0
(
h
)
:
(
1
2
)
Y
2
-
(
1
2
)
X
2
+
(
5
8
)
X
4
-
(
1
6
)
X
6
=
h
,
hhhhhh
-
(
11
96
)
<
h
<
0
.Having these notions in mind, we can obtain the Melnikov function for the system (15):
(20)
M
(
h
)
=
∮
Γ
0
(
h
)
(
q
(
X
,
Y
)
d
X
-
p
(
X
,
Y
)
d
Y
)
=
∮
Γ
0
(
h
)
(
-
1
+
b
~
X
2
)
Y
d
X
.Introducing the functions(21)
I
k
(
h
)
=
∮
Γ
0
(
h
)
X
2
k
Y
d
X
,
k
=
0
,
1
,
(22)
P
(
h
)
=
I
1
(
h
)
I
0
(
h
)
,the Melnikov function and its derivative with respect toh can be written as
(23)
M
(
h
)
=
-
I
0
(
h
)
+
b
~
I
1
(
h
)
=
I
0
(
h
)
[
-
1
+
b
~
P
(
h
)
]
,
(24)
M
′
(
h
)
=
-
I
0
′
(
h
)
[
-
1
+
b
~
P
(
h
)
]
+
b
~
I
0
(
h
)
P
′
(
h
)
.Note thatI
0
(
h
) and I
1
(
h
) are Abelian integrals.We will look for the zeros ofM
(
h
) with (
-
11
/
96
)
<
h
<
0. To find these zeros we need additional investigations of the functions I
0
(
h
), I
1
(
h
), and P
(
h
). The facts necessary for this purpose are collected in Lemma B.1 presented in Appendix B. A crucial factor in this consideration is the fact that the function P
(
h
) is positive and strictly monotone decreasing in the interval (
-
11
/
96
)
<
h
<
0. The graph of the function P
(
h
) obtained by numerical integration is shown in Figure 4.Figure 4
The graph of the functionP
(
h
)
=
I
1
(
h
)
/
I
0
(
h
) obtained by numerical integration of the Abelian integrals given by (
B
.
7
) with the help of Mathematica.A given zero of the Melnikov functionh
0 satisfies the equations
(25)
M
(
h
0
)
=
0
,
-
1
+
b
~
P
(
h
0
)
=
0
.Now we can make some conclusions. SinceP
(
h
) is a strongly monotone decreasing function, the second equation in (25) can have only one zero. The function P
(
h
) satisfies the inequality 0.4408
<
P
(
h
)
<
0.5. Therefore the Melnikov function has a single zero in the case when the following inequality is satisfied:
(26)
2
<
b
~
<
2.2686
.In this case the system (15) has a limit cycle Γ
ε
~
(
h
0
), which is localized in O
(
ε
~
)-neighborhood of the closed curve Γ
0
(
h
0
).According to the Melnikov theory, the stability of the limit cycle is determined by the sign of quantityε
~
M
′
(
h
0
): stable for ε
~
M
′
(
h
0
)
<
0 and unstable for ε
~
M
′
(
h
0
)
>
0 [27]. From (24) and Lemma B.1 in Appendix B, there follows that
(27)
M
′
(
h
0
)
=
b
~
I
0
(
h
0
)
P
′
(
h
0
)
<
0
.Taking into account the inequalityε
~
>
0, we conclude that the limit cycle Γ
ε
~
(
h
0
) is stable.Figure5 shows the limit cycle of system (15) obtained for h
0
=
-
0.095 and b
~
=
2.02.Figure 5
The phase portrait of system (15) obtained by numerical solution.Figure5 presents the results from the calculation of three different initial conditions for the time τ
=
1000. One of the initial conditions (
0.865
,
0
) coincides with the limit cycle and the corresponding phase trajectory stays the same during the calculation. Another initial condition (
0.815
,
0
) is inside the limit cycle and with the increase of time the trajectory approaches the limit cycle from the inside. The last initial condition (
0.915
,
0
) is outside the limit cycle and with the increase of time the trajectory approaches the limit cycle from the outside.The presented theory allows us to synthesize a system of the type (15), having in advance assigned limit cycle. Let h
=
h
0 be the Hamiltonian level around which we want to a limit cycle to arise. After computing P
(
h
0
), we get from (25)
(28)
b
~
=
1
P
(
h
0
)
,
(
-
11
96
)
<
h
0
<
0
.The obtained in this way value ofb
~ provides system (15) with the desired limit cycle.


## 5. The Melnikov Method for Analysis of Duffing-Van der Pol Equation

In the previous section it was shown that, under certain conditions, (11) has a simple stable limit cycle. This limit cycle arises around a closed phase trajectory of the unperturbed system, which is an oval curve. In the present section we will show that (11) can have a limit cycle that arises in a neighborhood of a homoclinic trajectory of the unperturbed system. In the first case, the closed oval curve, respectively, the emerging limit cycle, has a finite period. In the second case, as is well known, the homoclinic trajectory is described by the solution of the unperturbed equation for infinitely long time. On account of this, the limit cycles in the two cases have different properties and this leads to different behaviors of the initial wave equation.The process of arising of limit cycles in a perturbed Hamiltonian system around a homoclinic trajectory of the unperturbed Hamiltonian system is considered for the first time by Andronov; see Andronov et al. [26] (therefore this process is called Andronov mechanism). Later this theory was developed by Roussarie and shown in a form convenient for our applications in [43]. Useful and detailed information on these subjects is given in the book of Han and Yu [44].Formally the analysis below will be performed for the casec
1
<
0, c
3
>
0, and c
5
=
0. (Note however, that c
1, c
3 in this paragraph are short notations of the magnitudes c
^
1, c
^
3 defined in Appendix C.) This case has been studied earlier in [30] by means of the hyperbolic perturbation method [32] and the hyperbolic Lindstedt-Poincare perturbation method for homoclinic motion [33].After applying the scaling(29)
u
=
-
c
1
c
3
X
,
ξ
=
-
1
c
1
τ
,(10) takes the form
(30)
d
2
X
d
τ
2
-
X
+
X
3
=
ε
~
(
-
1
+
b
~
X
2
)
d
X
d
τ
,
where
(31)
ε
~
=
(
1
-
c
1
)
>
0
,
b
~
=
-
σ
(
c
1
c
3
)
>
0
.Further, we will assume thatε
~ is a small parameter. (Having in mind the definition of c
~
1
=
-
2
K in Appendix C this requires that the parameter K be sufficiently large.)Equation (30) can be written as a perturbed Hamiltonian system in the following way:
(32)
X
˙
=
Y
+
ε
~
[
-
X
+
(
1
3
)
b
~
X
3
]
,
Y
˙
=
X
-
X
3
,
where the dot means differentiation with respect to τ. The perturbation functions in this system are
(33)
p
(
X
,
Y
)
≡
-
X
+
(
1
3
)
b
~
X
3
,
q
(
X
,
Y
)
≡
0
.The unperturbed system (withε
~
=
0) is Hamiltonian with Hamiltonian function
(34)
H
(
X
,
Y
)
=
(
1
2
)
Y
2
-
(
1
2
)
X
2
+
(
1
4
)
X
4
=
h
.Moreover, the unperturbed system has three equilibrium points—a hyperbolic “saddle point”(
X
,
Y
)
=
(
0
,
0
) and two nondegenerate “centers” (
X
,
Y
)
=
(
±
1
,
0
). The value of the Hamiltonian function at the two “centers” is H
(
X
,
Y
)
=
h
=
-
1
/
4. It is necessary to note that the point (
0
,
0
) is also an equilibrium point “saddle point” for the perturbed system (32).The equationH
(
X
,
Y
)
=
h
=
0 defines a symmetric double homoclinic loop L
0 (figure eight-loop), consisting of two homoclinic orbits L
0
+ and L
0
- connecting the “saddle point” (
0
,
0
) to itself; that is, L
0
≡
L
0
-
∪
L
0
+. The homoclinic orbits L
0
+ and L
0
- are localized, respectively, in the half planes X
>
0 and X
<
0. The solution of the unperturbed system along the homoclinic orbits L
0
+ and L
0
- can be expressed in the time domain in the following way:(35a)
X
=
φ
0
(
τ
)
=
±
2
sech
τ
,
(35b)
Y
=
ψ
0
(
τ
)
=
∓
2
sech
τ
th
τ
.On the other hand, the equationH
(
x
,
y
)
=
h defines three continuous families of periodic orbits—{
Γ
-
(
h
)
:
(
-
1
/
4
)
<
h
<
0
,
X
<
0
}, {
Γ
+
(
h
)
:
(
-
1
/
4
)
<
h
<
0
,
X
>
0
}, and {
Γ
(
h
)
:
h
>
0
}. The periodic orbits Γ
-
(
h
) and Γ
+
(
h
) surround the centers (
±
1
,
0
) and tend, respectively, to L
0
- and L
0
+ as h
→
(
0
-
). The limit of Γ
(
h
) as h
→
(
0
+
) is the double homoclinic loop L
0.In general, the Melnikov function can be defined on each of the mentioned three families of periodic orbits and we have three Melnikov functions [44]. Let us introduce the Melnikov function M
(
h
) for the perturbed Hamiltonian system (32) defined along the periodic orbit Γ
+
(
h
). It is well known [43, 45] that the Melnikov function has the following asymptotic expansion near h
=
0 (0
<
|
h
|
≪
1, or h
→
(
0
-
))(36a)
M
(
h
)
=
s
0
+
s
1
h
ln
⁡
|
h
|
+
s
2
h
+
s
3
h
2
ln
⁡
⁡
|
h
|
+
s
4
h
2
+
⋯
,
where
(36b)
s
0
=
M
(
0
)
,
s
1
=
(
∂
p
∂
X
+
∂
q
∂
Y
)
|
(
X
,
Y
)
≡
saddle
point
.The quantityM
(
0
) is the Melnikov function computed along the separatrix L
0
+ (the computation along L
0
- is identical and leads to the same results) and for the system (32) it has the form
(37)
M
(
0
)
=
∫
-
∞
+
∞
[
q
(
φ
0
,
ψ
0
)
φ
˙
0
-
p
(
φ
0
,
ψ
0
)
ψ
˙
0
]
d
τ
.The coefficientss
0, s
1, s
2,… play an important role in the study of the limit cycle bifurcation. According to a theorem of Roussarie [43], if the following relations hold
(38)
s
0
=
s
1
=
⋯
=
s
k
-
1
=
0
,
s
k
≠
0
,then for smallε
~
≠
0 the system (32) can have the greatest number of limit cycles k in a small neighborhood of the homoclinic loop.In our case it is easy to computes
1
=
-
1
≠
0. Thus the equation s
0
=
M
(
0
)
=
0 gives the condition for which in system (32) a single limit cycle bifurcates from the homoclinic loop of the unperturbed system. Taking into account (32), (33), (35a), (35b), and (37) and having conducted some straightforward calculations, we obtain
(39)
M
(
0
)
=
∫
-
∞
+
∞
[
-
p
(
φ
0
,
ψ
0
)
ψ
˙
0
]
d
τ
=
∫
-
∞
+
∞
[
φ
0
-
(
1
3
)
b
~
φ
0
3
]
(
φ
0
-
φ
0
3
)
d
τ
=
(
1
45
)
(
-
60
+
48
b
~
)
=
0
.Hence(40)
b
~
=
-
σ
(
c
1
c
3
)
=
5
4
.Equation (40) expresses the condition for arising of two homoclinic orbits L
ε
+ and L
ε
- connecting the “saddle point” (
0
,
0
) to itself in system (32). It can be said, that under conditions (40), the perturbed Hamiltonian system (32) admits a single homoclinic limit cycle L
ε
~, consisting of the orbits L
ε
+ and L
ε
-; that is, L
ε
~
≡
L
ε
~
-
∪
L
ε
~
+. Let us remember that the perturbations in system (32) preserve the equilibrium point (
0
,
0
) unchanged. The limit cycle L
ε
~ is an isolated homoclinic loop located in a small neighborhood of the loop L
0 and passing through the point (
0
,
0
).Further, the stability of the limit cycleL
ε
~ is determined by the sign of the saddle point quantity [46]
(41)
σ
0
=
{
∂
∂
X
(
Y
+
ε
~
[
-
X
+
(
1
3
)
X
3
]
)
1
+
∂
∂
Y
(
X
-
X
3
)
}
|
(
X
,
Y
)
=
(
0
,
0
)
=
-
ε
~
<
0
.This means that the separatrix cyclesL
ε
+ and L
ε
-, as well as the limit cycle L
ε
~, are stable.Having in mind the derivation of(
C
.
6
) in Appendix C it is clear that the obtained here (40) precisely coincides with condition (
C
.
2
) derived in [30]. From this condition follows the relation between the critical value of velocity and the amplitude of the solitary wave solution related to the homoclinic loop given by (
C
.
3
). (Note that in (35a), (35b) a
0
=
2). However, as has been established in [30], by means of the soliton perturbation theory the relation between the equilibrium amplitude and velocity of the perturbed soliton solution (see (
C
.
4
)) in fact coincides with the same relation in (
C
.
3
). Consequently, the numerical confirmation of (40) and therefore of (
C
.
3
) coincides with the numerical confirmation of (
C
.
4
).


## 6. Discussion

Since (11) is derived from (1), we expect that the existence of this single limit cycle of (11) which arises around a closed phase trajectory of the unperturbed system will lead to the existence of periodical attractor of (1). As in this case the emerging limit cycle has a finite period, we expect that the corresponding periodic attractors of the generalized CCQGLE will represent the localized pulsating solutions of this equation. We should mention however that the dynamical system given by (11) possesses different equilibrium points depending on the values of the parameters c
j, we may expect the existence of a variety of limit cycles. Even in the case considered here the other limit cycles could be found in the region between homoclinic and heteroclinic trajectories, or in the neighborhood of a homoclinic or heteroclinic trajectory of the unperturbed system. It is clear that the question of how many limit cycles exist in (11) and the question of their type require further systematic mathematical investigation.Let us now discuss the question of numerical confirmation of the results obtained in Sections4 and 5. At this point we do not have any numerical confirmation of the results presented in Section 4. We have fixed the value of c
~
3
=
5
/
2, assuming, that c
1
=
-
1, c
3
=
5
/
2, c
5
=
-
1. As a result, we have obtained the following values of the parameters of (1): δ
=
-
0.371895, β
=
0.497519, ε
=
1.49504, μ
=
-
0.371895, ν
=
-
0.37375. The application of the Melnikov method leads to additional condition for b
~ given by (26). What we can say at this point is that the parameters β, ε, μ, ν are in regions for which stable localized fixed-shape solutions and localized pulsating solutions of (1) are identified in [14], while the values of δ and γ
=
1.005 are larger than the usual values [14, 16]. In order to confirm the limit cycles predicted by our approach in Section 4, we plan a further numerical investigation of the CCQGLE.In Section5 the situation is completely different. Having applied the Melnikov method to the equation of Duffing-Van der Pol oscillator we have proved the existence of a limit cycle that arises in a neighborhood of a homoclinic trajectory of the unperturbed system (see (40)). This result coincides with the condition (
C
.
2
) derived earlier in [30]. As the homoclinic trajectory is described by the solution of the unperturbed equation for infinitely long time, the condition of existence of the limit cycle in (40) gives the relation between the critical value of velocity and the amplitude of the solitary wave solution (see (
C
.
3
) as well as (17) in [30]). As has been mentioned in Section 5, there appears a problem in the numerical confirmation of (
C
.
3
) and (
C
.
4
).To prove the predicted relation between the critical value of velocity and the amplitude of the solitary wave solution(
C
.
3
) or of the equilibrium values of the amplitude and velocity of perturbed soliton solution (
C
.
4
), a detailed numerical investigation has been performed in [31]. Equation (1) with ν
=
μ
=
ε
=
0 has been numerically solved by means of the split-step Fourier method, applying the Blow-Wood RK4 scheme for the following values of parameters δ
⊂
(
5
·
1
0
-
3
÷
0.5
), β
⊂
(
1
·
1
0
-
3
÷
0.3
), and γ
⊂
(
1
·
1
0
-
4
÷
0.3
) [31]. The obtained results in [31] confirm with an excellent accuracy the relation between the equilibrium values of the amplitude and velocity of perturbed soliton solution (
C
.
4
) and therefore of (
C
.
3
) and finally (40).Taking into account the excellent numerical confirmation of results of Section 5 in [31], we could expect the numerical observation of the limit cycles for the strongly nonlinear Lienard-Van der Pol equation obtained for the first time here in Section 4.


## 7. Conclusion

We have studied the dynamics of the localized pulsating solutions of generalized cubic-quintic complex Ginzburg-Landau equation (CCQGLE) in the presence of intrapulse Raman scattering.The main result of this work is a proposal for an approach for identification of periodic attractors of the generalized CCQGLE. First, we use ansatz of the travelling wave and determine some conditions for the material parameters; then we derive the strongly nonlinear Lienard-Van der Pol equation for the amplitude of the nonlinear wave. Next, we apply the Melnikov method to this equation and we show that for a set of fixed material parameters a limit cycle arises around a closed phase trajectory of the unperturbed system. After that we prove its stability. Due to the complexity of the strongly nonlinear Lienard-Van der Pol equation, however, it is clear that the question of how many limit cycles exist as well as the question of their type will require further systematic mathematical investigation.We next have shown that the Melnikov method could be applied to the equation of Duffing-Van der Pol oscillator used for the investigation of the influence of the IRS on the bandwidth limited amplification [30]. We have proved the existence and stability of a limit cycle that arises in a neighborhood of the homoclinic trajectory of the corresponding unperturbed system. As this trajectory is described by the solution of the unperturbed equation for infinitely long time, the condition of existence of the limit cycle derived in (40) is equivalent to the relation between the critical value of velocity and the amplitude of the solitary wave solution found in [30] (see (17)).



---

*Source: 308947-2014-10-15.xml*