Pattern Formation in a Semi-Ratio-Dependent Predator-Prey System with Diffusion

We investigate spatiotemporal dynamics of a semi-ratio-dependent predator-prey system with reaction-diffusion and zero-flux boundary. We obtain the conditions for Hopf, Turing, and wave bifurcations of the system in a spatial domain by making use of the linear stability analysis and the bifurcation analysis. In addition, for an initial condition which is a small amplitude random perturbation around the steady state, we classify spatial pattern formations of the systemby using numerical simulations.The results of numerical simulations unveil that there are various spatiotemporal patterns including typical Turing patterns such as spotted, spot-stripelike mixtures and stripelike patterns thanks to the Turing instability, that an oscillatory wave pattern can be emerged due to the Hopf and wave instability, and that cooperations of Turing and Hopf instabilities can cause occurrence of spiral patterns instead of typical Turing patterns. Finally, we discuss spatiotemporal dynamics of the system for several different asymmetric initial conditions via numerical simulations.


Introduction
In recent years, pattern formations in nonlinear complex systems have been one of the central problems of the natural, social, technological sciences and ecological systems [1][2][3][4][5][6][7].Particularly, many researchers have studied the prey-predator system with reaction-diffusion.For example, in [8], Garvie and Trenchea presented the analysis of reaction-diffusion predator-prey systems with the Holling type II functional response and provided an  ∞ a priori estimate.Wang et al. in [9] investigated the spatial pattern formation of a predatorprey system with prey-dependent functional response of Ivlev type and reaction-diffusion.Also Zhang et al. in [10] studied a linear stability and bifurcation analysis including the Hopf and Turing bifurcation for a spatial Holling-type IV predator-prey model.In addition, Camara and Aziz-Alaoui, the authors of [11], considered a predator-prey system with a modified Leslie-Gower functional response modeled by a reaction-diffusion equation and derived the conditions for Hopf and Turing bifurcation in the spatial domain.
In this context, in this paper, we will focus on the following a semi-ratio-dependent predator-prey system with reaction-diffusion: where ,  stand for the population densities of prey and predator with time , respectively.The parameters , , , and  are all positive constants standing for the prey intrinsic growth rate, the carrying capacity of prey, the capturing rate, and the half-saturation constant, respectively.The predator grows logically with the growth rate  and the carrying capacity /ℎ.In fact, the parameter ℎ is a measure of the food quality that the prey provides for conversion into the predator.And  1 and  2 are diffusion coefficients.Here ∇ 2 =  2 / 2 + 2 / 2 is the usual Laplacian in two-dimensional and  and  stand for the space.
In fact, the authors in [12][13][14] have explored a similar model to system (1) without reaction-diffusion, which is a kind of a temporal predator-prey system, and figured out dynamic properties of the temporal system such as stabilities of equilibria and permanence of the system.However, as mentioned above, it is important to investigate spatiotemporal dynamical behaviors of a diffusive predator-prey system like system (1).Thus, in this paper, we will concentrate our concerns on studying the bifurcation analysis and the spatiotemporal pattern formation analysis of system (1).
Throughout this paper, we assume that no external input is imposed from the outside.Hence the boundary conditions are taken as / = / = 0 on Ω, Ω is the spatial domain in R 2 , and  is the outward unit normal vector of the boundary Ω.The main reason for choosing such boundary conditions is that we are interested in the self-organization of pattern.Also system (3) needs to be analyzed with the initial populations (, , 0) > 0 and (, , 0) > 0.
In order to minimize the number of parameters in system (1), we set For simplicity we will omit the bar notation in , ,  in the rest of the paper.Thus we can have the following system containing dimensionless quantities: where / = , / = , ℎ/ = ,  1 / =  1 , and  2 / =  2 .
The main object of this paper is to look into the spatial dynamic behaviors of system (3).For this, in Section 2, we investigate bifurcation phenomena, Hopf, Turing, and wave bifurcations.In particular, we find out sufficient conditions for the Hopf bifurcation, the Turing instability, and the wave bifurcation of system (3).In Section 3, we give numerical simulations of system (3) with random perturbation initial conditions of the stationary solution of the spatially homogeneous system to investigate various spatiotemporal pattern formations such as spotted, stripelike, and spiral patterns.Finally, we discuss dynamical behaviors of system (3) with several different initial conditions in Section 4.

