Discrete Dynamics in Nature and Society
Volume 2018 (2018), Article ID 8713651, 17 pages
https://doi.org/10.1155/2018/8713651

Research Article

Neimark-Sacker-Turing Instability and Pattern Formation in a Spatiotemporal Discrete Predator-Prey System with Allee Effect

Huayong Zhang, Xuebing Cong, Tousheng Huang, Shengnan Ma, and Ge Pan

Research Center for Engineering Ecology and Nonlinear Science, North China Electric Power University, China

Correspondence should be addressed to Huayong Zhang; rceens@ncepu.edu.cn

Received 20 January 2018; Accepted 19 April 2018; Published 26 June 2018

Academic Editor: Guang Zhang

Copyright © 2018 Huayong Zhang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

A spatiotemporal discrete predator-prey system with Allee effect is investigated to learn its Neimark-Sacker-Turing instability and pattern formation. Based on the occurrence of stable homogeneous stationary states, conditions for Neimark-Sacker bifurcation and Turing instability are determined. Numerical simulations reveal that Neimark-Sacker bifurcation triggers a route to chaos, with the emergence of invariant closed curves, periodic orbits, and chaotic attractors. The occurrence of Turing instability on these three typical dynamical behaviors leads to the formation of heterogeneous patterns. Under the effects of Neimark-Sacker-Turing instability, pattern evolution process is sensitive to tiny changes of initial conditions, suggesting the occurrence of spatiotemporal chaos. With application of deterministic initial conditions, transient symmetrical patterns are observed, demonstrating that ordered structures can exist in chaotic processes. Moreover, when local kinetics of the system goes further on the route to chaos, the speed of symmetry breaking becomes faster, leading to more fragmented and more disordered patterns at the same evolution time. The rich spatiotemporal complexity provides new comprehension on predator-prey coexistence in the ways of spatiotemporal chaos.

1. Introduction

Since the pioneering work of Lotka and Volterra on population dynamics [1, 2], predator-prey relationship has become a central topic for biological studies over the last one hundred years because of its significance in many problems [3]. To study the interactions between predator and prey in space, it is important to determine the specific forms of predator-prey functional response and spatial movement of the predator and prey populations [3]. Theoretically, there are a fair amount of works by using nonlinear reaction-diffusion equations to study the spatiotemporal dynamics of predator-prey systems [4, 5]. Many former studies demonstrated that spatial pattern formation is an important nonlinear phenomenon in the spatially extended predator-prey systems [6, 7]. Spatial pattern is a ubiquitous scenario in nature and often takes the role of modifying temporal dynamics and regulating population stability [8]. After the investigations in about four decades, pattern formation of predator-prey systems has been widely concerned and become one of timely central topics in biological science [9, 10]. It has been a great interest to understand spatial patterns and related mechanisms of interacting species in predator-prey systems.

The researchers have developed various reaction-diffusion predator-prey models in order to study the pattern formation of predator-prey systems in varying circumstances, such as different functional responses, growth forms, and population motions. In literature, it is found that many types of functional responses can contribute greatly to the spatiotemporal complexity of predator-prey systems, such as Beddington-DeAngelis type, Leslie-Gower type, Holling type IV, and ratio-dependent type functional response. Upadhyay et al. discovered that zooplankton-phytoplankton systems with Holling type IV functional response have a common mechanism as wave of chaos, which can result in the self-organization of complex patterns [11]. Banerjee and Petrovskii demonstrated that a ratio-dependent predator-prey system can develop patterns both inside and outside of the Turing parameter domain and Hopf bifurcation is essential for the onset of spatiotemporal chaos [6].

Different growth forms of prey can also exert influences on pattern formation. Generally, the most common growth forms include logistic form, exponential form, and the form with Allee effect. Many works have considered continuous-time models for the predator-prey systems with Allee effect. Numerical results in Rao and Kang [3] indicated that the spatial Michaelis-Menten type predation model with Allee effects exhibits rich dynamical transitions from wave growth not only to striped-like/spots patterns, but also to regular patterns, which implied that synergistic effects of diffusion and Allee effects can contribute to the patterns of the solutions in different ways. The dynamics of a reaction-diffusion predator-prey system with strong Allee effect in the prey population has been considered in Wang et al.’s research [12], which revealed that the impact of Allee effect essentially increases the system spatiotemporal complexity. More abundant investigations and results relating to the Allee effect influencing predator-prey pattern formation can be found in Petrovskii et al. [13], Sun et al. [14], and Rodrigues et al. [15].

Moreover, many researchers further found that the population motions play a key role in spatiotemporal pattern formation of predator-prey systems. As important natural phenomena, cross-diffusion and convection of the predator and prey are often considered and studied in spatially extended predator-prey systems. Ghorai and Poria found that the cross-diffusion supports the formation of a wide variety of spatial and spatiotemporal patterns [16]. Ling et al. explored the impacts of cross-diffusion on the formation of spatial patterns in a ratio-dependent predator-prey system with zero-flux boundary conditions and determined that the cross-diffusion can trigger the emergence of spatial patterns which is however impossible under the same conditions when cross-diffusion is absent [17]. In the predator-prey systems with predator and/or prey convection, the spatial patterns are often characterized by stripes [18].

