Spatiotemporal Complexity of a Leslie-Gower Predator-Prey Model with the Weak Allee Effect

Abstract

We investigate a diffusive Leslie-Gower predator-prey model with the additive Allee effect on prey subject to the zero-flux boundary conditions. Some results of solutions to this model and its corresponding steady-state problem are shown. More precisely, we give the stability of the positive constant steady-state solution, the refined a priori estimates of positive solution, and the nonexistence and existence of the positive nonconstant solutions. We carry out the analytical study for two-dimensional system in detail and find out the certain conditions for Turing instability. Furthermore, we perform numerical simulations and show that the model exhibits a transition from stripe-spot mixtures growth to isolated spots and also to stripes. These results show that the impact of the Allee effect essentially increases the model spatiotemporal complexity.

(ProQuest: ... denotes non-US-ASCII text omitted.)

Yongli Cai 1 and Caidi Zhao 2 and Weiming Wang 2

Academic Editor:Victor Kazantsev

1, School of Mathematics and Computational Science, Sun Yat-Sen University, Guangzhou 510275, China
2, College of Mathematics and Information Science, Wenzhou University, Wenzhou 325035, China

Received 23 April 2013; Revised 27 November 2013; Accepted 28 November 2013

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.

1. Introduction

The dynamics of a predator-prey model in a homogeneous environment can be described by the following reaction-diffusion equations: [figure omitted; refer to PDF] where u ( x , t ) and v ( x , t ) are the densities of prey and predator at time t and position x ∈ Ω ⊂ ^{... m} , respectively. The Laplace operator Δ describes the spatial dispersal with passive diffusion; _{D 1} > 0 and _{D 2} > 0 are the diffusion coefficients corresponding to species u and v . F ( u ) describes the per-capita growth rate of the prey; G ( u ) is the functional response of the predator, which corresponds to the saturation of their appetites and reproductive capacity; Q ( u , v ) , the so-called numerical response, is the per-capita growth rate of the predator [1-4].

Functions F ( u ) , G ( u ) , and Q ( u , v ) can be formulated in various specific situations. In general, F ( u ) is of the standard logistic growth: [figure omitted; refer to PDF] which was first created by Verhulst [5]. Here K > 0 is the prey carrying capacity and _{r 1} > 0 is the intrinsic growth rate of prey.

Some conventional functional response functions G ( u ) include Holling types I, II, and III (see [6-10]). Among many possible choices of G ( u ) , the Holling type-II functional response is most commonly used in the ecological literature, which is defined by [11]: [figure omitted; refer to PDF] where _{c 1} > 0 describes the maximum predation rate and _{K 1} ...5; 0 measures the extent to which environment provides protection to prey u . The Leslie-Gower type numerical response Q ( u , v ) is given by [figure omitted; refer to PDF] which was first proposed by Leslie [12], and has been discussed by Leslie and Gower [13] and Pielou [14]. A modified version of Leslie-Gower functional response [figure omitted; refer to PDF] is given by Aziz-Alaoui et al. [15, 16]. Here, _{r 2} > 0 describes the growth rate of the predator v ; _{c 2} > 0 has a similar meaning to _{c 1} ; δ > 0 takes on the role of the prey-dependent carrying capacity for the predator; _{K 2} ...5; 0 is the extent to which environment provides protection to predator v . Hence, we can rewrite model (1) as follows: [figure omitted; refer to PDF] The biological significance of all parameters in model (6) is as above.

For model (6), in the case of _{K 1} = _{K 2} = 0 , Du et al. [17, 18] mainly focused attention on the steady-state problem and observed some quite interesting phenomena of pattern formation. In the case of _{K 1} > 0 , _{K 2} = 0 , the so-called Holling-Tanner model, Peng and Wang [19, 20] analyzed the global stability of the unique positive constant steady-state and established the results for the existence and nonexistence of positive nonconstant steadystates; Shi and coworkers [21] studied the global attractor and persistence property, local and global asymptotic stability of the unique positive constant equilibrium, and the existence and nonexistence of nonconstant positive steady-states; Li et al. [22] considered the Turing and Hopf bifurcations of the equilibrium solutions; Liu and Xue [23] found the model exhibits the spotted, black-eye, and labyrinthine patterns. For model (6), that is, _{K 1} ...0; _{K 2} ...0; 0 , Camara and Aziz-Alaoui [24-26] paid more attention to pattern formation in the spatial domain and observed the labyrinth, chaos, and spiral wave patterns.

On the other hand, in population dynamics, any mechanism that can lead to a positive relationship between a component of individual fitness and either the number or density of conspecific can be termed an Allee effect [27-30], starting with the pioneer work of ecologist Allee [31]. In particular, theoretical investigations have shown that an Allee effect can greatly increase the likelihood of local and global extinction [32] and can lead to a rich variety of dynamical effects. As a consequence, it is necessary to explore the influence of Allee effect in the growth of a population.

The Allee effect has been modeled in different ways [33-37]. From an ecological viewpoint, the Allee effect has been modeled into strong and weak ones [33, 38-42]. In a recent analytic approach by Wang and Kot [38], the Allee effect is "strong" if the sign of the growth function u F ( u ) in the limit of law density is negative; that is, [figure omitted; refer to PDF] It is "weak" if the sign of the growth function u F ( u ) in the limit of law density is positive; that is, [figure omitted; refer to PDF] The strong Allee effect introduces a population threshold, and the population must surpass this threshold to grow. In contrast, the weak case has not any threshold [10, 35, 38, 42].

In particular, the growth function considering Allee effect is expressed by the equation: [figure omitted; refer to PDF] having an additive Allee effect, which was first deduced in [43] and applied in [34-36]. Where _{m 1} / ( u + _{b 1} ) is the term of additive Allee effect, _{m 1} and _{b 1} are the Allee effect constants. It should be noted that, if _{m 1} < _{b 1} , the Allee effect in (9) is weak, while if _{m 1} > _{b 1} , the Allee effect in (9) is strong.

Based on the above discussions, in this paper, we rigorously consider the spatiotemporal dynamics of the following modified Leslie-Gower predation model with the additive Allee effect on prey: [figure omitted; refer to PDF]

We make a change of variables: [figure omitted; refer to PDF] And for the sake of convenience, we still use variables u , v instead of u ~ , v ~ . Thus, the model to be studied is as follows: [figure omitted; refer to PDF] where [figure omitted; refer to PDF] Here, s is the growth rate of the predator v . m u / ( u + b ) is the term of additive Allee effect, and m and b are the Allee effect constants. Ω ⊂ ^{... 2} is a bounded domain with smooth boundary ∂ Ω , and ν is the outward unit normal vector on ∂ Ω . The initial data _{u 0} ( x ) and _{v 0} ( x ) are continuous functions on Ω ¯ , and the zero-flux boundary conditions mean that model (12) is self-contained and has no population flux across the boundary ∂ Ω [44, 45].

By the standard theory for semilinear parabolic systems (see, [46]), we have model (12) that admits a unique classical solution ( u ( x , t ) , v ( x , t ) ) for all time.