Hopf Bifurcation.
Hopf bifurcation is an instability induced by the transformation of the stability of a focus.In fact, the space-independent Hopf bifurcation breaks temporal symmetry of a system and gives rise to oscillations that are uniform in space and periodic in time.Mathematically speaking, the Hopf bifurcation occurs when Thus we can figure out that the critical value of Hopf bifurcation parameter  is equal to We assume that the value   is greater than 0 for the existence of Hopf bifurcation.In other words, from now on, we will consider the parameters  and  which are satisfied with the following condition: Thus the positive stationary solution ( * ,  * ) of system ( 3) is an unstable equilibrium solution if  <   and oscillates periodically with respect to time with the frequency The corresponding wavelength is 2.2.Turing Instability.Turing instability (or called Turing bifurcation) is a phenomenon that causes certain reactiondiffusion system to lead to spontaneous stationary configuration.That is why Turing instability is often called diffusiondriven instability.The Turing instability is not dependent on the geometry of the system but only on the reaction rates and diffusion.It can occur only when the inhibitor () diffuses faster than the activator () [3,8,9,20,21,28].
In fact, the Turing instability sets in when at least one of the solutions of (8) crosses the imaginary axis.In other words, the spatially homogeneous steady state will become unstable due to heterogeneous perturbation when at least one solution of ( 8) is positive.For this reason, at least one out of the following two inequalities is violated to occur the Turing instability phenomenon: Since perturbations of zero wavenumber are stable when diffusive instability sets in by definition (due to the stability for the nonspatial steady state), tr() =   +   < 0 and det() =     −     > 0 are satisfied.It is seen from these facts that the first condition in ( 16) always holds.Hence we  have only one for the instability condition; that is, ( 2 ) ≡ det(  ) < 0. Elementary calculations yield that Thus the minimum of ( 2 ) occurs at the critical wavenumber  2  , where By substituting  2 =  2  into ( 2 ), we can get a sufficient condition for the Turing instability as follows: where It is from replacing the inequality in (19) by the equality that the critical value of bifurcation parameter  for the Turing bifurcation equals At the Turing threshold   , the spatial symmetry of the system is broken and the patterns are stationary in time and oscillatory in space with the wavelength

Wave Bifurcation.
The wave instability caused by the wave bifurcation plays an important part in pattern formations in many areas [7].In fact, mathematically speaking, the wave bifurcation occurs when Thus, from elementary calculation, the critical value of wave bifurcation parameter  can be obtained as where It is well known that, at the wave threshold   , both spatial and temporal symmetries are broken and the patterns are oscillatory in space and time with the wavelength where If one takes  = 0.4,  = 5/6, and  1 = 1 in system (3), the bifurcation diagram, shown in Figure 1, can be obtained using the above linear stability analysis.Since Turing patterns occur when the inhibitor () diffuses faster than the activator (), we take into account the range of the value  2 greater than  1 = 1.According to Figure 1, the Hopf bifurcation line and the Turing bifurcation line intersect at a point (4.0611, 0.3587) and the Turing bifurcation line and the wave bifurcation line meet at a point (1, 0.0932).In addition, these bifurcation lines split the parametric space into five distinct domains.In domain I, located above two bifurcation lines, the steady state is the only stable solution of system (3).Domain II is the region of the pure Turing instability and domain III is the region of Hopf and wave instabilities.Hopf   and Turing instabilities take place in domain IV.On the other hand, all three instabilities occur in domain V. Furthermore, Figure 2 shows the dispersion relations corresponding to several values of two parameters  2 , .The transitions of Hopf, Turing, and wave modes from stable to unstable for system (3) can be ascertained clearly via Figure 2.While Hopf is subcritical, Turing and wave bifurcations are supercritical.

