Synchronized Cycles of Generalized Nicholson-Bailey Model

In this paper, we study a drive-response discrete-time dynamical system which has been coupled using convex functions and we introduce a synchronization threshold which is crucial for the synchronizing procedure. We provide one application of this type of coupling in synchronized cycles of a generalized Nicholson-Bailey model. This model demonstrates a rich cascade of complex dynamics from stable fixed point to periodic orbits, quasi periodic orbits and chaos. We explain how this way of coupling makes these two chaotic systems starting from very different initial conditions, quickly get synchronized. We investigate the qualitative behavior of GNB model and its synchronized model using time series analysis and its long time dynamics by the help of bifurcation diagram.


Introduction
Traditionally, the study of long time behavior of a dynamical systems was based on the examples of ordinary differential equations with regular solutions and those solutions which remained in a bounded region of the phase space could be divided into two different types based on their local behavior: first, a stable equilibrium point and second, a periodic (or quasi-periodic) oscillation. Edward Lorenz in 1961, by working on a simplified version of atmospheric transfer model which was consisting of three nonlinear ordinary differential equations, numerically observed that a very small changing in the initial conditions of the system equations makes a huge difference on the long term behavior of their solutions [1]. Indeed, his finding was due to one of the major properties of chaotic dynamical systems which later called as sensitive dependence on initial conditions or butterfly effect.
Chaos is a complex nonlinear phenomenon that has been increasingly studied in the last three decades. During those years, many fields of science and engineering have been affected by chaos studies. One of the most important achievements in nonlinear and complex dynamics is the discovery of synchronized chaos. Synchronization happens when two events take place in synchrony at the same time and when time approaches infinity, the error between solutions of the first system and its synchronized one vanishes and approaches to zero.
The synchronization between two dynamical systems is a well-known phenomena occurring in Physics, Biology or Engineering and refers to a phenomenon that may occur when two or more oscillators are coupled. For the first time, Christiaan Huygens discovered synchronization of coupled pendulum clocks in 1665 [2]. Occurrence of synchronization in coupled chaotic system composed of identical chaotic oscillators has been detected for the first time by Fujisaka and Yamada [3] and [4] and after that it has been reported by Pecora and Carroll [5].
The dynamics of coupled chaotic systems show properties which we cannot detect in the behavior of the individual elements [6]. Someone can find the same spatial synchronized fluctuations in biology, ecology and epidemiology [7] [8] [9] [10] [11]. Synchronization of complex population oscillations in natural systems has been examined widely by some researchers [12]. Bernd Blasius and Lewi Stone worked on a chaotic UPCA foodweb model and they claimed that the spatio-temporal structures associated with phase synchronization have important implications for conservation ecology. They proposed that even though perturbation of a local patch population can bring them to the brink of extinction, the periodicity of spatial phase synchronization can help to buffer the endangered population by colonizers. They also asserted that unlike this thought that population synchronization can cause global population extinction [13], however, phase synchronization can be useful for maintaining species persistence. Their finding indicated that synchronization can shape the distribution and abundance of species even in continental scale.
As we have already discussed, there are many varieties of synchronization. In this research, instead of exploring all of these different types of synchronization which have been proposed for different purposes and with various applications, we simply focus on the most fundamental case, and we will develop our new approach based on a basic mathematical concept. Indeed, the purpose of this paper is that after defining and setting the fundamental concepts which we need to establish the basis of our study, to demonstrate such configurations under a suitable coupling method is possible. Moreover, we explore the dynamical and ecological effects of synchronization of a host-parasitoid model which is a generalization of Nicholson-Bailey model (GNB). The GNB model demonstrates regular and irregular or chaotic oscillations. We define a lift function which is technically a convex function and maps the orbits of the drive system into the orbits of the response system. Using this convex function, we drive the response system which inherits all the complex qualitative dynamics of GNB model and mimics that certain properties of the motion which is shared between them. We investigate numerically that this method of coupling, synchronizes completely the stable and periodic cycles and even the chaotic motions of GNB model and to do that, we need to adjust the synchronization constant to be closer to zero. We demonstrate the complex dynamics of GNB model and its coupled system by conducting some time series and bifurcation analysis.