The stationary problem of model (12), which may display the dynamical behavior of solutions to model (12) as time goes to infinity, satisfies the following elliptic system: [figure omitted; refer to PDF] Unless otherwise specified, in this paper, we always assume that m < b ; that is, we only focus on the case of weak Allee effect.

The rest of the paper is organized as follows. In Section 2, we investigate the stability of nonnegative constant steady-state solutions. In Section 3, we mainly give a priori upper and lower bounds for positive solutions of model (14). In Section 4 we discuss existence and nonexistence of nonconstant positive solutions, which might give us some suggestions on the conditions under which the patterns may or may not occur. In Section 5 we first use the method of linearized stability analysis to deduce the conditions under which the Turing instability might occur, and next we perform a series of numerical simulations to show the occurrence of different patterns. Finally, in the last section we make a summary to our results and give some concluding remarks.

2. Dynamics Analysis of Model (12)

2.1. The Existence of the Constant Steady-State Solution

It is easy to verify that model (12) has the following nonnegative constant steady-state solutions:

(i) the trivial constant solution _{E 0} = ( 0,0 ) (extinction of two species);

(ii) the semitrivial constant solution _{E 1} = ( 0 , _{k 2} ) (extinction of the prey);

(iii): the semitrivial constant solution _{E 2} = ( ( 1 - b + ( 1 - b ^{) 2} - 4 ( m - b ) ) / 2 , 0 ) (extinction of the predator);

(iv) the unique positive constant solution _{E 3} = ( ^{u *} , ^{v *} ) (coexistence of two species), where ^{v *} = ^{u *} + _{k 2} , and ^{u *} is a real positive root of the cubic [figure omitted; refer to PDF] where 3 _{η 1} = b + c + _{k 1} - 1 , 3 _{η 2} = b c + b _{k 1} + c _{k 2} + m - b - _{k 1} , _{η 3} = b c _{k 2} + _{k 1} ( m - b ) .

By the transformation z = w + _{η 1} , (15) is reduced to [figure omitted; refer to PDF] where [figure omitted; refer to PDF]

It is worthy to note that if q = _{η 3} - 3 _{η 1 η 2} + 2 ^{η 1 3} < 0 holds, (16) has a real positive root. Considering (17) and (18), from q = 0 , one can determine m = ^{m *} , where [figure omitted; refer to PDF]

Hence, we have the following lemma regarding the existence of the positive constant steady-state solution of model (12).

Lemma 1.

If either of the following inequalities holds: [figure omitted; refer to PDF] model (12) has a unique positive constant steady-state solution _{E 3} = ( ^{u *} , ^{v *} ) .

2.2. Stability of the Constant Steady-State Solution

In this subsection, we will analyze the asymptotical stability of the nonnegative constant solutions for the corresponding reaction-diffusion dynamics (12).

For sake of simplicity, we rewrite model (12) in the vectorial form: [figure omitted; refer to PDF] where w = ( u , v ^{) T} , D = diag ( _{d 1} , _{d 2} ) , and [figure omitted; refer to PDF]

Let 0 = _{μ 0} < _{μ 1} < _{μ 2} < ... be the eigenvalues of the operator- Δ on Ω with the zero-flux boundary conditions. And set [figure omitted; refer to PDF] { _{[varphi] i j} |" j = 1 , ... , dim E ( _{μ i} ) } being an orthonormal basis of E ( _{μ i} ) , and _{X i j} = { c _{[varphi] i j} |" c ∈ ^{... 2} } ; then [figure omitted; refer to PDF] where _{X i} = ^{... j = 1 dim E ( _{μ i} )} _{X i j} .

Let E = ( u , v ) be any arbitrary constant steady-state solution of model (12). And the linearization of model (12) at the constant steady-state solution can be expressed by [figure omitted; refer to PDF] where _{J ( u , v )} = ( _{ξ 11 ξ 12 ξ 21 ξ 22} ) and [figure omitted; refer to PDF]

From [46], it is known that if all the eigenvalues of the operator £ have negative real parts, then E = ( u , v ) is asymptotically stable; if there is an eigenvalue with positive real part, then E = ( u , v ) is unstable; if all the eigenvalues have nonpositive real parts while some eigenvalues have zero real parts, then the stability of E = ( u , v ) cannot be determined by the linearization [10].

For each i ...5; 0 , _{X i} is invariant under the operator £ , and λ is an eigenvalue of £ if and only if λ is an eigenvalue of the matrix _{M i} = - _{μ i} D + _{J ( u , v )} for some i ...5; 0 .

In the following, we denote ^{ξ j j [ k ]} = _{ξ j j} evaluated at _{E k} , j = 1,2 and k = 0 , 1 , 2 , 3 . So, the local stability of the constant steady-state solution can be analyzed as follows.

Theorem 2.

For any positive parameters,

(a) the trivial constant solution _{E 0} = ( 0 , 0 ) is unstable;

(b) the semitrivial constant solution _{E 1} = ( 0 , _{k 2} ) is

(b1) locally asymptotically stable if c > ( _{k 1} / _{k 2} ) or ( _{k 1} / _{k 2} ) ( 1 - m / b ) < c ...4; ( _{k 1} / _{k 2} ) holds

(b2) unstable if c < ( _{k 1} / _{k 2} ) ( 1 - m / b ) holds;

(c) the semitrivial constant solution _{E 2} = ( ( 1 - b + ( 1 - b ^{) 2} - 4 ( m - b ) ) / 2 , 0 ) is unstable.

Proof.

The stability of the constant steady-state solution is reduced to consider the characteristic equation: [figure omitted; refer to PDF] with [figure omitted; refer to PDF]

(a) _{J E 0} = ( 1 - ( b / m ) 0 0 s ) , for i = 0 ; the eigenvalues are 1 - ( b / m ) > 0 and s > 0 , so _{E 0} = ( 0 , 0 ) is unstable.

(b) We can obtain _{J E 1} = ( 1 - ( b / m ) - ( c _{k 2} / _{k 1} ) 0 s - s ) .

(b1) If _{k 1} < c _{k 2} or _{k 1} ...5; c _{k 2} and b ( 1 - ( c _{k 2} / _{k 1} ) ) < m < b hold, then 1 - ( b / m ) - ( c _{k 2} / _{k 1} ) < 0 , so for i ...5; 0 , [figure omitted; refer to PDF]

: Hence, _{E 1} is locally asymptotically stable.

(b2) When _{k 1} ...5; c _{k 2} and 0 < m < b ( 1 - ( c _{k 2} / _{k 1} ) ) , then 1 - ( b / m ) - ( c _{k 2} / _{k 1} ) > 0 . For i ...5; 0 , det ( _{M i} ) = - s ( 1 - ( b / m ) - ( c _{k 2} / _{k 1} ) ) < 0 , which implies that (27) has at least one root with positive real part. Hence, _{E 1} is an unstable steady-state solution of model (12).