Spatiotemporal Pattern Analysis via Numerical Simulations
In this section, we will investigate spatiotemporal pattern formations of the spatially extended system (3) in twodimensional space via numerical examples.In fact, the nonuniform stationary states of system (3) that corresponds to spatial patterns cannot be found analytically.In other words, the analytical methods are not sufficient to fully understand the system, which is a reason why we need to use computer simulations.For this, first, take  = 0.4,  = 5/6, and  1 = 1 in system (3).Then the spatially homogeneous equilibrium ( * ,  * ) = (0.1892, 0.2270) is always the only nontrivial positive stationary state for the nonspatial system of system (3) regardless of the values  2 and .
In order to solve partial differential equations numerically, one has to discretize the space and time of the given problem.For this reason, the discrete domain for the continuous domain Ω in system ( 3) is assumed to be the 200 × 200 lattice sites and the spacing between the lattice points is assumed to be Δ = Δ = 1.The time evaluation is also discrete; that is, the time goes in steps of Δ.And all our numerical simulations employ the zero-flux boundary condition.
System ( 3) is solved numerically by using a finite difference scheme for the spatial derivatives and an explicit Euler method for the time integration with the time step Δ = 0.01.These numerical schemes for the diffusion equation ( 3) give stable solutions as long as the equalities  1 Δ/(Δ) 2 ,  2 Δ/(Δ) 2 < 1/2 are satisfied (see [31]).
It is well known that the spatiotemporal dynamics of a diffusion-reaction system depends on the choice of initial conditions [3,9,32,33].In this section, an initial condition is taken as a small amplitude random perturbation around the steady state ( * ,  * ) since it is very natural from the biological point of view.We stop the simulation when the numerical solutions either reach a stationary state or show oscillatory behaviors.
In the numerical simulations, different types of spatiotemporal dynamics are observed and we have found that the distributions of predator and prey are always of the same type.Consequently, we can restrict our analysis of pattern formations to one distribution.In this section, we show the distribution of predator, for instance.Now, we will classify spatiotemporal patter formations of system (3) for each bifurcation domain in Figure 1 by taking advantage of numerical simulations.
As mentioned in Section 2, the domain I in Figure 1 is the region where the steady state is the only stable solution of system (3).Thus we will omit to mention dynamical behaviors of system First, we will investigate dynamical behaviors of system (3) with a random small perturbation initial condition of the stationary solution ( * ,  * ) by taking into account the parameters in domain II.As mentioned in Section 2, the domain II is the region of the pure Turing instability.If we take the parameter values  2 = 7 and  = 0.6, we can see from Figure 3 that the random initial distribution (see Figure 3(a)) leads to the formation of a regular macroscopic spotted pattern which prevails over the whole domain at last, and after that the dynamics of the system does not change their behaviors anymore.
From now on, through this section, we will not display the snapshots of the initial pattern as Figure 3(a) since all following figures are obtained with a random small perturbation initial condition of the stationary solution as the previous case.
If we take the parameter values  2 = 7 and  = 0.535, we can observe from Figure 4 that the steady state of spotted pattern and the stripelike pattern coexist.In addition, for the parameters  2 = 7 and  = 0.4, Figure 5 shows that stripelike spatial patterns are prevalent in the whole domain eventually.In fact, we can see that the stripelike patterns grow steadily with time and finally they overwhelm the space domain (see Figure 5(c)) after the spotted pattern and the stripelike pattern coexist (see Figures 5(a Next, in order to study dynamical behaviors of system (3) in domain V, set the parameters  2 = 7 and  = 0.34.We can know from Figure 1 that, in this case, all three instabilities occur at the same time.In this context, we display Figure 6 to show snapshots of predator spatial patterns at 50000, 100000, and 200000 iterations.Since Turing instabilities take place in domain V, we can see the stripelike spatial patterns as shown in Figure 5.However, an oscillatory pattern in time can be observed even though we cannot observe this phenomenon via snapshots of Figure 6.To make it clear, we exhibit the local phase portrait of system (3) for a fixed point (0.1892, 0.2270) in Figure 7.As one sees, thanks to the Hopf instability, a stable cycle can be observed.Moreover, we can ascertain numerically that the frequency of the periodic oscillations in time  = 0.2685 and the corresponding wavelength  = 23.4.Also, it follows from ( 14) and ( 15) that the frequency Thus Figures 3,4,5,and 6 show typical Turing patterns of predator in system (3) for the fixed value  2 = 7 when the parameter  varies.In fact, as the value  decreases, the sequence spotted → spot-stripelike mixtures → stripelike pattern is observed.Now, consider the parameters  2 = 2 and  = 0.34 in domain III.It is from the linear stability analysis in Section 2 that Hopf and wave instabilities occur.As shown in Figure 8, a special pattern may not be formed such as spotted or stripelike patterns since Turing instabilities do not occur in this domain.However, in this case, an oscillatory wave pattern formation is observed.In order to deal with this phenomenon in more detail, we display the local phase portrait of system (3) for a fixed point (0.1892, 0.2270) in Figure 9(a), which shows a stable limit cycle, and the population density of predator in Figure 9(b), which is periodic in time and the space-time plots in Figure 9(c), obtained by piling up the lines  = 100 in each snapshot as time progresses.Thus, we can clearly figure out that Hopf and wave instabilities play an important role in forming a wave and periodic oscillatory pattern in space and time as mentioned in Section 2. In fact, numerical calculations yield that the frequency of the periodic oscillations in time  = 0.2474 and the corresponding wavelength  = 25.4.However, the theoretical values of the frequency and the corresponding wavelength are   = 0.2536 and   = 24.7738,respectively.However, the theoretical values   = 0.0789 and   = 79.6768can be obtained from ( 26) and (27).From these facts, we can surmise that, in this case, Hopf instabilities have an stronger effect on pattern formation in terms of time than wave instabilities.
On the contrary, domain IV is the region of Hopf and Turing instabilities according to Figure 1.To discuss dynamical behaviors of system (3) in this domain, we think about two parameter sets as  2 = 2,  = 0.15 and  2 = 7,  = 0.2.As one sees, Figures 10 and 11 show one of the typical spatial patterns called a spiral pattern instead of patterns shown in Turing domains II and V even if Turing instabilities occur in this domain.In fact, in many recent researches [8,9,28,33], spiral patterns have been reported.However, most of them used specific initial conditions which are not the same as the conditions used in this section to generate a spiral pattern.Thus a diffusive system starting with a small amplitude random perturbation around the steady state could have spiral patterns thanks to the cooperations of Turing and Hopf instabilities.

Conclusion and Discussion
Various spatiotemporal patterns generated by a semi-ratiodependent predator-prey system with reaction-diffusion are studied theoretically and numerically.The spatial parametric domain is divided into five regions by three bifurcation lines, Hopf, Turing, and wave bifurcation lines, shown in Figure 1.We show that in domains II and V, the sequence spotted → spot-stripelike mixtures → stripelike pattern is emerged as the growth rate of predator decreases since the Turing instability occurs in domain III, an oscillatory wave pattern emerges due to the Hopf and wave instabilities and in domain IV a typical spiral pattern is investigated because of Hopf and Turing instabilities.Moreover, we calculate numerically the values of the frequency and wavelength of the periodic oscillations in time and compare these values to the theoretical values obtained from ( 14), ( 15), (26), and (27).
The authors in [3,9,32,33] have investigated the patterns arising from special initial conditions.They indicated that spatiotemporal dynamics of a reaction-diffusion system depend on the choice of initial conditions.For discussing the evolution process of the spatial pattern formation of system (3) with some specific initial conditions, we employ two specific initial conditions which were considered in [3,9,32,33].
These figures beg the question that the spotted pattern could be a globally stable pattern even though we do not have any theoretical evidences.Similarly, we can figure out numerically that system (3) with the parameters  2 = 7 and  = 0.3, which lies on domain V, has the same phenomena even if we do not exhibit their snapshots.In this case, the system has a stripelike pattern as shown in Figure 5.