The nonlinear mechanism of pattern formation, Turing instability, can date back to the pioneering work of Turing who studied a reaction-diffusion model which included two chemical species, an activator and an inhibitor, in 1952 [19]. Levin and Segel originally extended the mechanism of Turing instability to ecology and theoretically investigate how diffusion leads to emergence of spatial patterns [9]. Their method becomes a routine framework in the studying of the spatial pattern formation of predator-prey systems. Much research effort has been devoted to the mechanisms behind the emergence of spatial patterns and a variety of pattern formation mechanisms have been found. These mechanisms can be basically classified into two categories. One category is that spatial patterns emerge from locally stable kinetic dynamics due to diffusion and/or convection, often called pure Turing instability, whereas the other category is that spatial patterns emerge from locally unstable kinetic dynamics (such as limit cycle and chaos) through spatially heterogeneous perturbations, and one of the mechanisms is known as Hopf-Turing instability [20].

Most of the reaction-diffusion predator-prey models developed in literature are time-and space-continuous. Recently, a new type of mathematical model was introduced into biology to study the pattern formation, i.e., coupled map lattices (CMLs). In the 1980s, the CML was used by Kaneko to explore the spatial-temporal structure of coupled logistic map, he found very complex spatiotemporal dynamics, including frozen chaos, defect turbulence, spatiotemporal intermittency, and fully developed spatiotemporal chaos [21–23]. Recently in 2012, Mistro et al. [24], Rodrigues et al. [15], and Punithan et al. [25] demonstrated that the CML model shows importance in understanding the predator-prey dynamics on a fragmented habitat, which can exhibit spatiotemporal chaos and spatiotemporal multistability. Huang and Zhang further revealed that considering the combining effect of pattern formation and bifurcations in the CML can result in rich and complicated dynamical behaviors and nonlinear characteristics [26]. Diverse predator-prey pattern self-organizations and complex pattern transitions can be induced by the coupled effects of Turing instability and Neimark-Sacker instability [10, 27, 28].

On the basis of former research works, we aim at studying the complex dynamics of a space- and time-discrete predator-prey system under the influence of Allee effect with the application of CML. For the discrete predator-prey system considered, four properties should be specified. First, the habitat where predator and prey populations dwell is distinctly patchy or even fragmented [24], and the individuals migrate among different patches. Second, the births and deaths of populations all occur at discrete times, or within certain intervals of time [29]. Third, the population numbers are small and in this case the application of CML model may be the most appropriate. Fourth, the Allee effect occurring in the predator-prey system represents a positive correlation between the population density and the per capita growth rate in population dynamics. It deserves to explore the spatiotemporal structures of such predator-prey system under the combining effect of Neimark-Sacker instability and Turing instability, which can help to further understand the spatiotemporal complexity of the discrete predator-prey system and are still few investigated in literature. The structure of the present research work is arranged as follows. Section 2 describes the CML model of the discrete predator-prey system with Allee effect and presents the stability analysis. In Section 3, we perform Neimark-Sacker bifurcation analysis and Turing instability analysis. Section 4 provides the numerical results of the complex dynamical behaviors, Neimark-Sacker bifurcation and pattern formation. In Section 5, we provide the discussion and conclusions of the research.

2. CML Model and Stability Analysis

2.1. CML Model of the Spatiotemporal Discrete Predator-Prey System

In this research, we investigate a predator-prey system with Allee effect and Michaelis-Menten type functional response. Allee effect, which was introduced by as well as named after Allee, refers to a positive correlation between population density and the per capita growth rate [3, 30]. Recently, Allee effect has been extensively explored to investigate its impact on dynamical features of predator-prey systems [12, 13, 15]. According to the research of Rao and Kang [3], the Michaelis-Menten type predator-prey system can exhibit great spatiotemporal complexity under the influence of Allee effect. Therefore, it holds significance and importance to further investigate the spatiotemporal dynamics of the discrete form of such predator-prey system, which may show more complicated nonlinear characteristics.

The Michaelis-Menten type reaction-diffusion predator-prey system with Allee effect in prey can be described as follows [3]:

where and are population densities of the prey and predator at time , respectively, is the maximum per capita growth rate of the prey, is the prey carrying capacity, and () is the “Allee threshold”. The value of is considered as a measure of the intensity of the Allee effect: the less the value of is, the less prominent is the Allee effect [13] (Petrovskii et al., 2002). More specifically, the Allee effect is called “strong” if (when the growth rate becomes negative for ), while it is called “weak” if , and no Allee effect exists when [31] (Morozov et al., 2004). is the predation rate, is the food utilization coefficient, and is the mortality rate of the predator.

Considering dimensionless variables with the following scaling:

and simultaneously employing the population diffusive terms on (1), then we obtain the following reaction-diffusion equations:

in which

In the above equations, is the relative intrinsic growth rate of the prey; describes the Allee threshold; represents the energy conversion rate from prey to predator; denotes the relative death rate of predator; denotes the Laplacian operator and describes the population diffusion in space; nonnegative constants and are the diffusion coefficients of and , respectively, measuring the speed of individual movements.