(c) Since ^{ξ 11 [ 2 ]} = - ( ( 2 ( 1 + b ^{) 2} - 4 m ( ( 1 + b ^{) 2} - 4 m + 1 + b - 2 m ) ) / ( b + 1 + ( 1 + b ^{) 2} - 4 m ^{) 2} ) , ^{ξ 22 [ 2 ]} = s , for i = 0 , one of the eigenvalues is s > 0 , so _{E 2} is unstable.

The proof is complete.

Straightforward calculations show that [figure omitted; refer to PDF] where ^{ξ 11 [ 3 ]} = ( ( c ^{u *} ( ^{u *} + _{k 2} ) ) / ^{( u * + _{k 1} ) 2} ) + ( m ^{u *} / ^{( u * + b ) 2} ) - ^{u *} .

The determinant of _{J E 3} is given by [figure omitted; refer to PDF] then, the sign of det ( _{J E 3} ) depends on the factor [varphi] ( ^{u *} ) : [figure omitted; refer to PDF] where ψ ( · ) is the same definition as (15), and [figure omitted; refer to PDF]

Therefore, we have the following.

Theorem 3.

Assume that ^{ξ 11 [ 3 ]} < s , _{ρ 1} ^{u * 2} + _{ρ 2} ^{u *} + _{ρ 3} > 0 , and the first eigenvalue _{μ 1} subject to the zero-flux boundary conditions satisfies [figure omitted; refer to PDF] Then the positive constant steady-state solution _{E 3} = ( ^{u *} , ^{u *} + _{k 2} ) of model (15) is uniformly asymptotically stable.

Proof.

When ^{ξ 11 [ 3 ]} < s , _{ρ 1} ^{u * 2} + _{ρ 2} ^{u *} + _{ρ 3} > 0 , then tr ( _{J E 3} ) < 0 and det ( _{J E 3} ) > 0 . So, for i ...5; 0 , [figure omitted; refer to PDF] Note that for any i ...5; 0 , we have det ( _{A i} ) > 0 > tr ( _{A i} ) . Therefore, the eigenvalues of the matrix - _{μ i} D + _{J E 3} have negative real parts. It thus follows from the Routh-Hurwitz criterion that, for each i ...5; 0 , the two roots _{λ i 1} and _{λ i 2} of _{[straight phi] i} ( λ ) = 0 all have negative real parts.

In the following, we prove that there exists δ > 0 such that [figure omitted; refer to PDF]

Let λ = _{μ i} ξ ; then [figure omitted; refer to PDF] Since _{μ i} [arrow right] ∞ as i [arrow right] ∞ , it follows that [figure omitted; refer to PDF]

By the Routh-Hurwitz criterion, it follows that the two roots _{ξ 1} , _{ξ 2} of _{[straight phi] ~ i} ( λ ) = 0 all have negative real parts. Thus, there exists a positive constant d ~ , such that Re { _{ξ 1} } , Re { _{ξ 2} } ...4; - d ~ . By continuity, we see that there exists _{i 0} such that the two roots _{ξ i 1} , _{ξ i 2} of _{[straight phi] ~ i} ( λ ) = 0 satisfy Re { _{ξ i 1} } ...4; - d ~ / 2 , Re { _{ξ i 2} } ...4; - d ~ / 2 , ∀ i ...5; _{i 0} . In turn, Re { _{λ i 1} } , Re { _{λ i 2} } ...4; - _{μ i} d ~ / 2 ...4; - d ~ / 2 , ∀ i ...5; _{i 0} . Let [figure omitted; refer to PDF] Then δ ~ > 0 and (36) holds for δ = min { δ ~ , ( d ~ / 2 ) } .

Consequently, the spectrum of £ , which consists of eigenvalues, lies in { Re λ ...4; - δ } . In the sense of [46], we obtain that the positive constant steady-state solution _{E 3} = ( ^{u *} , ^{u *} + _{k 2} ) of model (12) is uniformly asymptotically stable. This ends the proof.

3. A Priori Estimates

In this section, we give a priori estimates for the steady-state solutions of model (14). To prove that we recall the following maximum principle [47].

Lemma 4 (see [47, maximum principle]).

Let Ω be a bounded Lipschitz domain in ^{... n} and g ∈ C ( Ω ¯ × ... ) .

(a) Assume that w ∈ ^{C 2} ( Ω ) ∩ ^{C 1} ( Ω ¯ ) and satisfies [figure omitted; refer to PDF] If w ( _{x 0} ) = _{max Ω ¯} w ( x ) , then g ( _{x 0} , w ( _{x 0} ) ) ...5; 0 .

(b) Assume that w ∈ ^{C 2} ( Ω ) ∩ ^{C 1} ( Ω ¯ ) and satisfies [figure omitted; refer to PDF] If w ( _{x 0} ) = _{min Ω ¯} w ( x ) , then g ( _{x 0} , w ( _{x 0} ) ) ...4; 0 .

Theorem 5.

Suppose that c ( 1 + _{k 2} ) < _{k 1} ( 1 - ( m / b ) ) . Any positive solution ( u ( x ) , v ( x ) ) of model (14) satisfies [figure omitted; refer to PDF] where [figure omitted; refer to PDF]

Proof.

Assume that ( u ( x ) , v ( x ) ) is a positive solution of (14). Set [figure omitted; refer to PDF] Applying Lemma 4 to model (14), we obtain that [figure omitted; refer to PDF] By virtue of the definitions of ^{x i 0} ( i = 1,2 , 3,4 ), it follows from (45) that u ( ^{x 1 0} ) < 1 and v ( ^{x 2 0} ) < 1 + _{k 2} , and [figure omitted; refer to PDF] So, we have [figure omitted; refer to PDF] If c ( 1 + _{k 2} ) < _{k 1} ( 1 - ( m / b ) ) or c ( 1 + _{k 2} ) = _{k 1} ( 1 - ( m / b ) ) and c ( 1 + _{k 2} ) > ^{k 2 2} hold, from (48), we get that u ( ^{x 3 0} ) > B .

If c _{k 2} < _{k 1} ( 1 - ( m / ( 1 + b ) ) ) or c _{k 2} = _{k 1} ( 1 - ( m / ( 1 + b ) ) ) and c _{k 2} > ^{k 2 2} hold, from (47), we get that [figure omitted; refer to PDF] By simple computations, c ( 1 + _{k 2} ) < _{k 1} ( 1 - ( m / b ) ) indicates that c _{k 2} < _{k 1} ( 1 - ( m / ( 1 + b ) ) ) . So, if c ( 1 + _{k 2} ) < _{k 1} ( 1 - ( m / b ) ) holds, we can obtain B < u ( x ) < A and ( 1 + _{k 2} ) B < v ( x ) < ( 1 + _{k 2} ) A . The proof is complete.

In order to obtain the desired bounds, we need to use the following Harnack inequality due to [48].

Lemma 6 (see [48, Harnack inequality]).

