Pattern Formation in a Cross-Diffusive Holling-Tanner Model

We present a theoretical analysis of the processes of pattern formation that involves organisms distribution and their interaction of spatially distributed population with selfas well as crossdiffusion in a Holling-Tanner predator-prey model; the sufficient conditions for the Turing instability with zero-flux boundary conditions are obtained; Hopf and Turing bifurcation in a spatial domain is presented, too. Furthermore, we present novel numerical evidence of time evolution of patterns controlled by selfas well as cross-diffusion in the model, and find that the model dynamics exhibits a cross-diffusion controlled formation growth not only to spots, but also to strips, holes, and stripes-spots replication. And the methods and results in the present paper may be useful for the research of the pattern formation in the cross-diffusive model.


Introduction
A fundamental goal of theoretical ecology is to understand the interactions of individual organisms with each other and with the environment and to determine the distribution of populations and the structure of communities.Empirical evidence suggests that the spatial scale and structure of environment can influence population interactions 1 .The study of complex population dynamics is nearly as old as population ecology, starting with the pioneering work of Lotka and Volterra, a simple model of interacting species that still bears their joint names 2, 3 .
And the predator-prey system models such a phenomenon, pursuit-evasion, predators pursuing prey and prey escaping from the predators 4 .In other words, in nature, there is a tendency that the preys would keep away from predators and the escape velocity of the preys may be taken as proportional to the dispersive velocity of the predators.In the same manner, there is a tendency that the predators would get closer to the preys and the chase velocity of predators may be considered to be proportional to the dispersive velocity of the preys 5 .Keeping these in view, cross-diffusion arises, which was proposed first by Kerner 6 and first applied in competitive population system by Shigesada et al. 7 .From the pioneering work of Turing 8 , spatially continuous models formulated as reaction-diffusion equations have been intensively used to describe spatiotemporal dynamics and to investigate mechanisms for pattern formation 9 .And the appearance and evolution of these patterns have been a focus of recent research activity across several disciplines 10 .
In recent years, there has been considerable interest to investigate the stability behavior of a predator-prey system by taking into account the effect of self-as well as cross-diffusion 11-25 .Cross-diffusion expresses the population fluxes of one species due to the presence of the other species.
Furthermore, in 23 , the authors gave a numerical study of pattern formation in the Holling-Tanner model with self-and cross-diffusion as the form with the following nonzero initial conditions: and zero-flux boundary conditions: where N t , P t represent population densities of prey and predator, respectively, r is the intrinsic growth rate or biotic potential of the prey, a is the maximal predator per capita consumption rate, that is, the maximum number of prey that can be eaten by a predator in each time unit, and b is the number of preys necessary to achieve one-half of the maximum rate a 26 .And the nonnegative constants, D 12 and D 21 , called cross-diffusion coefficients, express the respective population fluxes of the prey and predators resulting from the presence of the other species, respectively.Lx and Ly give the size of the system in the directions of x and y, respectively.In 1.3 , n is the outward unit normal vector of the boundary ∂Ω, which we will assume is smooth.The main reason for choosing such boundary conditions is that we are interested in the self-organization of pattern; zero-flux conditions imply no external input 11, 27 .Biologically, cross diffusion implies countertransport 11 .And the induced crossdiffusion rate D 12 in 1.1 represents the tendency of the prey N to keep away from its predators P and D 21 represents the tendency of the predator to chase its prey.The crossdiffusion coefficients, D 12 and D 21 , may be positive or negative.Positive cross-diffusion coefficient denotes that one species tends to move in the direction of lower concentration of another species while negative cross diffusion expresses the population fluxes of one species in the direction of higher concentration of the other species 13 .
In 23 , the authors found that, for the equal self-diffusion coefficients i.e., the coefficients of ∇ 2 N and ∇ 2 P in 1.1 are both equal to 1 , only spots pattern could be obtained.In addition, they indicated that cross diffusion may have an effect on the distribution of the species, that is, it may lead the species to be isolated.
There comes a question: besides spots, does model 1.1 exhibit any other pattern replication controlled by cross-diffusion?
In this paper, based on the results of 23 , we mainly focus on the effect of self-as well as cross-diffusion on pattern formation in the two-species Holling-Tanner predator-prey model.In the next section, we give the sufficient conditions for the Turing instability with zero-flux boundary conditions.And, by using the bifurcation theory, we give the Hopf and Turing bifurcation analysis of the model.Then, we present and discuss the results of complex, not simple, pattern formations via numerical simulations, which are followed by the last section, that is, concluding remarks.