The CML model of spatiotemporal discrete predator-prey system is developed based on discretizing (3). For making the discretization, we give a time interval , a space interval , and a two-dimensional rectangular domain which includes lattices (the length of each lattice is ). Each lattice represents one site and is ascribed to two numbers, i.e., the prey density and the predator density. Correspondingly, two state variables, and ( and ), are newly defined in the discrete space and time, representing the prey density and the predator density in the site and at th iteration. If we give an initial time , then the time at th iteration is . According to the framework of coupled map lattice described in previous approach [10, 15, 24–28], the dynamics of the spatiotemporal discrete predator-prey at each discrete step from to iteration consists of two distinctly different stages, (a) the dispersal stage and (b) the “reaction” stage. The dispersal stage can be obtained by discretizing the spatial terms of (3), described by the following equations:

where and are the prey and predator densities after dispersal; denotes the discrete form of the Laplacian operator, i.e.,

Via discretizing the nonspatial part of (3), the predator-prey reaction stage described in the CML model is expressed by the following equations:

where and are the reaction functions determined by local inter- and intraspecific interactions, given by

Equations (5)–(8) describe the governing equations of the CML model of the spatiotemporal discrete predator-prey system. For applying the CML model, periodic boundary conditions are utilized in this research [10].

2.2. Stability Analysis on Homogeneous Stationary States

A homogeneous stationary state of the spatiotemporal discrete predator-prey system is a state where the prey and predator densities are homogeneous in space and simultaneously keep steady in time. Its stability determines whether the predator-prey system can stably stay at the corresponding homogeneous stationary state under some external disturbances. Generally, the stationary states of the discrete predator-prey system satisfy the following conditions: i.e.,

for all of , , and . Under such conditions, the predator-prey equations can be described by the following dynamic equations:

For convenience of the following analysis, (10) are rewritten in the form of a map,

Mathematically, the homogeneous stationary states of the discrete predator-prey system are exactly equivalent to the fixed points of the map (11). Therefore, the results for solving the fixed points of (11) and the corresponding stability analysis are consistent with that on the homogeneous stationary states. According to the definition of the fixed point of a map, direct calculation yields five nonnegative fixed points for map (11),

in which

When and are established, the two fixed points and are positive. Since we intend to study the predator-prey coexistence in this research, hereafter we always presume these two conditions of positive fixed points are satisfied.

In order to determine the stability of the fixed points, the method of Jacobian matrix is often applied. The Jacobian matrix associated with map (11) at any point is described as follows:

Substitute the value of each fixed point into (14) and then calculate the two eigenvalues of corresponding Jacobian matrix. We denote these two eigenvalues as and . If and , then the corresponding fixed point is either a stable node or a stable focus; if and , then the corresponding fixed point is either an unstable node or an unstable focus; if and or and , then the corresponding fixed point is a saddle point, which is also unstable.

Substituting into (14), then we obtain the Jacobian matrix as

Obviously the two eigenvalues of are , . If we want the fixed point to be stable, the absolute values of and must be less than one. According to this, the conditions for that is stable are determined as .

Similarly for the fixed points and , we have

It can be found that one of the eigenvalues of (16) is the same, described as . Because we presume that , we directly have , suggesting and are unstable.

As for the fixed points and , we can write their Jacobian matrix as

For convenience, we demote the four elements of matrix (17) as , , , and . The two eigenvalues of (17) are obtained as

in which

According to and , the conditions for that or is stable are determined as and . With the given parametric conditions, the stability of and will be further determined in the numerical simulations section. In the following theoretical analysis, we denote the stable coexistent fixed point as .

3. Analysis on Neimark-Sacker Bifurcation and Turing Instability

3.1. Neimark-Sacker Bifurcation Analysis

Besides the homogeneous stationary states, the discrete predator-prey can also generate many other complex spatially homogeneous states, which always keep oscillating with time. To determine the transition from homogeneous stationary states to homogeneous oscillating states, bifurcation analysis is one of the most reliable and efficient ways. It should be noticed that such transition in the spatiotemporal discrete predator-prey system is mathematically equivalent to the transition from the fixed point to other attractors in the map (11). Generally, discrete system may exhibit various types of bifurcations, including fold bifurcation, flip bifurcation, and Neimark-Sacker bifurcation [29]. The Neimark-Sacker bifurcation is the discrete analogue of the Hopf bifurcation that occurs in continuous systems and is highly relevant to the discrete dual congestion control algorithm [32]. Often in literature, the explicit Neimark-Sacker bifurcation is also called Hopf bifurcation for map, the calculations of which can be seen in Cheng and Cao [29].

When map (11) undergoes Neimark-Sacker bifurcation, the stable fixed point will change to an invariant closed curve. In the following bifurcation analysis, parameter is chosen as the bifurcation parameter. According to the Neimark-Sacker bifurcation theorem, the occurrence of such bifurcation needs the satisfaction of many conditions [33]. Firstly, Neimark-Sacker bifurcation takes place if the eigenvalues in (18) are conjugate complex numbers and their moduli are one. That results in

Calculating on (19) yields the following conditions:

Under the satisfaction of conditions (21), the fixed point of map (11) is then translated to the origin by the following translation:

Thus, the map (11) is transformed into a new map as

where and describes a function with order at least four in the variables , and the coefficients in (23) are described as follows:

The map (23) has a fixed point (0, 0) and the two eigenvalues at this fixed point can be written as

where , and are described in (19) with . Based on conditions (21), the moduli of and satisfy