Let w ∈ ^{C 2} ( Ω ) ∩ ^{C 1} ( Ω ¯ ) be a positive solution to Δ w ( x ) + c ( x ) w ( x ) = 0 , where c ∈ C ( Ω ¯ ) , satisfying zero-flux boundary conditions. Then there exists a positive constant ^{C *} = ^{C *} ( _{|| c || ∞} , Ω ) , such that [figure omitted; refer to PDF]

Theorem 7.

Let ^{d *} be an arbitrary fixed positive number. Then there exist positive C _ = C _ ( m , b , c , _{k 1} , _{k 2} , Ω ) , such that if _{d 1} ...5; ^{d *} and _{d 2} > 0 , any positive solution ( u , v ) of model (14) satisfies [figure omitted; refer to PDF]

Proof.

By Theorem 5, we note that u ( x ) < 1 and v ( x ) < 1 + _{k 2} . And it suffices to verify the lower bounds of ( u , v ) . We will verify the conclusion by a contradiction argument.

On the contrary, suppose that the conclusion is not true; then there exist sequences { _{d 1 , i} ^{} i = 1 ∞} and { _{d 2 , i} ^{} i = 1 ∞} with _{d 1 , i} ...5; d , _{d 2 , i} > 0 and the positive solution ( _{u i} , _{v i} ) of model (14) corresponding to ( _{d 1} , _{d 2} ) = ( _{d 1 , i} , _{d 2 , i} ) , such that [figure omitted; refer to PDF] and ( _{u i} , _{v i} ) satisfies [figure omitted; refer to PDF] We observe that Lemma 4 guarantees [figure omitted; refer to PDF] for all i ...5; 1 by use of the second equation of (53). On the other hand, by the Harnack inequality (Lemma 6), we know that there is a positive constant C independent of i , such that _{max Ω ¯ u i} ( x ) ...4; C _{min Ω ¯ u i} ( x ) for all i ...5; 1 . Consequently, [figure omitted; refer to PDF] which contradicts Theorem 5. We finish the proof.

4. Nonexistence and Existence of the Nonconstant Solutions

In this section, we apply the energy method and the topological degree theory [49] to establish some results on the nonexistence and existence of the positive nonconstant solutions of model (14), respectively.

4.1. Nonexistence of the Nonconstant Solutions

4.1.1. Nonexistence of Nonconstant Semitrivial Steady-State Solution

In the case that model (14) has the semitrivial steady-state solutions ( u ( x ) , 0 ) and ( 0 , v ( x ) ) , such ( u ( x ) , 0 ) and ( 0 , v ( x ) ) should satisfy, respectively, [figure omitted; refer to PDF]

Note that _{μ 1} is the smallest positive eigenvalue of the operator- Δ in Ω subject to the zero-flux boundary conditions. Now, using the energy estimates, for the semitrivial solution, we can claim the following.

Theorem 8.

(i) If _{d 1} > 1 / _{μ 1} , the only nonnegative solutions of model (56) are u = 0 and u = ( ( 1 - b + ( 1 - b ^{) 2} - 4 ( m - b ) ) / 2 ) .

(ii) If _{d 2} > s / _{μ 1} , the only nonnegative solution of model (57) is v = _{k 2} .

Proof.

Let u ( x ) , v ( x ) be a nonnegative solution of models (56) and (57), respectively. Denote u ¯ = ^{| Ω | - 1} _{∫ Ω} ... u d x and v ¯ = ^{| Ω | - 1} _{∫ Ω} ... v d x . Then _{∫ Ω} ... ( u - u ¯ ) d x = _{∫ Ω} ... ( v - v ¯ ) d x = 0 .

Multiplying the equation in (56) by u - u ¯ , we get [figure omitted; refer to PDF] Then with the Poincaré inequality [figure omitted; refer to PDF] we find that [figure omitted; refer to PDF] which implies that _{d 1} < 1 / _{μ 1} unless u - u ¯ ...1; 0 , and u must be a constant solution.

In a similar manner, we multiply the equation in (57) by ( v - v ¯ ) to have [figure omitted; refer to PDF] In view of _{d 2} > s / _{μ 1} , we have v - v ¯ ...1; 0 and v must be a constant solution. This ends the proof.

4.1.2. Nonexistence of the Positive Nonconstant Solutions

This subsection is devoted to the consideration of the nonexistence for the positive nonconstant solutions of model (12), and in the below results, the diffusion coefficients do play a significant role.

Theorem 9.

There exists a positive constant ^{d *} = ^{d *} ( b , m , c , s , _{k 1} , _{k 2} , Ω ) such that model (14) has no positive nonconstant solution provided that min { _{d 1} , _{d 2} } > ^{d *} .

Proof.

Let ( u , v ) be any positive solution of model (14). Then, multiplying the first equation of model (14) by ( u - u ¯ ) , integrating over Ω and using Theorem 7, we have that [figure omitted; refer to PDF] In a similar manner, we multiply the second equation in model (14) by ( v - v ¯ ) to have [figure omitted; refer to PDF] Using Theorem 7 again, we have that [figure omitted; refer to PDF] Similarly, [figure omitted; refer to PDF] From the well-known Poincaré inequality (62), (63), (64), and (65), we obtain that [figure omitted; refer to PDF] where _{α 1} = 1 + ( c ( 3 + 2 _{k 2} ) / 2 _{k 1} ) + ( s ( 1 + _{k 2} ^{) 2} / 2 ^{k 1 2} ) , _{α 2} = s + ( c / 2 _{k 1} ) + ( s ( 1 + _{k 2} ^{) 2} / 2 ^{k 1 2} ) .

This shows that if [figure omitted; refer to PDF] then ∇ ( u - u ¯ ) = ∇ ( v - v ¯ ) = 0 , which asserts our results.

4.2. Existence of the Nonconstant Positive Solutions

In this subsection, we will discuss the existence of the positive nonconstant solution of model (14).

Unless otherwise specified, in this subsection, we always require that Lemma 1 holds, which guarantees that model (14) has the unique positive constant solution _{E 3} = ( ^{u *} , ^{u *} + _{k 2} ) . From now on, we denote _{w 0} = _{E 3} .

Let X be the space defined in (23) and let [figure omitted; refer to PDF] We write model (14) in the form: [figure omitted; refer to PDF] where [figure omitted; refer to PDF] Then w is a positive solution of model (69) if and only if w satisfies [figure omitted; refer to PDF] where ( I - Δ ^{) - 1} is the inverse operator of I - Δ subject to the zero-flux boundary condition. Then [figure omitted; refer to PDF] where [figure omitted; refer to PDF]

If ∇ ... ( _{w 0} ) is invertible, by Theorem 2 . 8.1 of [50] the index of ... at _{w 0} is given by [figure omitted; refer to PDF] where γ is the multiplicity of negative eigenvalues of ∇ ... ( _{w 0} ) .

On the other hand, using the decomposition (24), we have that _{X i} is an invariant space under ∇ ... ( _{w 0} ) and ξ ∈ ... is an eigenvalue of ∇ ... ( _{w 0} ) in _{X i} , if and only if ξ is an eigenvalue of ( _{μ i} + 1 ^{) - 1} ( _{μ i} I - ...9C; ) . Therefore, ∇ ... ( _{w 0} ) is invertible, if and only if for any i ...5; 0 the matrix _{μ i} I - ...9C; is invertible.