Linear Stability Analysis
In the absence of diffusion, model 1.1 has two equilibrium solutions in the positive quadrant.One equilibrium point is given by N 0 , P 0 1, 0 .In the N, P phase plane, 1, 0 is a saddle 28 .Another equilibrium point E * N * , P * depends on the parameters a and b and is given by where M 1 − a − b 2 4b.And in the presence of diffusion, we will introduce small perturbations U 1 N − N * and U 2 P − P * , where |U 1 |, |U 2 | 1.The zero-flux boundary conditions 1.3 imply that no external input is imposed from outside.To study the effect of diffusion on the model system, we have considered the linearized form of the system as follows: where where r x, y and 0 < x < Lx, 0 < y < Ly. k k n , k m and k n nπ/Lx, k m mπ/Ly are the corresponding wavenumbers.
Having substituted u nm and v nm with 2.2 , we obtain where A general solution of 2.5 has the form C 1 exp λ 1 t C 2 exp λ 2 t , where the constants C 1 and C 2 are determined by the initial conditions 1.2 and the exponents λ 1 and λ 2 are the eigenvalues of the following matrix: Correspondingly, λ 1 and λ 2 arise as the solution of the following equation:

2.11
Noting that P * 1/a 1 − N * N * b , one can obtain

2.12
If b ≥ a/ 1 a holds, dV/dt ≤ 0. Next, we select the Liapunov function for model 1.1 two dimensional with diffusion case : where ∂N ∂x ∂P ∂x ∂N ∂y ∂P ∂y dx dy.
And from 19 , we know that the equilibrium E * is Turing unstable if it is an asymptotically stable equilibrium of model 1.1 without self and cross diffusion but is unstable with respect to solutions of model 1.1 .Hence, Turing instability sets in when the condition either tr D < 0 or det D > 0 is violated, which subject to the conditions J 11 J 22 < 0 and det J > 0. It is evident that the condition tr D < 0 is not violated when the requirement J 11 J 22 < 0 is met.Hence, only violation of condition det D > 0 gives rise to Turing instability.Then the condition for Turing instability is given by The critical wavenumber k c of the first perturbations to grow is found by evaluating k m from 2.22 .
Figure 2 shows that the linear stability analysis yields the bifurcation diagram with a 0.8, b 0.1, and D 21 0.8.Turing and Hopf lines intersect at the Turing-Hopf bifurcation point D 12 , r −0.00489, 0.12572 and separate the parametric space into four domains.On domain I, located above all two bifurcation lines, the steady state is the only stable solution of the model.Domain II is the region of pure Turing instability.In domain III, which is located above all two bifurcation lines, both Hopf and Turing instability occur.And domain IV is the region of pure Hopf instability.

Pattern Formation
In this section, we performed extensive numerical simulations of the spatially extended model 1.1 in 2-dimension spaces, and the qualitative results are shown here.All our numerical simulations employ the zero-flux boundary conditions with a system size of Lx × Ly, with Lx Ly 400 discretized through x → x 0 , x 1 , x 2 , . . ., x n and y → y 0 , y 1 , y 2 , . . ., y n , with n 200.Other parameters are fixed as a 0.8, b 0.1, and D 21 0.8.The numerical integration of 1.1 was performed by means of forward Euler integration, with a time step of τ 0.05 and spatial resolution h 2 and using the standard five-point approximation for the 2D Laplacian with the zero-flux boundary conditions 30, 31 .More precisely, the concentrations N n 1 i,j , P n 1 i,j at the moment n 1 τ at the mesh position i, j are given by τg N n i,j , P n i,j ,