When Neimark-Sacker bifurcation occurs, the modulus of the eigenvalues (25) should be varying with the bifurcation parameter. This requires

A direct calculation on leads to

Since , then we find that the value of is indeed not equal to zero. Another condition for the occurrence of Neimark-Sacker bifurcation is that the eigenvalues (25) are exactly neither real numbers nor imaginary numbers, implying

which is equivalent to . From (21), we find that . It means , i.e.,

On the basis of map (23), the canonical form is studied to obtain the last condition for the occurrence of Neimark-Sacker bifurcation. Applying the following transformation

on map (23), then we obtain

where

The second-order and third-order partial derivatives of and at (0, 0) are calculated as

According to the bifurcation theorem [33], if the map (32) undergoes Neimark-Sacker bifurcation at around , the following discriminatory quantity must not be zero, i.e.,

in which

A calculation on the quantity gives

where

Based on the above calculations, the occurrence of Neimark-Sacker bifurcation in the discrete predator-prey system needs the satisfaction of conditions (21), (30), and (38). Moreover, if and , then an attracting invariant closed curve bifurcates from the fixed point for ; otherwise, if and , a repelling invariant closed curve bifurcates from the fixed point for .

3.2. Turing Instability Analysis

When spatially heterogeneous perturbations occur at the homogeneous states, the spatiotemporal discrete predator-prey system may converge to heterogeneous states [10]. To determine the occurrence conditions for the stable heterogeneous states, Turing instability analysis is performed. Firstly, spatially heterogeneous perturbations are introduced to perturb the stable homogeneous stationary state , given by

It should be noticed that

Based on the calculations in former approach [10], the eigenvalues of the operator are determined as

in which

Substituting (40) and (41) into the spatiotemporal discrete predator-prey system leads to

where . The two-order terms in (44) can be ignored when the perturbations are small. Using the corresponding eigenfunction of the eigenvalues to multiply (44) gets

Summing (45) for all of and obtains

Let

Equations (47) can be transformed into the following form [10, 26]:

The Turing instability emerges when system (48) has unstable eigenvalues at the fixed point (0, 0) [10, 26]. Calculating these two eigenvalues of system (48), then we obtain

where

According to previous research works [10, 26–28], the occurrence of Turing instability demands that existing one group of to make the maximum value of larger than one. Define

Thus, when , Turing instability occurs; when , the spatiotemporal discrete predator-prey system returns to the homogeneous states under spatially heterogeneous perturbations. Moreover, the threshold condition for the occurrence of Turing bifurcation requires . Consequently, there exists a critical value satisfying .

4. Numerical Results

The purpose of numerical simulations is to demonstrate the Neimark-Sacker bifurcation as well as pattern formation under the influence of Neimark-Sacker-Turing instability. For such demonstration, the parametric conditions should be provided firstly. Based on former research work of Rao and Kang [3], a group of ecologically feasible parameter values can be given as = 1.63, = 0.1, = 0.32, = 0.22, = 0.6, = 1.2, and = 200. The values of parameters and should be set to ensure nonnegativity of the state variables and . In the following simulations, we let = 10 and as a varying parameter. Then according to previous research works [10, 27], the nonnegativity of and needs that and are both less than 0.5.

Under the above provided parametric conditions, the fixed point of the discrete predator-prey system is determined as (0.6538, 02972), and the critical point for Neimark-Sacker bifurcation is calculated as = 4.9444. When , the fixed point of the system is stable. Taking = 4.5 as an example, the two eigenvalues of the fixed point are = 0.8307 + 0.5287i, = 0.8307 − 0.5287i, the modulus of which is less than one, suggesting is a stable focus. When , we get the two eigenvalues of as = 0.8140 + 0.5809i, = 0.8140 − 0.5809i, and = 1. Then we check conditions (28), (30), and (38) for Neimark-Sacker bifurcation; we have = 0.0152 ≠ 0, = −1.6280 ≠ 0,1 and = −179.8092 ≠ 0. Hence, according to the Neimark-Sacker bifurcation theorem, the discrete predator-prey system indeed undergoes Neimark-Sacker bifurcation at . Furthermore, since < 0 and > 0, we know that an attracting invariant closed curve bifurcates from the fixed point for , and the stable point loses its stability.

Figure 1 shows the Neimark-Sacker bifurcation diagram of the discrete predator-prey system. From the bifurcation diagram and corresponding maximum Lyapunov exponent, it is known that the Neimark-Sacker bifurcation induces a route to chaos, leading to a dynamic transition from fixed point, to invariant curves, with periodic windows occurring in-between, to chaotic dynamics. Invariant curves suggest that the discrete predator-prey system follows dynamical behaviors which are homogeneous in space and quasiperiodically oscillating in time. In the periodic windows, we find period 11, 20, 31, 41, 51, 62, 72, and 82 orbits, revealing spatially homogeneous and periodically oscillating behaviors of the system. The maximum Lyapunov exponent diagram in Figure 1(b) says that the predator-prey system begins entering spatially homogeneous chaotic dynamics with the increase of the value of and finally enters chaotic dynamics at around = 5.778. Figure 2 explicitly displays the three typical dynamical behaviors in the phase portraits, i.e., invariant closed curve, periodic orbit (period 31), and chaotic attractor.