Let m ( _{u i} ) be the multiplicity of _{μ i} . For the sake of convenience, we denote [figure omitted; refer to PDF] Then, if _{μ i} I - ...9C; is invertible for any i ...5; 0 , with the same arguments as in [51], we can assert the following conclusion.

Lemma 10.

Assume that, for all i ...5; 0 , the matrix _{μ i} I - ...9C; is nonsingular; then [figure omitted; refer to PDF]

To compute index ( ... , _{w 0} ) , we have to consider the sign of H ( μ ) . A straightforward computation yields [figure omitted; refer to PDF] where _{θ 1} = ( ^{ξ 11 [ 3 ]} / _{d 1} ) - ( s / _{d 2} ) , _{θ 2} = ( s / _{d 1 d 2} ) ( ( c ^{u *} / ( ^{u *} + _{k 1} ) ) - ^{ξ 11 [ 3 ]} ) .

If ^{θ 1 2} - 4 _{θ 2} > 0 , then H ( _{d 2} , _{d 2} , μ ) = 0 has two positive solutions ^{μ ±} given by [figure omitted; refer to PDF]

Theorem 11.

Assume that ^{θ 1 2} - 4 _{θ 2} > 0 . If ^{μ -} ∈ ( _{μ i} , _{μ i + 1} ) and ^{μ +} ∈ ( _{μ j} , _{μ j + 1} ) for some 0 ...4; i < j , and ^{∑ k = i + 1 j} ... m ( _{u k} ) is odd, then model (14) has at least one nonconstant solution.

Proof.

By Theorem 9, we can fixed _{d ¯ 1} > _{d 1} and _{d ¯ 2} > _{d 2} such that

(i) model (12) with diffusion coefficients _{d ¯ 1} and _{d ¯ 2} has no nonconstant solutions;

(ii) H ( _{d ¯ 1} , _{d ¯ 2} , μ ) > 0 for all μ ...5; 0 .

By virtue of Theorem 7, there exists a positive constant M = M ( m , b , c , _{k 1} , _{k 2} ) such that ^{M - 1} < u , v < M .

Set [figure omitted; refer to PDF] and define [figure omitted; refer to PDF] by [figure omitted; refer to PDF] where [figure omitted; refer to PDF] It is clear that solving model (14) is equivalent to finding a fixed point of Ψ ( w , 1 ) in [physics M-matrix] . Further, by virtue of the definition of [physics M-matrix] , we have that Ψ ( w , 1 ) has no fixed point in ∂ [physics M-matrix] for all 0 ...4; t ...4; 1 .

Since Ψ ( w , t ) is compact, the Leray-Schauder topological degree deg ( I - Ψ ( w , t ) , [physics M-matrix] , 0 ) is well defined. From the invariance of Leray-Schauder degree at the homotopy, we deduce [figure omitted; refer to PDF]

Clearly, I - Ψ ( w , 1 ) = ... . Thus, if model (14) has no other solutions except the constant one _{w 0} , then Lemma 10 shows that [figure omitted; refer to PDF]

On the contrary, by the choice of _{d ¯ 1} and _{d ¯ 2} , we have that _{w 0} is the only solution of Ψ ( w , 0 ) = w . Furthermore, by (ii) above, we have [figure omitted; refer to PDF]

From (83)-(85), we get a contradiction, and the proof is completed.

Corollary 12.

Let _{d 2} be fixed. If ^{ξ 11 [ 3 ]} > 0 and all the eigenvalues _{μ i} have odd multiplicity, then, there exists a sequence of intervals { ( _{k n} , _{K n} ) } with 0 < _{k n} < _{K n} < _{k n + 1} [arrow right] 0 (as n [arrow right] ∞ ) such that the steady-state model (14) has at least one nonconstant solution for all _{d 1} ∈ _{... n ...5; 1} ... ( _{k n} , _{K n} ) .

Corollary 13.

Let _{d 1} be fixed. If ^{ξ 11 [ 3 ]} > 0 and _{∑ i ...5; 0,0 < μ i < ^{ξ 11 [ 3 ]}} ... m ( _{μ i} ) is odd, then there exists ^{D *} > 0 such that the steady-state model (14) has at least one nonconstant solution for any _{d 2} > ^{D *} .

We omit the proofs of Corollaries 12 and 13 here and refer the reader to more detailed proofs in [52].

Remark 14.

Our results suggest that if the parameters are properly chosen, both the general stationary pattern and more interesting Turing pattern can arise as a result of diffusion.

5. Turing Instability and Pattern Formation

5.1. Turing Instability

In this subsection, we mainly focus on the effects of reaction-diffusion on Turing instability for model (12).

Let us consider the spatially homogeneous system corresponding to model (12): [figure omitted; refer to PDF]

Mathematically speaking, a positive constant steady-state solution _{E 3} is Turing instability, which was emphasized by Turing in his pioneering work in 1952 [53], meaning that it is an asymptotically stable steady-state solution of model (86) but is unstable with respect to the solutions of spatial model (12). Since ^{ξ 11 [ 33 ]} - s < 0 , then tr ( _{A i} ) < 0 (the matrix _{A i} = - _{μ i} D + _{J E 3} ) is always true. Hence _{A i} has an eigenvalue with a positive real part; then it must be a real value and the other eigenvalue must be a negative real one. As a consequence, we can obtain the following results.

Theorem 15.

Assume that the following conditions are true:

(a) ^{ξ 11 [ 3 ]} < s ;

(b) _{η 1} ^{u * 2} + _{η 2} ^{u *} + _{μ 3} > 0 ;

(c) _{d 2} ^{ξ 11 [ 3 ]} - _{d 1} s > 2 _{d 1 d 2} det _{J E 3} ;

then the positive constant steady-state solution _{E 3} = ( ^{u *} , ^{u *} + _{k 2} ) of model (12) is Turing instability if 0 < _{k 1} < _{μ i} < _{k 2} for some _{μ i} , where [figure omitted; refer to PDF]

Proof.

Using the same approach as in Theorem 3, we have that the positive constant steady state solution _{E 3} of model (86) is asymptotically stable provided (a) and (b).

For the Turing instability, we must have det ( _{M i} ) < 0 for some _{μ i} . And we notice that det ( _{M i} ) achieves its minimum: [figure omitted; refer to PDF] at the critical value ^{μ *} > 0 when ^{μ *} = ( ( _{d 2} ^{ξ 11 [ 3 ]} - _{d 1} s ) / 2 _{d 1 d 2} ) .

We now discuss the conditions for Turing instability in the further calculation. The condition _{min μ i} det ( _{M i} ) < 0 is equivalent to 4 _{d 1 d 2} det ( _{J E e 2} ) - ( _{d 2} ^{ξ 11 [ 3 ]} - _{d 1} s ^{) 2} < 0 . In this case, det ( _{M i} ) = 0 has two positive roots _{k 1} and _{k 2} which satisfy [figure omitted; refer to PDF] Therefore, if we can find some _{μ i} such that _{k 1} < _{μ i} < _{k 2} , then det ( _{A i} ) < 0 . By [46] it follows that _{E 3} is uniformly asymptotically instable. This finishes the proof.

