Abstract
The nonlinear dynamics of predator-prey systems coupled into network is an important issue in recent biological advances. In this research, we consider each node of the coupled network represents a discrete predator-prey system, and the network dynamics is investigated. By applying Jacobian matrix, center manifold theorem and bifurcation theorems, stability of fixed points, flip bifurcation and Neimark-Sacker bifurcation of the discrete predator-prey system are analyzed. Via the method of Lyapunov exponents, the nonchaos-chaos transition of the coupled network along the routes to chaos induced by bifurcations is determined. Numerical simulations are performed to demonstrate the bifurcations, various attractors and dynamic transitions of the coupled network. Via comparison, we find that the coupled network exhibits far richer and more complex behaviors than single predator-prey system, including period-doubling cascades in orbits of period-2, period-4, period-8, invariant closed curves, dynamic windows for periodic orbits and invariant curves, quasiperiodic orbits, tori, and chaotic sets. Moreover, the attractors of the coupled network show more diverse and complicated structures. These results may provide a new perspective on the predator-prey dynamics in complex networks.
1. Introduction
Complex networks have been widely studied in the field of interdisciplinary, involving nonlinear science, complexity science, ecology, biology, and many other disciplines. The complex networks are topological abstractions of a large number of real complex systems, as well as topological foundations for the analysis of complexity. Via analyzing the characteristics of complex networks, many research works have played an important role in understanding origination and degree of complexity [1, 2]. Moreover, lots of complex networks investigations have contributed to solving practical problems, attracting attentions in the fields of science and engineering [3].
Complex networks were originally within the scope of graph theory and can be traced back to the problem of the famous Seven Bridges of Konigsberg raised by Euler. The simplest network is regular network, which is an important basis of the complex networks. Three common regular networks studied in literature are globally coupled network, star network, and nearest-neighbor coupled network [4]. In 1950s, mathematicians Erdos and Renyi came up with an approach to construct random networks [5]. Regular networks and random networks are two types of extreme complex networks. Later, the researchers further constructed two important types of networks, small world network, and scale-free network [6, 7]. The small world networks have two statistical characteristics, namely, large family coefficients and small average distance [6]. In the scale-free network, node degree obeys power-law distribution [7]. The above networks reflect important characteristics of real complex networks and show great significance for exploring complex networks [8, 9].
When each node in a complex network represents a nonlinear dynamic system, the complex network becomes a very complex dynamic system which demonstrates complicated nonlinear characteristics. The research on such dynamic complex networks mainly focused upon their diverse dynamic behaviors. For example, synchronization dynamics is one of the most interesting topics in research of dynamic complex networks [10, 11]. Wang et al. investigated the synchronization of coupled heterogeneous complex networks and derived sufficient criteria for quasi network synchronization [10]. The complex network control problems are attracted as well [12, 13]. Wang and Chen studied the control of a scale-free dynamical network by applying local feedback injections to a fraction of network nodes [13]. The researchers also concerned many other interesting dynamics of complex networks, such as the emergence of avalanche dynamics and cascading failures.
Among the rich dynamics of complex networks, the transition of dynamics in coupled dynamical networks is a commonly concerned issue [4, 14]. Li et al. studied the transition to chaos in complex dynamical networks and verified that, during the transition from nonchaotic to chaotic states, if the network topology is more heterogeneous, the coupling strength required to achieve the transition can be decreased [14]. Zhang et al. researched the emergence of chaos in complex dynamical networks and found that the chaotic state will emerge without changing each node’s parameter if these nodes are connected through a certain type of network [4]. Moreover, a possible bifurcation route from chaos to hyperchaos was discovered in coupled dynamical networks [15]. How does network topology influence the dynamical transition of the complex networks has also become one focus of the research [14]. These approaches greatly improved the understanding of nonchaos-chaos transition in coupled dynamical networks.
In recent decades, mathematical biology has developed to be a hot research field, and mathematical models have been established to analyze various biological phenomena. Due to ubiquity and importance of the predator-prey relationship, the dynamics of predator-prey systems has attracted the attentions of many scholars [16–20]. However, the dynamics of coupled networks of predator-prey systems is still seldom investigated and documented in literature.
Early studies of predator-prey systems mainly utilized continuous models, via which the stability of equilibrium and existence and attractivity of limit cycles are often studied [21, 22]. When the populations have short life time and generations are nonoverlapping, it is more suitable applying discrete models to illustrate predator-prey systems [23]. Generally, discrete predator-prey models can often demonstrate more abundant dynamical behaviors than the continuous ones and therefore are particularly important for revealing nonlinear characteristics of discrete predator-prey systems [24, 25]. A great deal of research works manifested that the discrete predator-prey systems can exhibit various bifurcations and complex dynamic behaviors [16–20, 24, 25]. Through the research of He and Lai [16], Liu and Xiao [17], Jing and Yang [18], Agiza et al. [19], Hu et al. [20], Huang and Zhang [25], and etc., they found flip bifurcation and Neimark-Sacker bifurcation are basic bifurcations occurring in discrete predator-prey systems. The two bifurcations often start a route to chaos, revealing the dynamic transition from nonchaotic states to chaotic states. Moreover, complex dynamic behaviors can emerge on the routes to chaos, including cascades of period-doubling, periodic windows, quasiperiodic orbits and chaotic sets [17, 18, 25]. Via calculating the Lyapunov exponents and fractal dimension, Agiza et al. even determined diverse strange attractors of the predator-prey systems [19].
The former research works on discrete predator-prey systems normally studied the dynamics of single system, but rarely took into account of the coupled networks of predator-prey system. Since natural predator-prey systems are usually not isolated and closely connected with each other, the coupled networks of predator-prey systems can be more accurate to reflect the characteristics of real ecological systems and deserve further investigations. On the other hand, the coupled networks of predator-prey systems may exhibit more complex dynamic behaviors. An investigation on coupled networks may provide new phenomena of predator-prey dynamics.
In this research, we will study the bifurcations and dynamic transition between nonchaos and chaos in a coupled network of discrete predator-prey system. The exploration is arranged as follows. Section 2 will give the model for the coupled dynamical network of discrete predator-prey system. Section 3 will analyze the stability of fixed points, period-doubling bifurcation, and Neimark-Sacker bifurcation when the discrete predator-prey system is not coupled. In Section 4, we will determine the sufficient conditions for emergence of chaotic state in the coupled dynamical network. In Section 5, numerical simulations will be carried out to demonstrate the complex dynamic behaviors. And finally, Section 6 will provide the conclusions.
2. Coupled Dynamical Network of Discrete Predator-Prey System
In this research, the complex dynamical network is considered to be modeled on the basis of a discrete Leslie-Gower predator-prey model. Leslie-Gower model is an important extension of Lotka-Volterra model. The Leslie-Gower predator-prey model assumes that interacting species grow according to the logistic law and that the environmental carrying capacity for the predator is not a constant but proportional to the population size of the prey [26]. However, due to the rarity of the prey, the predator can switch over to other food, but its growth is still limited by the fact that its favorite prey is not available in abundance [26].
The predation relationship studied in the predator-prey system is considered as Crowley-Martin type functional response, which can accommodate interference among predators [27]. The Crowley-Martin type functional response is used for data sets that indicate feeding rate that is affected by predator density [28]. It is assumed that predator-feeding rate decreases by higher predator density even when prey density is high, and therefore the effects of predator interference in feeding rate remain important all the time whether an individual predator is handling or searching for a prey at a given instant of time [28].
With Crowley-Martin type functional response considered, the continuous form of the Leslie-Gower predator-prey model can be described as follows: in which and represent population densities of prey and predator at the time of t respectively; R1, R2, A, B, m, D, E, and K are positive constants; R1 and are the intrinsic growth rate and the carrying capacity of the prey; the constants , , and are the saturating Crowley-Martin type functional response parameters, where measures the magnitude of interference among prey, expresses the interference among the predator, and represents the maximum consumption rate; R2 describes the growth rate of predator, is the maximum value which per capita reduction rate of can attain, and measures the extent to which environment provides protection to the predator.
In order to reduce the parameter numbers in the predator-prey model, we change the variables and parameters as follows, . Then (1a) and (1b) can be transformed intowhere , , , , and .
Applying the forward Euler scheme to (2a) and (2b), the discrete-time Leslie-Gower predator-prey system can be expressed by the following equations:in which represents the sequence number of iterations and describes the time step size. It should be noticed that .
Consider a coupled dynamical network of linearly and diffusively coupled identical nodes, with a full and diagonal coupling, where each node is a two-dimensional dynamical system described by (3a) and (3b). Therefore, the state equations of the coupled dynamical network of discrete predator-prey system are specified by [4, 14]where and denote the sequence number of the nodes in the coupled dynamical network, represents the coupling strength of the network, and denotes the coupling configuration of the network. If there is a connection between node and node , then ; otherwise, . Simultaneously, let , , where is the degree of node and can be defined to be the number of its outreaching connections, which means that the network is fully connected in the sense of having no isolated clusters.
To investigate the dynamics of coupled dynamical network (4a) and (4b), we need to study the dynamics of one single node firstly, i.e., system (3a) and (3b). Via Jacobian determinant, bifurcation theorem, and center manifold theorem, stability analysis and bifurcation analysis are performed on system (3a) and (3b) to learn how the system dynamics responds to parameter variation. On the basis of investigations on system (3a) and (3b), then the dynamics of the coupled dynamical network is analyzed via applying the method of maximum Lyapunov exponent.
3. Stability Analysis and Bifurcation Analysis on the Discrete Predator-Prey System
3.1. Stability Analysis
For the convenience of analysis, we rewrite the expression of system (3a) and (3b) as the form of a map, i.e.,
According to the definition of a fixed point for the map, the fixed points of map (6) can be calculated by solving the following equations:
Direct calculation obtains four nonnegative fixed points, i.e.,where is a positive real root of the following cubic equation:
If the inequalities are satisfied then the map (6) has an unique positive equilibrium ,
The stability of the fixed points is determined by the method of Jacobian matrix. The Jacobian matrix associated with any point of map (6) is
The four fixed points , , , and are substituted into the Jacobian matrix (12), respectively, and the eigenvalues and are calculated. If and are found, the corresponding fixed point is stable; if or , an unstable fixed point emerges. The stability of the four fixed points can be determined as follows.(1)The two eigenvalues of are and . Obviously, we have and , then we know is unstable.(2)The two eigenvalues of are and . Hence, the conditions for that is stable can be determined as(3)The two eigenvalues of are and . Since , the fixed point is also unstable.(4)For , the corresponding Jacobian matrix is where Calculating the two eigenvalues of as where
By calculating and , the stable conditions for are determined as
3.2. Flip Bifurcation Analysis
Flip bifurcation is performed around the fixed point when parameter variation occurs in the discrete predator-prey system. Parameter is chosen as the bifurcation parameter. According to the flip bifurcation theorem [29], the flip bifurcation occurs when one of the two eigenvalues of is equal to −1. In such case, we have , and direct calculation yieldswhere
When condition (19) establishes, the other eigenvalue of is obtained as . For the occurrence of flip bifurcation, we need , namely,
Let and we can translate the fixed point to the origin. Simultaneously, via doing Taylor expansion near the fixed point, the map (6) can be rewritten as where , , , and are described in (15) andIn (24), all the τ values satisfy .
On the basis of map (23), we make a reversible transformation as follows: Then the map (23) is transformed to a norm form in which
Then we determine the central manifold of map (26) at the fixed point (0, 0, 0). By using the central manifold theorem [29], there is a central manifold existing at the fixed point (0,0,0), which can be approximated byIn (28), is assumed to beAccording to , we can get
Using the balance of (30) to determine the values of , , , and , we get the following results:
Consider the dynamics of map (26) restricted in the central manifold . Hence, a one-dimensional map is obtained, i.e.,
The coefficients in map (32) are expressed as follows:
If the map (32) experiences flip bifurcation, two determinant values, and , should be nonzero, i.e., It can be obtained by direct calculation,
According to the flip bifurcation theorem [29], the conditions of the flip bifurcation of the discrete predator-prey system are described as follows: if conditions (19) and (21) as well as (35) establish, the discrete predator-prey system undergoes flip bifurcation at the fixed point when the parameter varies in a small neighborhood of . Moreover, if , the period-2 orbit bifurcates from which is stable; if , the period-2 orbit bifurcates from which is unstable.
3.3. Neimark-Sacker Bifurcation Analysis
According to the Neimark-Sacker bifurcation theorem [29], the occurrence of Neimark-Sacker bifurcation first requires that the two eigenvalues in (16) are a pair of conjugate complex numbers whose module is one. This means and , which lead to
When parametric conditions satisfy (36a) and (36b), we translate the fixed point to the origin by the translation Applying the Taylor expansion to map (6), we obtain a new mapwhere the coefficients in map (38) are described in (24), but with . Under the conditions of (36a) and (36b), the two eigenvalues of the Jacobian matrix of map (38) at the fixed point (0, 0) are also a pair of conjugate complex with module as one. We here write these two eigenvalues aswhere , and explicitly, .
At the same time, the occurrence of Neimark-Sacker bifurcation needs that the derivation of to parameter is not zero at . With a calculation, we haveOn the other hand, to avoid that is a real number or a pure imaginary number, we also need that Condition (41) is equivalent to . Since , we know . Hence, condition (41) becomes , i.e.,
Then we study the norm form of map (38). With the invertible transformationthe map (38) can be transformed into the following new map:where
Calculate the two- and three-order derivatives of and at and . The following results are obtained:
In order for map (44) experiencing Neimark-Sacker bifurcation at , the following determinant value must not be zero:where By using the above formulas, condition (47) becomesin which
Taking the above calculations together, the occurrence of Neimark-Sacker bifurcation at the fixed point in the discrete predator-prey system needs the satisfaction of conditions (36a) and (36b), (42), and (49). Moreover, if and , then an attracting invariant closed curve bifurcates from the fixed point for ; otherwise, if and , an repelling invariant closed curve bifurcates from the fixed point for .
4. Analysis on Nonchaos-Chaos Transition in the Coupled Dynamical Network
For studying the dynamic transition between nonchaos and chaos in the coupled dynamical network, Lyapunov exponent is an important tool. As well known, the Lyapunov exponent can determine whether chaos exists or not. If the maximum Lyapunov exponent is positive, chaotic dynamics exists in nonlinear systems. In this research, the analyzing method as described in Zhang et al. [4] and Li et al. [14] is applied. Via studying the eigenvalues of coupling matrix of the coupled dynamics network, the positivity and negativity of the maximum Lyapunov exponent can be determined. According to the relationship between maximum Lyapunov exponent and chaos, the conditions for nonchaotic and chaotic network dynamics can be found.
Before analyzing the dynamics of the coupled dynamical network, we need firstly describe the properties of the coupling matrix . According to the research and description in Wu and Chua [30] and Zhang et al. [4], if the coupling matrix is an irreducible matrix and satisfies the conditions defined in Section 2, then has the following properties:
(1) 0 is the eigenvalue of the matrix Λ, and the corresponding eigenvector is , where means the transpose of a matrix.
(2) All other eigenvalues of the matrix Λ are real and strictly negative.
On the basis of the above properties, we define and rank the eigenvalues of the coupling matrix Λ as the following:
For convenience of analysis, we rewrite the equations of the coupled dynamical network (see (4a) and (4b)) in a vector form, in which and . Generally, the Lyapunov exponents of the coupled dynamical network can be calculated aswhere is the Jacobian matrix of the iterated map starting from a random initial state and is a set of orthonormal vectors in the tangent space of the map (52).
We denote the maximum Lyapunov exponent of the discrete predator-prey system (3a) and (3b) as , which can be described as [31]where , i.e., the functions defined in system (3a) and (3b). Choose parameter as varying parameter. Then from (54), we know that there exists a relationship between and when the value keeps in variation. We here express this relationship as the following: Generally, can be determined by numerical method.
Furthermore, according to Zhang et al. [4] and Li et al. [14], there is also a relation between and which can be described as where is the th eigenvalue of the coupling matrix Λ as described in (51). Because the coupling strength is positive, then the Lyapunov exponents of the coupled dynamical network can be ranked as follows:
The purpose here is to investigate the aggregative dynamics of the coupled dynamical network and determine the nonchaotic and chaotic regions with the variation parameter . In this research, we consider the case that the parameters of all coupling nodes are the same; i.e., the dynamics of the discrete predator-prey system at any node is attracted to the same attractor if the initial conditions are provided the same. Then, if the coupled dynamical network exhibits chaotic dynamics, there is at least one positive Lyapunov exponent in (57), implying . On this basis, we can classify the Lyapunov exponents of the coupled dynamical network into two types, and each type of Lyapunov exponents satisfies the following: Here, is a positive integer and is the largest negative Lyapunov exponent.
Combining (53)-(58a) and (58b), we can obtain the condition where at least one positive Lyapunov exponent emerges in the coupled dynamical network, i.e., Via solving inequality (59), the range of parameter τ in which the coupled dynamical network exhibits chaotic dynamics is determined. Therefore, with the increase of τ value, the distribution of nonchaotic and chaotic regions is known. Combining with the bifurcation analysis of the discrete predator-prey system, we can investigate the dynamic transition between nonchaos and chaos in the coupled network when flip bifurcation or Neimark-Sacker bifurcation take place. In the following numerical simulations, we will also demonstrate how the coupled dynamical network influences the bifurcations of the discrete predator-prey system.
5. Numerical Simulations
In this section, numerical simulations are performed to verify the previous theoretical results and exhibit complex and rich dynamical behaviors of the discrete predator-prey system and the coupled dynamical network. Bifurcation diagrams, phase portraits, and maximum Lyapunov exponents are plotted. Via comparing the dynamics of single predator-prey system and the coupled network, how the coupled network influences the predator-prey dynamics is investigated.
5.1. Bifurcations and Complex Dynamics of the Discrete Predator-Prey System
To exhibit the flip bifurcation of the discrete predator-prey system, we give a group of parameter values as , , , , and and choose parameter as bifurcation parameter. At the beginning, the occurrence of the flip bifurcation under such parametric conditions is demonstrated. The fixed point of the discrete predator-prey system is determined as (0.0125, 0.9897), and the critical point for flip bifurcation is calculated as . When , the fixed point of the system is stable. Taking for example, the two eigenvalues of the fixed point are , , the modulus of which are both less than one, suggesting is a stable node. When , we get the two eigenvalues of as , . Then map (26) restricted to the center manifold is got as . On the basis of the one-dimensional map , the two discriminatory quantities and are got as , . Hence, according to the flip bifurcation theorem, the discrete predator-prey system indeed undergoes flip bifurcation at . Furthermore, because , we know that a period-2 orbit bifurcating from is stable and the fixed point becomes unstable. Taking , for example, the two eigenvalues of are .0463, , respectively, therefore is a saddle point. Simultaneously, the discrete predator-prey system is attracted to the two stable periodic points (0.0128, 0.8841) and (0.0131, 1.0923).
Figure 1(a) shows the flip bifurcation diagram of the discrete predator-prey system with parameter ranging in under the above given parametric conditions. Explicitly, a period-doubling cascade is observed and the dynamic transition from periodic orbits to chaotic behaviors is found. Figures 1(b) and 1(c) are local amplifications of Figure 1(a), demonstrating the detail of periodic windows with period-6 and period-5 orbits. In the periodic windows, period-doubling cascade also occurs, resulting to period-10, period-20, period-40 orbits and period-12, period-24, period-48 orbits, etc. Moreover, other periodic windows can also be found in the chaotic region. Figure 2 shows the maximum Lyapunov exponents corresponding to the three bifurcation diagrams in Figure 1. It can be seen in Figure 2 that the discrete predator-prey system first goes into chaotic region at about .
(a)
(b)
(c)
(a)
(b)
(c)
Figure 3 shows different attractors in the phase portraits corresponding to the flip bifurcation diagram in Figure 1(a). Figures 3(a)–3(d) display the period-1, period-2, period-4, period-8 orbits. The doubling of periodic points in these orbits coincides with the period-doubling process in Figure 1(a). Figures 3(e), 3(g), and 3(i) show a cloud of random points, indicating chaotic behaviors. The maximum Lyapunov exponents of these dynamical behaviors are 0.0756, 0.1334, and 0.3265, respectively, confirming that these behaviors are indeed chaotic. Figures 3(f) and 3(h) display the period-6 and period-5 orbits; in fact, they correspond to the periodic windows as shown in Figures 1(b) and 1(c). The dynamic transition from period to chaos demonstrated in Figure 3 suggests a disordering process occurring in the discrete predator-prey system.
(a) τ = 0.2
(b) τ = 0.24
(c) τ = 0.26
(d) τ = 0.263
(e) τ = 0.267
(f) τ = 0.2703
(g) τ =0.2715
(h) τ = 0.2814
(i) τ = 0.295
Figures 4–6 demonstrate the dynamics of Neimark-Sacker bifurcation of the discrete predator-prey system. Before describing these figures, an example is provided to illustrate the occurrence of Neimark-Sacker bifurcation. Given parameter values as , , , , and , the fixed point is obtained as (0.8313, 0.2518), and the critical point for Neimark-Sacker bifurcation is determined as . When , taking , for instance, we get the two eigenvalues of are and ; therefore is a stable focus. When we take , the two eigenvalues are and . Besides, through calculation, we have and . It is verified that the discrete predator-prey system undergoes Neimark-Sacker bifurcation at . Moreover, on account of and , we have an attracting invariant closed curve bifurcating from the fixed point for .
(a)
(b)
(c)
(a)
(b)
(c)
(a) τ = 2.815
(b) τ = 2.9
(c) τ = 3.015
(d) τ = 3.037
(e) τ = 3.06
(f) τ = 3.07
(g) τ = 3.08
(h) τ = 3.15
(i) τ = 3.19
Figure 4 shows the Neimark-Sacker bifurcation diagram with parameter ranging in [2.8, 3.2]. Under the above given parametric conditions, when , a locally asymptotically stable fixed point dominates in the phase space, whereas when and the value of is at the nearby of , the system dynamics converges to invariant circles. Via amplifying the bifurcation diagram in Figure 4(a), periodic windows with period-26 and period-7 orbits are found, as shown in Figures 4(b) and 4(c). Moreover, as the value of parameter continuously increases, the discrete predator-prey system will present chaotic behaviors. As demonstrated by the maximum Lyapunov exponent in Figure 5, the occurrence of chaos begins at about . And in the chaotic region, periodic windows can repeatedly emerge (see Figures 5(b) and 5(c)).
Figure 6 shows a series of representative dynamical behaviors corresponding to the Neimark-Sacker bifurcation diagram. Figure 6(a) shows a locally asymptotically stable focus, Figures 6(b) and 6(d) display two different invariant circles, Figures 6(f), 6(h), and 6(i) all exhibit a cloud of random points, which suggest the emergence of chaos, and Figures 6(c), 6(e), and 6(g) demonstrate a few periodic windows occurring in-between the regions of invariant circles and chaos. The variation of these dynamical behaviors actually reveals the dynamic transition on the route to chaos induced by Neimark-Sacker bifurcation.
5.2. Bifurcations and Complex Dynamics of the Coupled Dynamical Network
To illustrate the nonchaos-chaos transition in the coupled dynamical network, here we present one example via applying the globally coupled network. Considering a globally coupled network with nodes, the coupling matrix of which can be described as
Through solving the solutions of , the eigenvalues of the coupling matrix are determined as and , respectively. According to (57), the maximum Lyapunov exponent of the coupled dynamical network is . When , the coupled dynamical network exhibits chaotic behaviors. This condition can be further converted via (59) to determine the range of parameter τ where chaotic dynamics emerges in the coupled dynamical network (4a) and (4b).
Figure 7 exhibits the flip bifurcation diagram of the coupled dynamical network when . Comparing with Figure 1, the bifurcation diagram here shows a few variations under the coupling effects of the discrete predator-prey system. First, the chaotic region extends. Following the calculation in previous section, the coupled network first experiences chaotic dynamics at about , which comes earlier than that in the single discrete predator-prey system. Likewise, the periodic windows in Figure 1 also shrink or even disappear in the bifurcation diagrams of the coupled dynamical network. Second, the dynamics at the bifurcation critical point is greatly changed. As shown in Figures 7 and 8(b), invariant circles emerge when the value of parameter stays around . Third, new periodic window may appear in the chaotic region. As demonstrated in Figures 7(b) and 8(f), the new periodic window shows the existence of period-8 orbit, the occurrence of which results from that maximum Lyapunov exponent of the discrete predator-prey system is very small, and therefore the network dynamics is not chaotic and stays at the periodic orbit of the discrete predator-prey system.
(a) δ = 0.01
(b) δ = 0.015
(c) δ = 0.02
(a) τ = 0.2
(b) τ = 0.2052
(c) τ = 0.252
(d) τ = 0.256
(e) τ = 0.262
(f) τ = 0.263
(g) τ = 0.265
(h) τ = 0.27
Figure 8 displays the predator-prey dynamics in the phase diagrams at one node of the coupled network. As demonstrated in the figures, it can be found that the coupled dynamical network shows more complex dynamical behaviors than single discrete predator-prey system. At the bifurcation critical point and around, complex invariant circles occur and replace periodic behaviors (Figure 8(b)). Like the discrete predator-prey system, the coupled network also experiences period-doubling process, generating period-2, period-4, period-8 orbits. When the network dynamics goes into period-8 orbit zone, chaotic dynamics comes sequently, exhibiting various chaotic attractors (Figures 8(c)-8(d) and 8(f)–8(h)). In the chaotic region, the network dynamics once shortly falls back to periodic behaviors (Figure 8(e)).
When the value of N is given as 20, the flip bifurcation diagram of the coupled dynamical network remains similar with that N = 5 (see Figure 9). However, with the increase of nodes, the phase diagrams of the coupled dynamical network may present more complex and interesting configurations (Figure 10). Around the bifurcation critical point , the phase diagrams of one node just show two simple invariant circles (Figures 10(a)-10(b)).
(a) δ = 0.001
(b) δ = 0.002
(c) δ = 0.003
(a) τ = 0.206
(b) τ = 0.207
(c) τ = 0.25
(d) τ = 0.254
When the discrete predator-prey system undergoes Neimark-Sacker bifurcation, the bifurcation dynamics of the coupled network will show another situation. Comparing with Figure 4, the change of the Neimark-Sacker bifurcation diagram under network coupling effect mainly reflects in the following aspects. First, the bifurcation critical point where dynamics transits from fixed point to attracted invariant circles occurs at a smaller value of parameter . In Figure 11, for , ; for , ; for , . It can be found that the stronger the coupling effect, the smaller the bifurcation critical point. Second, the τ value where the coupled network first goes into chaotic region comes earlier. Via the calculation in (59), the chaotic dynamics first takes place at about , 2.7918, and 2.7847 for the three cases in Figure 11. Third, even though the dynamics of coupled network goes into the chaotic region, periodic windows can also emerge. Moreover, in Figure 11, there exists another dynamic window, where the network dynamics returns to invariant circles.
(a)
(b)
(c)
The Neimark-Sacker bifurcation of the coupled dynamical network can induce rich and complex dynamical behaviors as demonstrated in Figures 12. When the value of parameter comes across the bifurcation critical point, the predator-prey dynamics at each node of the coupled network follows the behavior of invariant circle (Figure 12(a)). As the value grows, the network dynamics will suddenly become chaotic (Figure 12(b)). Following this chaotic region, the network dynamics returns to be nonchaotic, since the maximum Lyapunov exponent is smaller than zero. As shown in Figures 12(c)–12(g), different tori dominate in the phase diagrams. With further increasing value, dynamic windows of invariant circle and periodic orbit sequently take place (Figures 12(h)-12(i)). Finally, after the periodic window, the network dynamics again transits to chaos. Figures 12(j)–12(l) exhibit several chaotic attractors of the coupled dynamical network.
(a) τ= 2.805
(b) τ= 2.81
(c) τ = 2.815
(d) τ= 2.819
(e) τ= 2.82
(f) τ= 2.825
(g) τ= 2.84
(h) τ= 2.92
(i) τ = 3.015
(j) τ = 3.02
(k) τ = 3.05
(l) τ = 3.08
Our numerical analysis reveals that the dynamic transition exhibits similar trend under other coupling strengths when the coupled network experiences Neimark-Sacker bifurcation. Figure 13 demonstrates various network dynamical behaviors in the phase diagrams of one node. Like the nonchaos-chaos dynamic transition in Figure 12, the dynamics in Figure 13 sequently undergoes the alternation between invariant circles and tori, then the periodic window, and finally the chaos. However, the attractors of invariant circles, tori, and chaos exhibit different appearance with that in previous figures. This suggests the diversity and complexity of the predator-prey dynamics in the coupled network.
(a) τ = 2.72
(b) τ = 2.729
(c) τ = 2.73
(d) τ = 2.74
(e) τ = 2.755
(f) τ = 2.78
(g) τ = 2.785
(h) τ = 2.8
Figure 14 shows the Neimark-Sacker bifurcation diagram of the coupled dynamical network when the value of is given as 20. Compared with Figure 4, the bifurcation critical point and first chaotic dynamics point here both move slightly toward the origin. Compared with Figure 11, the chaotic region in Figure 14 shows no windows for invariant circle and period orbit. Generally, the dynamic transition along the increasing direction in the bifurcation diagrams presents similar but simpler trend. Three regions in these diagrams can be divided. In the region of straight line, the network dynamics is attracted to fixed point. In the region of bullet shape, invariant circles prevail. In the third region, the coupled network exhibits tori and chaotic dynamics. At the intersection between the second and the third regions, new attractors for invariant circles or tori can often appear, such as the examples given in Figure 15. These attractors are also a reflection of rich dynamics of the coupled network.
(a) δ = 0.0015
(b) δ = 0.002
(c) δ = 0.005
(a) δ = 0.002, τ =2.82
(b) δ = 0.002, τ =2.823
(c) δ = 0.005, τ =2.783
(d) δ = 0.005, τ =2.784
6. Conclusions
Generally, the dynamics of coupled network is far richer and more complex than single dynamic system considered without system coupling. Previous research has been made, revealing that even though the dynamic system at each node holds parametric conditions which generate nonchaotic dynamic behaviors, the coupled dynamical network can still exhibit chaotic dynamics [4, 14]. On the basis of former research works, this research investigates the complex dynamic behaviors and nonlinear characteristics in a coupled network of discrete predator-prey system. New improvements are obtained on the bifurcations and nonchaos-chaos transition of the coupled dynamical network. Via comparing the dynamics between single discrete predator-prey system and the coupled dynamical network, how the coupled network influences the predator-prey dynamics is determined.
The complex behaviors of the discrete predator-prey system are described firstly. The system may have four fixed points, among which the stable unique positive fixed point reflects the predator and prey can stably coexist. Moreover, the discrete predator-prey system can undergo flip bifurcation and Neimark-Sacker bifurcation at this fixed point. The two bifurcations both trigger a route to chaos, which is also a disordering process from period to chaos. On the routes to chaos, the system displays complex behaviors and various attractors, including period-doubling cascades in orbits of period-2, 4, 8, periodic windows with period-5, 6, 7, 10, 12, 20, 24, 26, 40, 48 orbits, invariant circles, and chaotic sets, which imply that the predator and prey can coexist in many possible complicated ways.
In comparing with the discrete predator-prey system, the coupled dynamical network can exhibit much more interesting behaviors. To illustrate the complex dynamics and nonchaos-chaos transition in the coupled network, globally coupled network is taken for investigation. The flip bifurcation and Neimark-Sacker bifurcation also emerge and induce dynamic transitions from fixed point to chaotic dynamics. Very complicated behaviors are found with the increase of varying parameter. With the findings of this research, the network dynamics can be concluded as follows:(1)Compared with the flip bifurcation of single discrete predator-prey system, the flip bifurcation of the coupled network also triggers period-doubling cascade with period-2, 4, 8 orbits. However, the network dynamics at the bifurcation critical point turns to be invariant circles. Furthermore, the coupled network goes into chaotic dynamics earlier and the chaotic region is extended. New periodic window (such as period-8 orbit) may appear in the chaotic region.(2)Compared with the Neimark-Sacker bifurcation of single discrete predator-prey system, the bifurcation critical point where dynamic behaviors transit from fixed point to attracted invariant circle and the point where the dynamics firstly becomes chaotic both comes earlier in the coupled dynamical network.(3)In the coupled network, the dynamic transition on the routes to chaos exhibits complex process. When flip bifurcation occurs, the network dynamics may sequently experience invariant circles, period-doubling cascade, chaotic sets, and periodic windows and finally stay at the chaotic region. When Neimark-Sacker bifurcation takes place, the coupled network may sequently undergo invariant circles, chaotic sets, tori, invariant circles, periodic orbits, and eventually chaotic dynamics.
Data Availability
The data of numerical results are generated during the study.
Conflicts of Interest
The authors declare no conflicts of interest.
Acknowledgments
This research was financed by the National Water Pollution Control and Treatment Science and Technology Major Project (No. 2017ZX07101-002, No. 2015ZX07204-007), the National Natural Science Foundation of China (no. 11802093), and the Fundamental Research Funds for the Central Universities (no. JB2017069).