(a)

(b)

Figure 1: (a) Neimark-Sacker bifurcation diagram and (b) maximum Lyapunov exponent corresponding to (a). a = 1.63, b = 0.1, c = 0.32, and d = 0.22.

(a)τ = 5.1

(b)τ = 5.525

(c)τ = 5.83

Figure 2: Three typical dynamical behaviors on the route to chaos induced by the Neimark-Sacker bifurcation, (a) invariant closed curve, (b) period-31 orbit, and (c) chaotic attractor. = 1.63, = 0.1, = 0.32, and = 0.22.

When heterogeneous perturbations occur on the spatially homogeneous states, the discrete predator-prey system may experience Turing instability and the system dynamics converges to stable spatially heterogeneous states. Under the previous given parametric conditions, the critical point for Turing bifurcation is overlapping with the Neimark-Sacker bifurcation critical point. Via calculation, we find the values corresponding to the three cases in Figure 2 are (5.1) = 1.0060, (5.525) = 1.0241, and (5.83) = 1.0385, respectively. Therefore, it is verified that Turing instability can indeed take place on the route to chaos induced by Neimark-Sacker bifurcation. For the three cases in Figure 2, Turing instability is spatial symmetry breaking on spatially homogeneous periodically, quasiperiodically, and chaotically oscillating states, resulting in different pattern formation. For simulating the pattern formation, two types of initial conditions are applied, i.e., stochastic initial condition and deterministic initial condition. Stochastic initial condition is given by perturbing the homogeneous stationary state (, ) with small random perturbations. And deterministic initial condition is given by the following:

Figures 3–5 show the evolution of spatiotemporal patterns of prey population (the predator patterns exhibit similar configuration) corresponding to spatial symmetry breaking occurring on the three dynamical behaviors in Figure 2, starting from the stochastic initial condition. At the beginning, the prey patterns evolve to a spatial configuration which is composed of large patches of higher and lower prey densities (such as Figures 3(b), 4(c), and 5(c)). As the time progresses, the large prey patches gradually fragment and change to smaller prey patches, the shape of which may present as stripes, arcs, spots, rings, etc. These smaller patches display disordered distribution in space and may combine together to form complex structure. After long-term evolution of patterns, we find that the prey patterns are always oscillating similarly in time and never reach equilibrium configuration.

(a) = 1

(b) = 500

(c) = 1000

(d) = 1500

(e) = 2000

(f) = 5000

Figure 3: Pattern formation resulting from spatially symmetry breaking on invariant closed curves and starting with stochastic initial condition. = 1.63, = 0.1, = 0.32, = 0.22, = 0.6, = 1.2, = 5.1, = 10, and = 200.

(a) = 1

(b) = 10

(c) = 100

(d) = 500

(e) = 1000

(f) = 5000

Figure 4: Pattern formation resulting from spatially symmetry breaking on period-31 orbit and starting with stochastic initial condition. = 1.63, = 0.1, = 0.32, = 0.22, = 0.6, = 1.2, = 5.525, = 10, and = 200.

(a) = 1

(b) = 10

(c) = 100

(d) = 500

(e) = 2000

(f) = 5000

Figure 5: Pattern formation resulting from spatially symmetry breaking on chaotic attractor and starting with stochastic initial condition. = 1.63, = 0.1, = 0.32, = 0.22, = 0.6, = 1.2, = 5.83, = 10, and = 200.

As exhibited in Figures 3–5, pattern formation in the three different cases presents similar evolution trend and process. However, two main differences are found in these pattern evolutions. When the local kinetics of the discrete predator-prey system goes further on the route to chaos, first, the pattern evolution speed becomes faster, and second, the prey patterns evolving to long time shows more fragmented, more disordered characteristics with smaller patches.

In order to exhibit the dynamical trajectory of pattern evolution, Figure 6 is plotted using the change of state variables at one lattice in Figures 3–5, which here is chosen as lattice (100, 100), with the increase of time. As shown in Figures 6(a)–6(c), three new attractors emerge in the phase portraits. Comparing with Figure 2, it can be found that, under the influence of Neimark-Sacker-Turing instability, the invariant cycle and periodic orbit break and transform to more complex attractors. The clouds of random points in Figures 6(b)-6(c) suggest that and may experience almost all states in the black area. Moreover, the evolution process of prey patterns in the three cases as time progresses presents disordered waves, as displayed by Figures 6(d)–6(f).

(a)

(b)

(c)

(d) = 5.1

(e) = 5.525

(f) = 5.83

Figure 6: Phase portraits and wave diagrams corresponding to the change of state variables in lattice (100, 100) of Figures 3–5 with the increase of time.

To explore the nonlinear characteristics of pattern evolution process in Figures 3–5, sensitivity analysis is carried out with the variations of initial conditions and parameter and the results are shown in Figure 7. In Figures 7(a)–7(c), we apply two close initial conditions, and , the difference of which is given as in merely one lattice, and then we observe the change of in another lattice. At the beginning, stays at zero; after a certain time of pattern evolution, the amplitude of suddenly bursts and follows chaotic oscillations. For the pattern formation in the three cases of Figures 3–5, this phenomenon suggests the occurrence of spatiotemporal chaos, which can result in the self-organization of numerous different spatiotemporal patterns with the variation of initial conditions. Simultaneously, it also reveals that the discrete predator-prey system holds the property of multistability when Neimark-Sacker instability and Turing instability both occur. Figures 7(d)–7(f) show the sensitivity of pattern evolution to the variation of parameter , which is controlled as 0.00001. Even with such tiny variation of parameter , the pattern evolution process in Figures 3–5 will follow totally different trajectories.