5.2. Pattern Formation

In this section, we perform extensive numerical simulations of the spatially extended model (12) in two-dimensional space, and the qualitative results are shown here. All our numerical simulations employ the zero-flux boundary conditions with a system size of 100 × 100 .

The numerical integration of model (12) is performed by using a finite difference approximation for the spatial derivatives and an explicit Euler method for the time integration [54] with a time stepsize of 1 / 100 . The initial condition is always a small amplitude random perturbation around the positive constant steady-state solution _{E 3} = ( ^{u *} , ^{u *} + _{k 2} ) . After the initial period during which the perturbation spread, the model goes either into a time-dependent state, or to an essentially steady-state solution (time independent).

We use the standard five-point approximation [55] for the 2D Laplacian with the zero-flux boundary conditions. More precisely, the concentrations ( ^{u i , j n + 1} , ^{v i , j n + 1} ) at the moment ( n + 1 ) τ at the mesh position ( _{x i} , _{y j} ) are given by [figure omitted; refer to PDF] with the Laplacian defined by [figure omitted; refer to PDF] where f ( u , v ) = u ( 1 - u - ( m / ( u + b ) ) ) - ( c u v / ( u + _{k 1} ) ) and g ( u , v ) = s v ( 1 - ( v / ( u + _{k 2} ) ) ) , and the space stepsize h = 1 / 3 .

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 u , for instance. We have taken some snapshots with red (blue) corresponding to the high (low) value of prey u .

Now, we show the Turing pattern for the different values of the parameters. Via numerical simulation, one can see that the model dynamics exhibits spatiotemporal complexity of pattern formation, including spots, stripes, and holes Turing patterns.

In Figure 1, we show the time process of spots pattern formation of the prey u at t = 0,300,1000,3000 for b = 0.495 , c = 0.5 , _{k 1} = 0.3 , _{k 2} = 0.2 , m = 0.32 , s = 0.25 , _{d 1} = 0.1 , and _{d 2} = 2 . In this case, there appears a competition between stripes and spots. The pattern takes a long time to settle down; starting with a homogeneous state _{E 3} = ( 0.066566 , 0.266566 ) (cf., Figure 1(a)), the random perturbations lead to the formation of stripes and spots (cf., Figure 1(b)), and the later random perturbations make these stripes decay and end with the time-independent regular spots pattern (cf., Figure 1(d)), which is isolated zones with low prey densities. Ecologically, spots pattern shows that the prey population is driven by predators to a very high level in those regions. The final result is the formation of patches of high prey density surrounded by areas of high prey densities; that is to say, under the control of these parameters, the prey is predominant in the area.

The formation of spots patterns of prey u with parameters b = 0.495 , c = 0.5 , _{k 1} = 0.3 , _{k 2} = 0.2 , m = 0.32 , s = 0.25 , _{d 1} = 0.1 , and _{d 2} = 2 . Time: (a) t = 0 , (b) t = 300 , (c) t = 1000 , and (c) t = 3000 .

(a) [figure omitted; refer to PDF]

(b) [figure omitted; refer to PDF]

(c) [figure omitted; refer to PDF]

(d) [figure omitted; refer to PDF]

In Figure 2, we show the time process of stripes pattern formation of the prey u at t = 0,300,1000,3000 for b = 0.495 , c = 0.5 , _{k 1} = 0.3 , _{k 2} = 0.2 , m = 0.3 , s = 0.25 , _{d 1} = 0.1 , _{d 2} = 2 . In this case, starting with a homogeneous state _{E 3} = ( 0.142318 , 0.342318 ) (cf., Figure 2(a)), the random perturbations lead to the formation of stripes and spots (cf., Figure 2(b)), and the later random perturbations make these stripes decay and end with the time-independent regular spots pattern (cf., Figure 2(d)).

The formation of stripes patterns of prey u with parameters b = 0.495 , c = 0.5 , _{k 1} = 0.3 , _{k 2} = 0.2 , m = 0.3 , s = 0.25 , _{d 1} = 0.1 , and _{d 2} = 2 . Time: (a) t = 0 , (b) t = 300 , (c) t = 1000 , and (c) t = 3000 .

(a) [figure omitted; refer to PDF]

(b) [figure omitted; refer to PDF]

(c) [figure omitted; refer to PDF]

(d) [figure omitted; refer to PDF]

In Figure 3, we show the time process of holes pattern formation of the prey u at t = 0,300,1000,3000 for b = 0.4 , c = 0.5 , _{k 1} = 0.3 , _{k 2} = 0.2 , m = 0.24 , s = 0.25 , _{d 1} = 0.1 , and _{d 2} = 2 . In this case, starting with a homogeneous state _{E 3} = ( 0.2 , 0.4 ) (cf., Figure 3(a)), the random perturbations lead to the formation of independent regular holes pattern (cf., Figure 1(d)), which is isolated zones with low prey densities.

The formation of holes patterns of prey u with parameters b = 0.4 , c = 0.5 , _{k 1} = 0.3 , _{k 2} = 0.2 , m = 0.24 , s = 0.25 , _{d 1} = 0.1 , and _{d 2} = 2 . Time: (a) t = 0 , (b) t = 300 , (c) t = 1000 , and (c) t = 3000 .

(a) [figure omitted; refer to PDF]

(b) [figure omitted; refer to PDF]

(c) [figure omitted; refer to PDF]

(d) [figure omitted; refer to PDF]

Comparing Figure 1(d) with Figure 3(d), we find that they share similarities. Figure 1(d) consists of red (maximum density of u ) spots on a blue (minimum density of u ) background; that is, the preys are isolated zones with high population density. Figure 3(d) consists of blue (minimum density of u ) spots on a red (maximum density of u ) background; that is, the preys are isolated zones with low population density. For the sake of learning the pattern dynamics of model (12) further, we illustrate the three-dimensional patterns in the space ( u , x , y ) . From Figure 4, one can realize the relations of the patterns (e.g., Figures 4(a) and 4(c)) with their corresponding numerical solutions (e.g., Figures 4(b) and 4(d)). In fact, the patterns (e.g., Figures 1-3) are the projections in x y plane of the numerical solutions u ( x , y , t ) to the model (12). In the software MATLAB, imagesc(u) displays u ( x , y ; t ) as an image; each numerical solution of u ( x , y ; t ) corresponds to a rectangular area in the image; that is, the values of the numerical solutions u ( x , y ; t ) to the model (12) are indices into the current colormap that determine the color of each patch.

Pattern dynamics of model (12) in two- and three-dimensional space. (a) spots pattern; (b) spots pattern in three-dimensional space; (c) holes pattern; (d) holes pattern in three-dimensional space.

(a) [figure omitted; refer to PDF]

(b) [figure omitted; refer to PDF]