Drive-Response System Derivation
In this section, we derive the coupled system corresponding to the drive system by defining a convex link function. We consider the following drive system: and it is called response system. Now, for the following non-linear discrete-time dynamical system, we are going to develop a theorem which helps us to find the synchronization threshold. We consider the following drive-response system: where , n n n X Y ∈  are the state vectors of drive system (2.1) and response system (2.2) respectively, g is a mapping from n  to itself and a constant 0 1 s < ≤ . The Jacobian matrix for drive-response system (2.3) has the following form: Definition 2.2. We say that the drive system (2.1) and response system (2.2) are in complete synchronization if Here, we imply to a general definition in synchronization theory which is crucial for proposed coupling method.
where α is called a contraction constant of g on E [15]. Now, we have the following result for drive-response system (2.3). Proof. For drive-response system (2.3), we have:  To obtain our rigorous results for complete synchronization, we need to find the normal form for drive-response system (2.3), we need to perform a few linear coordinate transformations that will put (2.3) into a form which is easier to work with [16]. First we transform the fixed point ( ) X Y of the system (2.3) to the origin by the translation The Jacobian matrix for drive-response system (2.6) has the following form: T. Azizi, G. Kerr American Journal of Computational Mathematics Then, we split off the linear part of the system (2.6) and write ( ) Let Q be the matrix that transforms the matrix J  into (real) Jordan canonical form which has the following form Then, under the transformation  . We remark that the transformation (2.10) has simplified the linear part of (2.8) as much as possible. One can continue the task of simplifying the nonlinear part. However, for our purpose, we only need to focus on the linear part of the system. The schematic representation of the procedure of complete synchronization in a general discrete-time drive-response dynamical system has been demonstrated in Figure 1.  and Grobman (1959) proved that the orbit structure near a hyperbolic fixed point has the same qualitative structure of associated linearized system [17] [18] [19]. According to Hartman-Grobman theorem, the dynamical systems behave similar to their linearization part around the fixed point. However, this theorem needs the linearization part without eigenvalues with real part zero for continuous-time system and for discrete-time systems, it needs the absolute value of the eigenvalues of linearized part not become one. Remark 1. Because of the nature of the contraction mapping theorem and therefore the proposed coupling method, the Hartman-Grobman theorem can be applied directly into our problem. Lemma 2.5. Using Hartman-Grobman theorem, to find the condition under which the drive-response system (2.11) achieves complete synchronization, we only need to look at the linear part of (2.11) Theorem 2.6. Given Jacobian matrix 1 J Q JQ − =  , for which the following inequality holds: where α is the contraction constant, Ĵ ρ is the spectral radius of Ĵ and i λ for 1, , i n =  are the eigenvalues of Jacobian matrix Ĵ . Then, the mapping G is a contraction and the drive-response system (2.11) satisfies the complete synchronization properties, i.e. Proof. Using the contraction mapping theorem [20].
In dynamical system point of view, it is possible for a point arbitrarily close to fixed point ( ) of system (2.11) to generate an orbit which stays arbitrarily close to * Z . An orbit which could circle around the equilibrium * Z staying within the proposed bounded region, for initial conditions sufficiently close to * Z , could eventually approach * Z . In this case, * Z and all invariant set of points which demonstrating the same attractive property called attractor. Using this statement, we define an attractor of drive-response system (2.11): i.e., ( ) G Λ ⊆ Λ . We define the distance between Λ and a point Z, as , then Λ is called an attractor for drive-response system (2.11). For a stable fixed point, it is of special interest to determine the set of initial conditions whose subsequent orbits end up at this fixed point and we call this set as the basin of attraction of stable fixed point which achieves by the following definition:  a continuously differentiable map, then, the compact and invariant set Γ is called the attracting set of drive-response system (2.11) if there exists an open neighborhood B of Γ such that ( ) The largest such B is called the basin of attraction for system.