(a)

(b)

(c)

(d) = 5.1

(e) = 5.525

(f) = 5.83

Figure 7: Sensitive analysis of pattern evolution process to variations of initial conditions (graphs a–c) and parameter (graphs d–f), with the application of state variables changing process in Figure 6.

The above sensitivity analysis suggests that change of initial conditions may lead to the self-organization of different patterns. Figures 8–10 demonstrate the evolution of prey patterns with the application of deterministic initial condition as described in (52). In these three cases, the pattern evolution starts from a square in the background. Then, a target wave pattern emerges replacing the square. As the rings of target wave pattern gradually expand, a new different structure may emerge at the core of the target and the target wave breaks from the corners, leading to the symmetry breaking increasing with time. In the pattern evolution process, symmetrical irregular transient patterns dominate the pattern evolution process all along, such as in Figure 8(g), from which one can see that the degree of symmetry is still strong as vertical, horizontal, and central symmetry. However, at large iteration time (e.g., = 200000), the patterns for the three cases all become disordered and irregular and the spatial symmetry of the patterns reduces to be weak, like the horizontal symmetry observed in Figures 8(i), 9(f), and 10(f).

(a) = 1

(b) = 500

(c) = 2000

(d) = 3000

(e) = 5000

(f) = 7000

(g) = 10000

(h) = 100000

(i) = 200000

Figure 8: Pattern formation resulting from spatially symmetry breaking on invariant closed curves and starting with deterministic initial condition. = 1.63, = 0.1, = 0.32, = 0.22, = 0.6, = 1.2, = 5.1, = 10, and = 200.

(a) = 1

(b) = 500

(c) = 1000

(d) = 1500

(e) = 2000

(f) = 200000

Figure 9: Pattern formation resulting from spatially symmetry breaking on period-31 orbit and starting with deterministic initial condition. = 1.63, = 0.1, = 0.32, = 0.22, = 0.6, = 1.2, = 5.525, = 10, and = 200.

(a) = 1

(b) = 100

(c) = 500

(d) = 1000

(e) = 1500

(f) = 200000

Figure 10: Pattern formation resulting from spatially symmetry breaking on chaotic attractor and starting with deterministic initial condition. = 1.63, = 0.1, = 0.32, = 0.22, = 0.6, = 1.2, = 5.83, = 10, and = 200.

The decrease of symmetry in Figures 8–10 as the time progresses results from the combined effects of Turing instability and Neimark-Sacker instability. A space-independent Neimark-Sacker bifurcation breaks the temporal symmetry of the system and gives rise to oscillations that are uniform in space and periodic in time. On the contrary, the diffusion-driven Turing instability breaks the spatial symmetry and manifests itself through the formation of patterns that are stationary in time and oscillatory in space. As explored in literature, the combined effects of both may lead to the emergence of spatiotemporal chaos, which presents disordered characteristics [10, 27]. The deterministic initial conditions provide order states which the predator-prey system starts and keeps at the beginning of the evolution process, but with the increase of time the spatiotemporal chaos will dominate and result in disordered states with lessening spatial symmetry.

Figure 11 shows the dynamical trajectories for the three cases of pattern evolution in Figures 8–10. Compared with Figure 6, we find that although the initial conditions are different, the attractors in phase portraits present similar configurations. Further sensitivity analysis demonstrates that the pattern evolution processes in these cases are also sensitive to tiny variation of parameter .

(a) = 5.1

(b) = 5.525

(c) = 5.83

Figure 11: Phase portraits corresponding to the change of state variables in lattice (100, 100) of Figures 8–10 with the increase of time.

The simulation results in Figures 8–11 demonstrate that on the front part of the route to chaos (such as < < 5.2), patterns can still exhibit strongly symmetrical property at large evolution time (see = 10000 in Figure 8), including horizontal, vertical, and central symmetry. Figure 12 displays a few transient strongly symmetrical patterns with the variation of parameter . These patterns show rich regular spatial structures. It reveals the spatiotemporal complexity of the discrete predator-prey system and the regularity in pattern self-organization.

(a) = 5.0

(b) = 5.02

(c) = 5.04

(d) = 5.06

(e) = 5.08

(f) = 5.12

(g) = 5.14

(h) = 5.16

(i) = 5.18

Figure 12: Various transient strongly symmetrical patterns when Neimark-Sacker-Turing instability occurs on the front part of route to chaos. = 1.63, = 0.1, = 0.32, = 0.22, = 0.6, = 1.2, = 10, = 200, and = 10000.

5. Discussion and Conclusions