(c) [figure omitted; refer to PDF]

(d) [figure omitted; refer to PDF]

6. Concluding Remarks

In the current investigation, we propose and analyze the dynamics of a reaction-diffusive Leslie-Gower predator-prey model with the additive Allee effect on prey. We are mainly concerned with the coexistence of the predator and prey and focus on the case of weak Allee effect (i.e., m < b ). The value of this study is twofold. First, it shows the nonexistence and existence of the nonconstant positive steady-states, which guarantees the existence of Turing patterns. Second, it rigorously proves Turing instability by linear stability analysis and illustrates all three categories of Turing patterns close to the onset of Turing bifurcation via numerical simulations which indicates that the model dynamics exhibits complex pattern replication.

We summarize our findings as well as their related biological implications as follows.

(1) Theorems 2 and 3 provide us with a full picture on the dynamics of the model with weak Allee effect. The dynamics of the model introduced can be very complicated due to Allee effects. In Theorem 2(a), we find that the trivial constant solution _{E 0} = ( 0,0 ) of the model which is subject to an Allee effect is unstable, which is exactly consistent with the model without Allee effect [24-26]. This demonstrates that there is no extinction of a species in the present of Allee effects. By comparing them to their corresponding models without Allee effect, we can conclude that Allee effect can make the extinction of the prey although the maximum predation rate c is small (see Theorem 2(b)). In this case, the predator species intends to change its food habits as predator has sufficient resources for alternative foods [56]. Furthermore, from Theorem 3, one can obtain that _{E 3} = ( ^{u *} , ^{u *} + _{k 2} ) is locally uniformly asymptotically stable, which means that nonconstant positive solution (stationary pattern) of model (12) unlikely exists.

(2) Theorems 9 and 11 indicate the existence and nonexistence of nonconstant steady-states with respect to various parameters. Roughly speaking, we can state that there is no pattern if the diffusion coefficients are suitably chosen. While pattern occurs provided that ^{θ 1 1} - 4 _{θ 2} > 0 and all the eigenvalues _{μ i} have odd multiplicity (see Theorem 11). Hence, we can conclude that the multiplicity of patterns seems very interesting from the viewpoint of mathematics.

(3) Theorem 15 and numerical simulations give us the existence of conditions of Turing instability and the types of Turing pattern. From the numerical results, one can see that our model has rich and complex spatiotemporal behavior. We find three typical Turing patterns, that is, stripes pattern, spot-stripe mixtures pattern, and spots pattern. That is to say, the effect of the Allee effect for pattern formation is tremendous. Therefore, the results of the present paper and [24, 25] show that the types of Turing pattern in the biological models depend on the effect of the Allee effect. In other words, the Allee effect may be one of the determining factors in producing spots and spot-stripe mixtures Turing patterns.

It is believed that our results made in this investigation related to predator-prey interactions due to the effect of Allee effect would certainly be of some help to theoretical mathematicians and the ecologists who are engaged in performing experimental work. Further studies are necessary to analyze the behaviour of a reaction-diffusion predator-prey model with the strong Allee effect (i.e., m > b ).

Acknowledgments

The authors thank the anonymous referee and Professor Shangbin Cui for very helpful suggestions and comments which led to improvement of their original paper. This research was supported by the National Science Foundation of China (61373005, 11171357, and 11271290) and Zhejiang Provincial Natural Science Foundation (LY12A01014).

References

[1] R. M. May Stability and Complexity in Model Ecosystems , Princeton University Press, Princeton, NJ, USA, 1973.

[2] F. Albrecht, H. Gatzke, A. Haddad, N. Wax, "The dynamics of two interacting populations," Journal of Mathematical Analysis and Applications , vol. 46, pp. 658-670, 1974.

[3] R. M. May, R. M. May, "Models for two interacting populations," Theoretical Ecology: Principles and Applications , pp. 78-104, Blackwell Scientific Publications, London, UK, 1981.

[4] P. A. Abrams, L. R. Ginzburg, "The nature of predation: prey dependent, ratio dependent or neither?," Trends in Ecology & Evolution , vol. 15, no. 8, pp. 337-341, 2000.

[5] P. F. Verhulst, "Notice sur la loi que la population suit dans son accroissement," Correspondances Mathmatiques et Physiques , vol. 10, pp. 113-121, 1838.

[6] C. S. Holling, "The components of predation as revealed by a study of small mammal predation of the European Pine Sawfly," The Canadian Entomologist , vol. 91, no. 5, pp. 293-320, 1959.

[7] N. D. Kazarinoff, P. van den Driessche, "A model predator-prey system with functional response," Mathematical Biosciences , vol. 39, no. 1-2, pp. 125-134, 1978.

[8] M. L. Rosenzweig, "Paradox of enrichment: destabilization of exploitation ecosystems in ecological time," Science , vol. 171, no. 3969, pp. 385-387, 1971.

[9] J. Wang, J. Shi, J. Wei, "Predator-prey system with strong Allee effect in prey," Journal of Mathematical Biology , vol. 62, no. 3, pp. 291-331, 2011.

[10] J. Wang, J. Shi, 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.

[11] R. Peng, J. Shi, "Non-existence of non-constant positive steady states of two Holling type-II predator-prey systems: strong interaction case," Journal of Differential Equations , vol. 247, no. 3, pp. 866-886, 2009.

[12] P. H. Leslie, "Some further notes on the use of matrices in population mathematics," Biometrika , vol. 35, no. 3-4, pp. 213-245, 1948.

[13] P. H. Leslie, J. C. Gower, "The properties of a stochastic model for the predator-prey type of interaction between two species," Biometrika , vol. 47, pp. 219-234, 1960.

[14] E. C. Pielou An Introduction to Mathematical Ecology , Wiley-Interscience, New York, NY, USA, 1969.

[15] M. A. Aziz-Alaoui, "Study of a Leslie-Gower-type tritrophic population model," Chaos, Solitons & Fractals , vol. 14, no. 8, pp. 1275-1293, 2002.

[16] M. A. Aziz-Alaoui, M. Daher Okiye, "Boundedness and global stability for a predator-prey model with modified Leslie-Gower and Holling-type II schemes," Applied Mathematics Letters , vol. 16, no. 7, pp. 1069-1075, 2003.

[17] Y. Du, S.-B. Hsu, "A diffusive predator-prey model in heterogeneous environment," Journal of Differential Equations , vol. 203, no. 2, pp. 331-364, 2004.

[18] Y. Du, M. Wang, "Asymptotic behaviour of positive steady states to a predator-prey model," Proceedings of the Royal Society of Edinburgh A , vol. 136, no. 4, pp. 759-778, 2006.

[19] R. Peng, M. Wang, "Positive steady states of the Holling-Tanner prey-predator model with diffusion," Proceedings of the Royal Society of Edinburgh A , vol. 135, no. 1, pp. 149-164, 2005.

[20] R. Peng, M. Wang, "Global stability of the equilibrium of a diffusive Holling-Tanner prey-predator model," Applied Mathematics Letters , vol. 20, no. 6, pp. 664-670, 2007.