with the Laplacian defined by
Initially, the entire system is placed in the stationary state N * , P * 0.370156, 0.370156 , and the propagation velocity of the initial perturbation is thus on the order of 5 × 10 −4 space units per time unit.And the system is then integrated for 100000 or 300 000 time steps and some images saved.After the initial period during which the perturbation spreads, the system goes either into a time dependent state, or to an essentially steady state time independent .
In the numerical simulations, different types of dynamics are observed, and it is found that the distributions of predator and prey are always of the same type.Consequently, we can restrict our analysis of pattern formation to one distribution.In this section, we show the distribution of prey N, for instance.
Firstly, we show the pattern formation for the parameters D 12 , r located in domain II cf. Figure 2 ; the region of pure Turing instability occurs while Hopf stability occurs.We have performed a large number of simulations and found that there only exhibits spots pattern in this domain.As an example, we show the time evolution of spots pattern of prey N at 0, 10000, 30000, and 300000 iteration for D 12 , r 0.12, 0.1278 in Figure 3.In this case, one can see that for model 1.1 , the random initial distribution cf. Figure 3 a leads to the formation of regular spots except for apparently stable defects cf. Figure 3 d .
Next, we show the pattern for the parameters D 12 , r located in domain III cf. Figure 2 , both Hopf and Turing instability occur.The model dynamics exhibits

Concluding Remarks
In summary, we have investigated a cross-diffusive Holling-Tanner predator-prey model with equal self-diffusive coefficients.Based on the bifurcation analysis Hopf and Turing , we give the spatial pattern formation via numerical simulation, that is, the evolution process of the system near the coexistence equilibrium point N * , P * .
In contrast to the results in 23 , we find that the model dynamics exhibits a crossdiffusion controlled formation growth not only to spots in 23 , Sun et al. claimed that the spots pattern is the only pattern of the model , but also to stripes, holes, and stripes-spots replication.That is to say, the pattern formation of the Holling-Tanner predator-prey model is not simple, but rich and complex.
On the other hand, in the predator-prey model, predators will tend to gravitate toward higher concentrations of prey while prey will preferentially move toward regions where predators are rare.Models of predator-prey systems with cross diffusion have been extensively analyzed in the literature, though often with respect to their mathematical properties rather than to provide insight into the kinds of patterns that can emerge 24 .And the methods and results in the present paper may be useful for the research of the pattern formation in the cross-diffusive predator-prey model.
D 12 D 21 .Following Malchow et al. 29 , we can know that any solution of the system 2.2 can be expanded into a Fourier series so that D 12 D 21 k 4 −J 11 − J 22 D 12 J 21 D 21 J 12 k 2 det J .If 0 < D 12 , D 21 < 1, and b ≥ a/ 1 a , then the uniform steady state E * of model 1.1 is globally asymptotically stable.Proof.For global stability of nonspatial model of 1.1 , we select a Liapunov function: Substituting the value of dN/dt and dP/dt from the nonspatial model of 1.1 , we obtained Equation 2.20 leads to the following final criterion for Turing instability: J 11 J 22 − D 12 J 21 − D 21 J 12 2 > 4 1 − D 12 D 21 det J .
min det J − J 11 J 22 − D 12 J 21 − D 21 J 12 2 4 1 − D 12 D 21 Theorem 2.2.The equilibrium E * of model 1.1 is Turing instability if J 11 J 22 > D 12 J 21 D 21 J 12 and J 11 J 22 − D 12 J 21 − D 21 J 12 > 2 1 − D 12 D 21 det J .It is easy to see that the minimum of H k 2 occurs at k 2 k 2 m , where