This research investigates the spatiotemporal dynamics of a space- and time-discrete predator-prey system with Allee effect with the application of CML model. The development of the CML model in the present research is based on two types of discrete models in literature. The first is the type derived by directly discretizing partial differential equations with finite difference scheme, for example, the exponential discrete Lotka-Volterra model proposed by Li et al. [34]. The second type is the CML model as given by Mistro et al. [24], which abandoned the structure of partial differential equations. Taking the advantages of the two types of models, the present CML model is still developed from discretization of partial differential equations and simultaneously utilizes the structure of the CML model in Mistro et al. [24]. Therefore, the same ecological mechanisms with the corresponding continuous model as well as the advantages of CML model are both kept. The present CML model further bridges the continuous models and the discrete models [10].

In the research, we find new nonlinear characteristics for the system dynamics and pattern formation under Neimark-Sacker-Turing instability. Summarizing the results obtained in theoretical analysis and numerical simulations, the concluding remarks of this research can be made as the following.

(1)The system can undergo Neimark-Sacker bifurcation inducing a route to chaos, on which three typical dynamical behaviors emerge, including invariant closed curves, periodic orbits, and chaotic attractors.(2)Neimark-Sacker-Turing instability triggers spatial symmetry breaking on invariant closed curves, periodic orbits, and chaotic attractors and leads to the formation of heterogeneous patterns, which keep oscillating property all along.(3)Different initial conditions result in diverse pattern evolution processes, revealing the multistability of the predator-prey system.(4)Under the effect of Turing instability, the invariant closed curves and periodic orbits change to more complex chaotic attractors. Moreover, spatiotemporal chaos occurs when Neimark-Sacker instability and Turing instability both take place.(5)With the application of deterministic initial conditions, transient symmetrical patterns are observed and the symmetry breaking increases as time continuously progresses. Moreover, when local kinetics of the system goes further on the route to chaos, the speed of symmetry breaking becomes faster, leading to more fragmented and more disordered patterns at the same evolution time.(6)Under Neimark-Sacker-Turing instability, the pattern evolution process is disordered and chaotic. However, at transient times, the predator-prey pattern may exhibit some symmetrical property. This demonstrates that ordered structures can exist in chaotic processes.(7)Since pattern evolution processes are sensitive to variations of initial conditions and parameter , the system generates rich spatiotemporal complexity under the coupled effects of Neimark-Sacker instability and Turing instability. This reveals numerous possible ways for predator-prey interactions and coexistence.

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that they have no financial and nonfinancial conflicts of interest for this research work. They confirm that the received funding mentioned in Acknowledgments will not lead to any conflicts of interest regarding the publication of this manuscript.

Acknowledgments

The authors would like to acknowledge with great gratitude the support of the National Major Science and Technology Program for Water Pollution Control and Treatment (nos. 2017ZX07101-002 and 2015ZX07203-011) and the Fundamental Research Funds for the Central Universities (no. JB2017069).

References

A. J. Lotka, The Elements of Physical Biology, Williams & Wilkins, Baltimore, Maryland, USA, 1925.

V. Volterra, “Fluctuations in the abundance of a species considered mathematically,” Nature, vol. 118, no. 2972, pp. 558–560, 1926.

F. Rao and Y. Kang, “The complex dynamics of a diffusive prey–predator model with an Allee effect in prey,” Ecological Complexity, vol. 28, pp. 123–144, 2016.

M. R. Garvie, “Finite-difference schemes for reaction-diffusion equations modeling predator-prey interactions in MATLAB,” Bulletin of Mathematical Biology, vol. 69, no. 3, pp. 931–956, 2007.

E. Hernandezmartinez, H. Puebla, T. Perezmunoz et al., “Spatiotemporal dynamics of telegraph reaction-diffusion predator-prey models,” International Symposium on Mathematical and Computational Biology, pp. 268–281, 2012.

M. Banerjee and S. Petrovskii, “Self-organised spatial patterns and chaos in a ratio-dependent predator-prey system,” Theoretical Ecology, vol. 4, no. 1, pp. 37–53, 2011.

G.-Q. Sun, J. Zhang, L.-P. Song, Z. Jin, and B.-L. Li, “Pattern formation of a spatial predator-prey system,” Applied Mathematics and Computation, vol. 218, no. 22, pp. 11151–11162, 2012.

L. Zhang, J. Liu, and M. Banerjee, “Hopf and steady state bifurcation analysis in a ratio-dependent predator-prey model,” Communications in Nonlinear Science and Numerical Simulation, vol. 44, pp. 52–73, 2017.

S. A. Levin and L. A. Segel, “Hypothesis for origin of planktonic patchiness,” Nature, vol. 259, no. 5545, p. 659, 1976.

T. Huang, H. Zhang, H. Yang, N. Wang, and F. Zhang, “Complex patterns in a space- and time-discrete predator-prey model with Beddington-DeAngelis functional response,” Communications in Nonlinear Science and Numerical Simulation, vol. 43, pp. 182–199, 2017.

R. K. Upadhyay, N. Kumari, and V. Rai, “Wave of chaos and pattern formation in spatial predator-prey systems with Holling type IV predator response,” Mathematical Modelling of Natural Phenomena, vol. 3, no. 4, pp. 71–95, 2008.

J. Wang, J. Shi, and J. Wei, “Dynamics and pattern formation in a diffusive predator-prey system with strong Allee effect in prey,” Journal of Differential Equations, vol. 251, no. 4-5, pp. 1276–1304, 2011.