[21] H.-B. Shi, W.-T. Li, G. Lin, "Positive steady states of a diffusive predator-prey system with modified Holling-Tanner functional response," Nonlinear Analysis: Real World Applications , vol. 11, no. 5, pp. 3711-3721, 2010.

[22] X. Li, W. Jiang, J. Shi, "Hopf bifurcation and Turing instability in the reaction-diffusion Holling-Tanner predator-prey model," IMA Journal of Applied Mathematics , vol. 78, no. 2, pp. 287-306, 2013.

[23] P.-P. Liu, Y. Xue, "Spatiotemporal dynamics of a predator-prey model," Nonlinear Dynamics , vol. 69, no. 1-2, pp. 71-77, 2012.

[24] B. I. Camara, M. A. Aziz-Alaoui, "Dynamics of a predator-prey model with diffusion," Dynamics of Continuous, Discrete &; Impulsive Systems A , vol. 15, no. 6, pp. 897-906, 2008.

[25] B. I. Camara, M. A. Aziz-Alaoui, "Turing and Hopf patterns formation in a predator-prey model with Leslie-Gower-type functional response," Dynamics of Continuous, Discrete & Impulsive Systems B , vol. 16, no. 4, pp. 479-488, 2009.

[26] B. Camara Complexite de dynamiques de modéles proie-prédateur avec diffusion et applications [Ph.D. thesis] , Universite du Havré, Lorient, France, 2009.

[27] F. Courchamp, T. Clutton-Brock, B. Grenfell, "Inverse density dependence and the Allee effect," Trends in Ecology & Evolution , vol. 14, no. 10, pp. 405-410, 1999.

[28] P. A. Stephens, W. J. Sutherland, R. P. Freckleton, "What is the Allee effect?," Oikos , vol. 87, no. 1, pp. 185-190, 1999.

[29] P. A. Stephens, W. J. Sutherland, "Consequences of the Allee effect for behaviour, ecology and conservation," Trends in Ecology & Evolution , vol. 14, no. 10, pp. 401-405, 1999.

[30] F. Courchamp, L. Berec, J. Gascoigne Allee Effects in Ecology and Conservation , Oxford University Press, Oxford, UK, 2008.

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

[32] S.-R. Zhou, Y.-F. Liu, G. Wang, "The stability of predator-prey systems subject to the Allee effects," Theoretical Population Biology , vol. 67, no. 1, pp. 23-31, 2005.

[33] L. Berec, E. Angulo, F. Courchamp, "Multiple Allee effects and population management," Trends in Ecology & Evolution , vol. 22, no. 4, pp. 185-191, 2007.

[34] P. Aguirre, E. González-Olivares, E. Sáez, "Two limit cycles in a Leslie-Gower predator-prey model with additive Allee effect," Nonlinear Analysis: Real World Applications , vol. 10, no. 3, pp. 1401-1416, 2009.

[35] P. Aguirre, E. González-Olivares, E. Sáez, "Three limit cycles in a Leslie-Gower predator-prey model with additive Allee effect," SIAM Journal on Applied Mathematics , vol. 69, no. 5, pp. 1244-1262, 2009.

[36] Y. Cai, W. Wang, J. Wang, "Dynamics of a diffusive predatorprey model with additive Allee effect," International Journal of Biomathematics , vol. 5, no. 2, 2012.

[37] E. González-Olivares, J. Mena-Lorca, A. Rojas-Palma, J. D. Flores, "Dynamical complexities in the Leslie-Gower predator-prey model as consequences of the Allee effect on prey," Applied Mathematical Modelling , vol. 35, no. 1, pp. 366-381, 2011.

[38] M.-H. Wang, M. Kot, "Speeds of invasion in a model with strong or weak Allee effects," Mathematical Biosciences , vol. 171, no. 1, pp. 83-97, 2001.

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

[40] C. M. Taylor, A. Hastings, "Allee effects in biological invasions," Ecology Letters , vol. 8, no. 8, pp. 895-908, 2005.

[41] J. Shi, R. Shivaji, "Persistence in reaction diffusion models with weak Allee effect," Journal of Mathematical Biology , vol. 52, no. 6, pp. 807-829, 2006.

[42] G. A. K. van Voorn, L. Hemerik, M. P. Boer, B. W. Kooi, "Heteroclinic orbits indicate overexploitation in predator-prey systems with a strong Allee effect," Mathematical Biosciences , vol. 209, no. 2, pp. 451-469, 2007.

[43] B. Dennis, "Allee effects: population growth, critical density, and the chance of extinction," Natural Resource Modeling , vol. 3, no. 4, pp. 481-538, 1989.

[44] J. D. Murray, "Discussion: Turing's theory of morphogenesis - Its influence on modelling biological pattern and form," Bulletin of Mathematical Biology , vol. 52, no. 1-2, pp. 119-152, 1990.

[45] J. Chattopadhyay, P. K. Tapaswi, "Effect of cross-diffusion on pattern formation--a nonlinear analysis," Acta Applicandae Mathematicae , vol. 48, no. 1, pp. 1-12, 1997.

[46] D. Henry Geometric Theory of Semilinear Parabolic Equations , vol. 840, of Lecture Notes in Mathematics, Springer, Berlin, Germany, 1981.

[47] Y. Lou, W. Ni, "Diffusion, self-diffusion and cross-diffusion," Journal of Differential Equations , vol. 131, no. 1, pp. 79-131, 1996.

[48] C.-S. Lin, W.-M. Ni, I. Takagi, "Large amplitude stationary solutions to a chemotaxis system," Journal of Differential Equations , vol. 72, no. 1, pp. 1-27, 1988.

[49] J. Smoller Shock Waves and Reaction-Diffusion Equations , vol. 258, of Fundamental Principles of Mathematical Science, Springer, New York, NY, USA, 1983.

[50] L. Nirenberg Topics in Nonlinear Functional Analysis , Courant Institute, New York, NY, USA, 2001., 2nd.

[51] R. Peng, J. Shi, M. Wang, "On stationary patterns of a reaction-diffusion model with autocatalysis and saturation law," Nonlinearity , vol. 21, no. 7, pp. 1471-1488, 2008.

[52] M. Ghergu, V. Radulescu, "Turing patterns in general reaction-diffusion systems of Brusselator type," Communications in Contemporary Mathematics , vol. 12, no. 4, pp. 661-679, 2010.

[53] A. Turing, "The chemical basis of morphogenesis," Philosophical Transactions of the Royal Society of London B , vol. 237, pp. 37-72, 1952.

[54] 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.

[55] A. Munteanu, R. V. Solé, "Pattern formation in noisy self-replicating spots," International Journal of Bifurcation and Chaos , vol. 16, no. 12, pp. 3679-3685, 2006.

[56] L. N. Guin, M. Haque, P. K. Mandal, "The spatial patterns through diffusion-driven instability in a predator-prey model," Applied Mathematical Modelling , vol. 36, no. 5, pp. 1825-1841, 2012.