Stable Threshold for Synchronization in Discrete-Time Dynamical Systems
We continue this section by defining a new concept in synchronization theory of discrete-time dynamical systems. The concept of stability in this study is similar to the one that we have in contraction mapping theorem which is different from the relative stability of equilibrium point and some nominal motion. We say that a system is stable if the final state of the system is independent on initial conditions and we call a system is attracting if the orbits of that system get pulled in or converge towards each other [21]. In general, stability can be interpreted as a property of solutions of a dynamical system which means all solutions converge towards each other [22].

Synchronized Cycles in Generalized Nicholson Bailey (GNB) Model: Description of the Model
Without loss of generality, we assume b r = and follow the same model assumptions as Asheghi in 2014 [25]. Therefore, (3.2) can be written as the form: x n x n y n x n Now, we apply this coupling method on the system (3.4). Consider the following drive-response system:

Synchronized Cycles in GNB Model without Parasitoid
Here, we show that in system (3.4), when the parasitoid populations go extinct T. Azizi, G. Kerr : where p is a function of 1 x and 2 x .
We have demonstrated the time-series corresponding to the solution of the drive-response system (3.11) -(3.12) in Figure 2. As we can see, with increasing the growth rate r, the behavior of system changes from stable equilibrium point to periodic behavior and then to irregular and chaotic dynamics which was expected since the system (3.11) -(3.12) has the same form and so dynamics of Ricker model.
Also, we performed a one co-dimensional bifurcation analysis for system (3.11) -(3.12) with respect to growth rate r in Figure 3(a) to discover the long term behavior of the system and we have compared the solution of drive system (3.11) with the response system (3.12) by showing the error between the solutions in Figure 3(b). As we will discuss also in Section 4, when synchronization constant s is larger, the drive-response system cannot be synchronized completely. Moreover, we have shown the Lyapunov Exponent corresponding the drive-response system (3.11) -(3.12) in Figure 3(c) which is the best place to  investigate the stability or chaotic behavior of system (3.11) -(3.12). As we know, the negative Lyapunov Exponent implies to stable behavior and when it is positive, we expect to see chaotic behavior.