S. V. Petrovskii, A. Y. Morozov, and E. Venturino, “Allee effect makes possible patchy invasion in a predator-prey system,” Ecology Letters, vol. 5, no. 3, pp. 345–352, 2002.

G.-Q. Sun, L. Li, Z. Jin, Z.-K. Zhang, and T. Zhou, “Pattern dynamics in a spatial predator-prey system with allee effect,” Abstract and Applied Analysis, vol. 2, Article ID 921879, pp. 18–29, 2013.

L. A. D. Rodrigues, D. C. Mistro, and S. Petrovskii, “Pattern formation in a space- and time-discrete predator-prey system with a strong Allee effect,” Theoretical Ecology, vol. 5, no. 3, pp. 341–362, 2012.

S. Ghorai and S. Poria, “Turing patterns induced by cross-diffusion in a predator-prey system in presence of habitat complexity,” Chaos Solitons & Fractals the Interdisciplinary Journal of Nonlinear Science & Nonequilibrium & Complex Phenomena, vol. 91, pp. 421–429, 2016.

Z. Ling, L. Zhang, and Z. Lin, “Turing pattern formation in a predator-prey system with cross diffusion,” Applied Mathematical Modelling: Simulation and Computation for Engineering and Environmental Systems, vol. 38, no. 21-22, pp. 5022–5032, 2014.

M. Dam, M. Brøns, J. Juul Rasmussen, V. Naulin, and J. S. Hesthaven, “Sparse identification of a predator-prey system from simulation data of a convection model,” Physics of Plasmas, vol. 24, no. 2, Article ID 293001, 2017.

A. M. Turing, “The chemical basis of morphogenesis,” Philosophical Transactions of the Royal Society of London. Series B. Biological Sciences, vol. 237, no. 641, pp. 37–72, 1952.

C. Tian, Z. Ling, and L. Zhang, “Nonlocal interaction driven pattern formation in a prey-predator model,” Applied Mathematics and Computation, vol. 308, pp. 73–83, 2017.

K. Kaneko, “Period-doubling of kink-antikink patterns, quasiperiodicity in antiferro-like structures and spatial intermittency in coupled logistic lattice. Towards a prelude of a “field theory of chaos”,” Progress of Theoretical and Experimental Physics, vol. 72, no. 3, pp. 480–486, 1984.

K. Kaneko, “Transition from torus to chaos accompanied by frequency lockings with symmetry breaking. In connection with the coupled-logistic map,” Progress of Theoretical and Experimental Physics, vol. 69, no. 5, pp. 1427–1442, 1983.

K. Kaneko, “Spatiotemporal intermittency in coupled map lattices,” Progress of Theoretical and Experimental Physics, vol. 74, no. 5, pp. 1033–1044, 1985.

D. C. Mistro, L. A. D. Rodrigues, and S. Petrovskii, “Spatiotemporal complexity of biological invasion in a space- and time-discrete predator-prey system with the strong Allee effect,” Ecological Complexity, vol. 9, pp. 16–32, 2012.

D. Punithan, D.-K. Kim, and R. I. McKay, “Spatio-temporal dynamics and quantification of daisyworld in two-dimensional coupled map lattices,” Ecological Complexity, vol. 12, pp. 43–57, 2012.

T. Huang and H. Zhang, “Bifurcation, chaos and pattern formation in a space- and time-discrete predator-prey system,” Chaos Solitons & Fractals the Interdisciplinary Journal of Nonlinear Science & Nonequilibrium & Complex Phenomena, vol. 91, pp. 92–107, 2016.

T. Huang, H. Zhang, and H. Yang, “Spatiotemporal complexity of a discrete space-time predator-prey system with self- and cross-diffusion,” Applied Mathematical Modelling: Simulation and Computation for Engineering and Environmental Systems, vol. 47, pp. 637–655, 2017.

F. Zhang, H. Zhang, S. Ma, T. Meng, T. Huang, and H. Yang, “Self-organized patterns induced by Neimark-Sacker, flip and turing bifurcations in a discrete predator-prey model with Lesie-Gower functional response,” Entropy, vol. 19, no. 6, p. 258, 2017.

L. Cheng and H. Cao, “Bifurcation analysis of a discrete-time ratio-dependent predator-prey model with Allee effect,” Communications in Nonlinear Science and Numerical Simulation, vol. 38, pp. 288–302, 2016.

W. C. Allee, Animal Aggregations: A Study in General Sociology, University of Chicago Press, Chicago, Ill, USA, 1931.

A. Morozov, S. Petrovskii, and B.-L. Li, “Spatiotemporal complexity of patchy invasion in a predator-prey system with the Allee effect,” Proceedings Biological Sciences, vol. 271, pp. 1407–1414, 2006.

Y. Li, “Oscillation deduced from the Neimark-Sacker bifurcation in a discrete dual Internet congestion control algorithm,” Advances in Difference Equations, vol. 1, pp. 1–11, 2013.

J. Guckenheimer and P. Holmes, “Nonlinear oscillations, dynamical systems and bifurcations of vector fields,” Physics Today, vol. 38, no. 11, pp. 102–105, 1985.

M. Li, B. Han, L. Xu, and G. Zhang, “Spiral patterns near Turing instability in a discrete reaction diffusion system,” Chaos, Solitons & Fractals, vol. 49, pp. 1–6, 2013.