Assessing the Non-Linear Dynamics of a Hopf–Langford Type System

In this paper, the non-linear dynamical behavior of a 3D autonomous dissipative system of Hopf–Langford type is investigated. Through the help of a mode transformation (as the system’s energy is included) it is shown that the 3D nonlinear system can be separated of two coupled subsystems in the master (drive)-slave (response) synchronization type. After that, based on the computing first and second Lyapunov values for master system, we have attempted to give a general framework (from bifurcation theory point of view) for understanding the structural stability and bifurcation behavior of original system. Moreover, a family of exact solutions of the master system is obtained and discussed. The effect of synchronization on the dynamic behavior of original system is also studied by numerical simulations.


Introduction
Non-linear chaotic behavior has been observed in systems of different nature as this motion is based on homoclinic (heteroclinic) structures which instability accompanied by local divergence and global contraction [1][2][3]. It is well-known that autonomous nonlinear differential system of the form .
where n ≥ 3, λ is the vector of parameters and f : R n → R n is a smooth vector function (i.e., continuously differentiable) in some domain Ω, can display a rich diversity of periodic, multiple periodic, chaotic and hyperchaotic flows dependent upon the specific values of one or more bifurcation (control) parameters [4,5]. The investigation of dynamical processes in coupled nonlinear systems is an interesting problem from both theoretical (mathematical) and applied (engineering) points of view. Phenomena such as stability in interacting subsystems can be observed in nature and science. Usually, that phenomenon is called synchronization [6,7]. There are four basic types of synchronization known: complete, generalized, phase and lag synchronization [8]. Phase synchronization is the phenomenon of the onset of balance between the phases of the subsystems state variables oscillations, which is caused by an onset of the energy balance.
A principal problem toward complete understanding of nonlinear interactions is to identify where in its phase space one dynamical system is structurally stable. For example, in a small neighborhood of a structurally stable Poincaré homoclinic orbit lie only saddletype periodic orbits. On the contrary, both structurally unstable and attractive periodic orbits, in addition to saddle ones, may exist near a structurally unstable homoclinic orbit [9]. Note that after Smale's works [10,11] these systems are said to be Morse-Smale systems.
The structural stability (roughness) investigation of steady state and of limit cycles or other types of trajectories is a main problem in bifurcation theory. It is well-known that there is critical dependence of the stability conditions of limit cycles on the stability conditions of its steady states. Based on classical works [12][13][14], it was defined that by knowing the sign of Lyapunov values, L k (λ i ), (also called focus values, Lyapunov quantities (coefficients)) we can efficiently study the structure of complicated nonlinear system trajectories. In other words, the type of: (1) stability loss of equilibrium and (2) winding/unwinding of system trajectories in small neighborhoods of equilibrium depend on the sign of Lyapunov value [15,16]. Furthermore, the (sequential) number of the last non zero Lyapunov value determines the number of possible limit cycles (stable/unstable).
In the present paper, we focus our study on the following system .
where µ, α, β, γ and l are positive control parameters. The infinite structure from ordinary differential equations of system (1) was originally introduced by Hopf [17] in order to describe a fluid's turbulence dynamics. Later, firstly in a private communication and after that in a paper [18], Langford introduced (1) for α = β = γ = 1 and l = 0. In our next considerations, when α = 1, β = 1 , γ = 1 and l = 0, we will call system (1)-a generalized Hopf-Langford system (GHL). After the works of Hassard et al. [19] and Nikolov et al. [20], even today, the non-linear dynamics of GHL system has been intensively investigated with the help of both theoretical and numerical approaches [21][22][23].
Here we investigate the previously unexplored parameter regions of the generalized Hopf-Langford system (GHL). A new qualitative picture of behavior near bifurcation points can be obtained when the GHL system (1) has the following modified form .
where E = 1 2 x 2 1 + x 2 2 + x 2 3 is the energy. Note that in Appendix B, we show how . E in (2) is obtained. Moreover, we qualitatively and numerically investigate the structural stability of (2) (which is equivalent to the original system (1)), whereas for this goal a specific version of bifurcation theory, based on the computing of Lyapunov values (not exponents), L k (λ i ), was used.
The equilibrium (steady state) points (FPs) of system (2) are found by equating the right-hand sides of (2) to zero (see Appendix D). Thus, we obtain that equilibrium points of the system are O 2,3,4 : second, third and fourth EPs.
Here we note that the polynomial function for x 3 in (4) (according to Cardan's method [24], p. 23, the substitution x 3 = z − m 1 3 is valid, i.e., z 3 + K 1 z + K 2 = 0) has three different real roots if where K 1 = − l . Moreover, according to [24], if: (i) Q > 0, then the equation for z (the third equation in (4), respectively) has one real root and a pair of imaginary roots and (ii) Q = 0, then for K 1 = K 2 = 0 all the roots of the same equation are zero and for K 1 = 0 , K 2 = 0 two real roots.
While dealing with modified system (2), the main simplification is that it can be separated of two coupled subsystems in the master (drive)-slave (response) synchronization type. As we can see, the master system has the following form .
Since it is obvious that the corresponding dynamic behavior of the original system (2) principally depends on the behavior of the master system [7,25], below we investigate only the bifurcation dynamic of this system. It is seen that master system describes the evolution in time of variable x 3 and energy E.
The paper consists of three sections apart from the Introduction. In Section 2 we present analytical results concerning the system (6). Section 3 gives an overview of numerical simulations, and finally in Section 4 we summarize and conclude our results.

Qualitative Analysis
Our objective here is to consider the local (near equilibrium points) behavior of system (6) which presents an autonomous nonlinear two-dimensional dynamical model. According to the general theory of ordinary differential equations [26], the equilibrium (steady state) values of system (6) are as those in (3) and (4)-the part for x 3 and E. In order to determine the type of these equilibrium points, we make the following substitutions into (6) x Hence, after accomplishing some transformations, system (6) (in local coordinates) is given through where The stability of equilibrium points (3) and (4) is defined by the following Routh-Hurwitz conditions We remark that the notations p, q and R in (10) and (11) are taken from [14], as they are the coefficients of the characteristic equation If (12) (i.e., the Routh-Hurwitz conditions) is not valid, then O 1 is: (i) always unstable or (ii) is degenerate (critical case). Note that in case (ii) (i.e., q = 0) the implicit function theorem may no longer be applied, and the existence of O 1 in a neighboring system is not necessarily guaranteed. For investigation of this critical case, the method of reduction to the center manifold and the method of normal forms can be used [3]. By definition, when condition (10) is not valid, the steady states (3) and (4) become unstable, as in this case according to [3,14] "soft" (reversible) or "hard" (un-reversible) stability loss takes place. To define the type of stability loss of steady states (3) and (4) it is necessary to calculate the so-called first Lyapunov value (L 1 (λ 0 )) on the boundary of stability R = 0. For instance, in [3], Sections 9.2 and 9.3 and [27] the analytical results (in the form of Theorems) for stability and stability loss regarding limit cycle are shown. In the case of two first order differential equations, this value can by determined analytically by the formula in [14] (for details see Appendix C). For system (8) we have Thus, in this case, L 1 (λ 0 ) has the form: It can be seen that: √ q > 0 and (iii) if l = 0, x 3 = 0 then L 1 can be positive, negative or equal to zero. In accordance with the Lyapunov-Andronov theory, if L 1 = 0, then the second one (L 2 ) must be calculated. Notice that sign of L 2 determines the stability/instability of the external limit cycle. For L 2 > 0, the external limit cycle is unstable. On the other hand, for L 2 < 0, the external limit cycle is stable-for a detailed discussion see Appendix A or [3,[12][13][14]28].
Let us now calculate L 2 , when L 1 = R = 0. We do the calculation in the above case (iii) where we can use formula (A3) in the Appendix A. The reader can easily complete the details or find them in the following literature- [3,12,14]. In this case, we can introduce b 04 = l and (13) for system (8), and the expression of second Lyapunov value becomes: If we consider (15), it follows that in dependence of the system's parameters L 2 can have a different sign or to be zero.
In order to simplify our numerical calculations for L 1 and L 2 , some of the parameters are kept constant. The parameter values are: α = 1 and β = 0.9. The first Lyapunov value obtained for system (6) as a function of bifurcation parameters µ ∈ [0.53, 0.533] and l ∈ [0.208, 0.24] is plotted in Figure 1. It is seen that for smaller values of µ and larger values of l, that L 1 has negative values, i.e., a soft stability loss takes place. On the contrary, the first Lyapunov value is positive (a hard stability loss occurs) for larger values of µ and smaller values of l. Indeed, using this fact, we conclude that for some combinations of parameter's values L 1 can be zero. Thus, when L 1 = 0, we are interested in knowing the sign of the second Lyapunov value.
In Figure 2, L 2 (calculated at R = L 1 = 0 in (15)) is shown for different values of the bifurcation parameters µ and l, i.e., µ ∈ [0.53, 0.533] and l ∈ [0.208, 0.24]. It can be seen that L 2 is always negative. Therefore, following the theory introduced in [12,14] we have soft stability loss for the external limit cycle. Hence, one of the situations shown in Figures A3 and A4 (see Appendix A) is valid.
It will be shown in Section 3 that the situation in Figure A4 is valid for system (6).

A Family of Exact Solutions
Let us return to system (6) in which x 3 and E are denoted by u and v, respectively, that is Solving the first one of the above equations with respect to v one obtains Then, substituting expressions (17) into the second equation of Equation (16) one arrives at the following second-order nonlinear ordinary differential equation for the function u(t), namely, where Equation (18) is of Liénard type [29] (see also, e.g., [30,31] and the references therein). Equations of this family describe the behavior of a large number of physical systems. For instance, the famous classical and generalized force-free Duffing and van del Pol oscillator equations, the so-called modified Painlevé-Ince equations and many other important equations belong to this family (see, e.g., [32][33][34]).
Fortunately, it turned out that Equation (18) in which the functions g(u) and h(u) are given by expressions (19) can be factorized, without any restriction on the involved parameters α, µ and l, in the following explicit form where This result, which can be easily verified by direct computation, opens a way of constructing a family of particular solutions for the regarded Equation (18) since each solution of the simpler equation is, evidently, a particular solution of Equation (18) as well. On account of the expression (21) for the function φ 1 (u), Equation (23) above takes the form of the following generalized Riccati equation (cf. [35], Sec. 2.15) whose general solution reads where t 0 and λ are arbitrary constants.
In this way, substituting Equation (25) into the first equation of Equation (17), we obtain a two-parameter family of exact solutions of the regarded system (16), namely, Actually, in the (u, v)-plane each such solution lies on the parabola v = 1 2 u 2 , see Figure 3. The right-hand sides of formulas (26) simply show how this parabola is traversed over time. They imply: (i) when λ > 0 the function u is positive for each t ∈ (−∞, +∞) and tends to zero or µ when the time tends to −∞ or +∞, respectively; (ii) if λ < 0 is such that λ + e (t 1 +t 0 )µ = 0, then the function u is negative when t ∈ (−∞, t 1 ) and tends to zero or −∞ when the time tends to −∞ or t 1 (from below), respectively, while, when t ∈ (t 1 , +∞) the function u is positive then we have the stationary solution u = µ, v = 1 2 µ 2 .  (26). Here, (0, 0) and µ, 1 2 µ 2 are equilibrium points of this system.
In Figures 4 and 5, we numerically demonstrate that when master subsystem (6) is stable (i.e., the energy E of the system tends to a constant) or has a stable limit cycle (i.e., L 1 < 0) then the slave subsystem (x 1 , x 2 ) into (2) has periodic solutions with period one, or quasi-periodic solutions (oscillations) with period different to one. Note here that those data are in accordance with the results obtained in our analytical study.  At µ = 0.533 and l = 0.21907, another qualitative behavior of master (x 3 , E) subsystem and slave (x 1 , x 2 ) subsystem can be observed-see Figure 6. Here, subsystem (x 3 , E) is (structurally) unstable, whereas subsystem (x 1 , x 2 ) has chaotic solutions. Remarkably, this is valid only for x 3 > 0 and energy E has a still smaller amplitude than one borderline case-limit value. It should be remarked here that, in a strict sense, the energy of original system (2) is a major regulator of its qualitative behavior. For clarity, we briefly discuss the dynamics of master subsystem (6) for different values of either L 1 , L 2 , or both. If L 1 < 0, then reversible (soft) stability loss takes place, and (6) is structurally stable-see Figures 4 and 5. In the case of L 1 > 0, the behavior of the system is expected to be structurally unstable-see Figure 6. Finally, if L 1 = 0, then for L 2 < 0, system (6) is structurally stable similar to the cases in Figures 4 and 5-not shown here.

Conclusions
Assume a three-dimensional autonomous nonlinear dissipative system of Hopf-Langford type. The introduction of energy as a new variable in the original system then leads to an equivalent four-dimensional system which can be separated into two coupled subsystems in the master (drive) (x 3 , E)slave (response) (x 1 , x 2 ) synchronization type.
Here we qualitatively and numerically investigated the nonlinear behavior (such as effect of synchronization) of an equivalent 4D dissipative system. As our study has shown, if the master system (energy) solutions: (1) tend to a constant, then the slave system (x 1 , x 2 ) has periodic solutions; (2) oscillate, then the slave system has quasi-periodic behavior and (3) are (structurally) unstable, then the slave system has chaotic solutions. Our results for the first Lyapunov value, L 1 and second Lyapunov value, L 2 , (when R = L 1 = 0) presented in Section 2 suggest that system (2) is structurally stable for some intervals of its parameters, and therefore soft stability loss takes place. It is also evident from the analytical calculations that master system (6) can be presented as the famous Liénard type. Hence, a family of exact solutions is obtained.
In Figures A1-A4, the possible types of qualitative behavior of system (Equation (A1)) are shown, when L 2 (λ 0 ) = 0. If L 2 (λ 0 ) > 0 (at R = L 1 = 0), then in the neighborhood of the fixed point the following structures are valid: (i) an unstable limit cycle exists, and the fixed point is stable focus-see Figure A2a (ii) an external unstable limit cycle and a stable (with small amplitude) limit cycle occur, when the safe boundary of stability is crossed. In this case, if the initial perturbations are smaller from the amplitude of the external cycle, then the system will be stay in the neighborhood of the fixed point/stable limit cycle-see Figure A2b; and (iii) irreversible behavior, as after very small perturbation the system leave forever vicinity of the fixed point-see Figure A2c.   Figure A1. Legend: solid circle line-stable limit cycle; dashed circle line-unstable limit cycle. (a) an unstable limit cycle exists, and the fixed point is stable focus; (b) an external unstable limit cycle and a stable limit cycle occur; (c) irreversible behavior, as after very small perturbation the system leave forever vicinity of the fixed point.
In Figure A3, the qualitative behavior of system (Equation (A1)) is shown, when L 2 (λ 0 ) < 0 (at R = L 1 = 0). For small variation of the system's parameters, in vicinity of a fixed point three different structures are possible-see Figure A4. In this case, the system (Equation (A1)) will stay in the neighborhood of fixed point independently from the perturbation size, and the external cycle is stable. Figure A3. Vicinity of a point from boundary of stability R = 0, in which first Lyapunov value L 1 is zero, and second Lyapunov value L 2 is negative. Legend: U-unstable zone (R < 0); S-stable zone (R > 0); P-particular zone. Figure A4. Qualitative structures in the vicinity of a fixed point, when second Lyapunov value L 2 is negative. Note that different zones U, S and P are those in Figure A3. Legend: solid circle line-stable limit cycle; dashed circle line-unstable limit cycle. (a) structure in zone S; (b) structure in zone P; (c) structure in zone U.

Appendix B
Derivation of .