Numerical Simulations
In chaotic systems, it seems to not being possible to reproduce exactly the same initial conditions and parameters and force the orbits converge. However, in this section, we will numerically show that by using a sufficiently strong coupling method, we can change the track of the orbits to converge. Therefore, there exists a possible way to get a complete synchronization in chaotic systems whereas they have been coupled by a suitable coupling method.
In this section, we demonstrate some numerical simulation to describe the qualitative behavior of drive-response system (3.5) - (3.8). The orbits of the system (3.5) - (3.8) in chaotic regime can be considered as chaotic oscillations. Now we want to study the evolution of the dynamic variables 1 2 , x x corresponding to the host population of drive system and response system, and 1 2 , y y which are corresponding to the parasitoid population of drive system and response system respectively. All analysis and numerical simulations which have been conducted in this section, are expected to reveal the type of attractor from equilibrium point, periodic and quasi-periodic orbits, and chaotic attractors for which the dynamics will eventually settle down and remain forever.
In Figure 4, we can see the evolution of the attractors of drive system (red color) and response system (black color) when we change the growth rate r. As we see, as long as we are increasing r, the dynamics of the system change from the stable equilibrium point which loses stability and arises to a limit cycle. This figure demonstrates that the drive and response system with different initial conditions, and for the smaller value of the synchronization constant s, become completely synchronized.
T. Azizi, G. Kerr In Figure 5, we increase the value of the synchronization constant s and we notice that two systems keep the same qualitative behavior from stable equilibrium point to limit cycle but they are not completely synchronized.
In Figure 6, we demonstrated the evolution of host population for drive-response system (3.5) -(3.8) (The interested parameters values have been selected from the given bifurcation diagram). For fixed parameter values 0.5 s = , 40 a = , 10 k = , and different growth rate r, we can easily see how the time series for 1 2 , x x change from stable and periodic oscillations to chaotic motions. However, as we discussed before, for smaller synchronization constant s, both 1 x and its coupled 2 x are completely synchronized. With previous illustration about the chaotic synchronization, we conclude that using this method of coupling, we can expect to get a complete synchronization in chaotic regime for smaller synchronization constant s which this has been demonstrated for growth rate 3.8 r = in Figure 6.
However, when we compare Figure 6 with Figure 7, the coupling method which has been explained in Section 2, is successful when the synchronization constant s has smaller values and is closer to zero. With increasing the synchronization constant s, we noticed that two convex functions (3.9) and (3.10) which we introduced in Section 2 cannot make a complete synchronization between the drive and response systems (3.5) -(3.8). For different growth rate r, bifurcation diagrams for drive system (3.5), (3.6) and response system (3.7), (3.8), have been demonstrated for fixed synchronization constant 0.5 s = in Figure 8 and for fixed synchronization constant 0.95 s = in Figure 9. The one-co-dimensional bifurcation diagram helps us to know about the dependence of the drive-response dynamical systems (3.5) -(3.8) on the certain parameter which here is the growth rate r. Here, we are expecting to get completely synchronization for 0.5 s = ( Figure 8) and we do not expect to get a complete synchronization phase for 0.95 s = ( Figure 9).
As it can be easily seen in Figure 8 and Figure 9, the dynamics of host population for drive-response system (3.5) -(3.8), for some range of parameter values     , 10 k = , red (drive system) and blue (response system). Down: the error between the solutions of drive system and response system receptively. r and when ( ) ( ) 1 2 , 0,0 y y = (without parasitoid), is similar to classical bifurcation diagram of Ricker model where the routes to chaos happen through the cascade of period-doubling bifurcations and crisis corresponding to the extinction of parasitoid for drive and response system. However, for smaller values of growth rate r, the host and parasitoid in both systems can coexist through the periodic and quasi-periodic cycles of the Neimark-Sacker bifurcation of interior steady states. We also can easily observe that the sudden changes of attractors (crisis) happen frequently when we increase the parameter values of r.

Conclusion
In this paper, we developed a drive-response system by defining a convex continuous link function which maps the orbits of the drive system into the orbits of its coupled system and keeps the same qualitative dynamics. We found an appropriate normal form for drive-response system and we obtained the condi-tions under which the solutions of drive and response system become completely synchronized. We provided a new concept in chaos synchronization, called, synchronization threshold, which means that the solutions of drive and response system diverge from each other and lose the complete synchronization properties when they pass the threshold. Also, we studied a coupled discrete-time two-dimensional host-parasitoid model which is a generalization of famous Nicholson-Bailey model. One of our objectives in this paper was to investigate the rich dynamics of drive-response system (3.5) -(3.8) around its equilibrium points and achieving the chaos synchronization. We developed a new drive-response system by defining a convex continuous link function which maps the orbits of the drive system keeping the same qualitative properties such as stability and periodicity into the orbits of its coupled system. We observed that this coupling method can be successful for drive-response system (3.5) -(3.8) to get a complete synchronization when the synchronization constant has smaller values, closer to zero. Moreover, numerical verification is performed to show the existence of wide range of dynamics of drive-response system (3.5) -(3.8) around the positive equilibrium point. We also changed the values of synchronization constant s in its range between (0, 1) and we observed that the response system (3.7) -(3.8) for smaller values of synchronization constant s is completely synchronized with its original drive system (3.5) -(3.6) and when we increased the values of synchronization constant s, we noticed that the qualitative behaviors of both systems remain the same, however, we do not get a complete synchronization between the solutions of drive and response system (3.5) - (3.8). In chaotic regime, for larger values of synchronization constant s, closer to one, we could not get a complete synchronization. But, for smaller synchronization constant s, closer to zero, we have shown that two systems are in complete synchronization when the dynamic is chaotic.