Dynamics of populations with individual variation in dispersal on bounded domains

ABSTRACT Most classical models for the movement of organisms assume that all individuals have the same patterns and rates of movement (for example, diffusion with a fixed diffusion coefficient) but there is empirical evidence that movement rates and patterns may vary among different individuals. A simple way to capture variation in dispersal that has been suggested in the ecological literature is to allow individuals to switch between two distinct dispersal modes. We study models for populations whose members can switch between two different nonzero rates of diffusion and whose local population dynamics are subject to density dependence of logistic type. The resulting models are reaction–diffusion systems that can be cooperative at some population densities and competitive at others. We assume that the focal population inhabits a bounded region and study how its overall dynamics depend on the parameters describing switching rates and local population dynamics. (Traveling waves and spread rates have been studied for similar models in the context of biological invasions.) The analytic methods include ideas and results from reaction–diffusion theory, semi-dynamical systems, and bifurcation/continuation theory.


Introduction
The movement of organisms plays an important role in determining the spatial distributions and interactions of populations, and thus influences many ecological processes. Most classical models for the movement of organisms assume that all individuals in a given population have the same patterns and rates of movement (for example, diffusion with a fixed diffusion coefficient) but there is empirical evidence that movement rates and patterns may vary among different individuals, or for the same individual under different conditions; see [5,6,13,19,35,36,40]. Furthermore, phenotypic variation in movement patterns can have important effects on population dynamics and the spatial structure of populations; see [10].
A simple way to capture variation in dispersal patterns is to allow individuals to switch between two dispersal modes. There has been some modelling of populations that can switch between different dispersal modes. The case where organisms can switch between moving and non-moving (quiescent) states has been studied a considerable amount, especially in the context of organisms living in advective environments. See, for example [17,24] and the references in those papers. The papers [37,41] focus on how switching between discrete movement modes influences the distribution of populations on the time scale of foraging or other movement, but do not connect the dispersal process directly to population dynamics. There has also been work showing that if individuals move by ordinary diffusion but with rates that are drawn from a continuum of possible rates according to a suitable probability distribution, the resulting movement pattern at the population level can mimic what would be expected from anomalous diffusion, specifically Lévy flights; see [18,32,33].
In a related but distinct line of research, researchers have studied systems arising as models for the evolution of dispersal, where mutations allow populations to change their dispersal rates. This idea was introduced in [11] in the case of a discrete set of dispersal rates, and has been explored by several researchers for a continuous set of dispersal rates; see [26] and the references in that paper. Another related but distinct source of interest in populations with multiple dispersal modes is the observation that the presence of dispersal polymorphism can affect the spread rates of biological invasions and that dispersal traits may evolve during invasions. Some empirical results on that topic are described in [34]. Theoretical results on travelling waves, spread rates and related topics for systems with multiple dispersal modes in the context of biological invasions are given in [7,12,15,16,31]. An extensive set of references to related work is given in [15]. Finally, the question of existence of positive equilibria for models with a discrete set of dispersal modes on bounded domains was addressed in [20].
In this paper, we will consider the system in a smooth bounded domain ⊂ R n (n ≥ 1) with classical (Dirichlet, Robin or Neumann) boundary conditions. Here all parameters are positive. We will describe in some detail how the equilibria and dynamics of the model depend on the parameters. Our results complement those of Girardin and Hei and Wu [15,20]. The results of Girardin [15] show that for a class of models including system (1) but allowing arbitrarily many possible movement modes, if the zero equilibrium is unstable then there exist travelling waves connecting the zero equilibrium to some positive equilibrium. They also include formulas for spreading speeds. However, the precise structure of the set of positive equilibria for the models is not studied in detail in [15], and the models are set on R rather than on bounded domains in R n . Some more detailed results on travelling waves for certain specific cases of models with two possible dispersal modes, again set on R, are derived in [12,31]. The results of Hei and Wu [20] give sufficient conditions for the existence of a positive equilibrium and some results on stability of small equilibria for a rather general class of models along the lines of system (1) on bounded domains. The class of models for which the results of Hei and Wu [20] hold contains some cases of system (1), and allows arbitrarily many movement modes and general competitive nonlinearities, but the hypotheses of the main existence theorem of Hei and Wu [20] impose significant restrictions on the parameters when applied to system (1). Roughly speaking, the conditions in [20] require α and β to be larger than certain expressions involving combinations of the other parameters in system (1) and the quantities ±(α − β). The paper is structured as follows: In Section 2, we show that the system is dissipative, examine the stability of the equilibrium (0, 0), calculate the minimal patch size needed for the instability of (0, 0) under Dirichlet conditions, and derive some other general properties of the system. In Section 3 we study in detail the equilibria and dynamics of the system of ordinary differential equations arising from setting the diffusion rates equal to zero in system (1). In Section 4 we study the dynamics of the system with Neumann boundary conditions, partly by using the results from the previous sections. We give a summary of the results and discuss their biological significance in Section 5.

Dissipativity and instability of (0, 0)
The system (1) can be written as follows: where We will also assume that u and v satisfy homogeneous Neumann, Dirichlet, or Robin boundary conditions. The local existence of classical solutions follows from standard results, see, for example, the discussion and references in [9, Sections 1.6.5 and 1.6.6], or [3,21,22]. Global existence follows if solutions are bounded by some finite B(T) in [L ∞ ( )] 2 on any time interval (0, T) with T > 0; see, for example [1,2]. Using the so-called method of contracting rectangles (see, e.g. [39,Chapter 14.E]), we have the following result on the uniform boundedness of solutions of system (1).
Thus, any solution of system (2) with nonnegative bounded initial data exists for all t ≥ 0, and eventually lies in the rectangle

Proof:
For system (2) with all fixed positive parameters, it is easy to see that g 1 (0, v) = βv > 0 and g 2 (u, 0) = αu > 0 for any u,v > 0. Moreover, there exist M 1 and N 1 such that a − α − bM 1 < 0, β − cM 1 < 0, d − β − fN 1 < 0 and α − eN 1 < 0. Then there exists some small > 0, such that for any . Also, the system obtained by replacing g 1 and g 2 in (2) with g + 1 and g + 2 is cooperative, and solutions to the original system are sub-solutions to the modified system. Let (U(t), V(t)) be the solution of ODE System This will be either a solution (in the Neumann case) or a super-solution (in the Dirichlet or Robin case) for the system obtained by replacing g 1 , g 2 with g + 1 , g + 2 in (2). By the comparison theorem (see, e.g. [38,Chapter 7 Linearizing system (1) at (0, 0), we have where B denotes a Dirichlet, Neumann, or Robin boundary operator. The system (4) is cooperative and irreducible, so the operator on the right-hand side of the partial differential equations has a compact positive resolvent by standard elliptic theory. By the celebrated Krein-Rutman theorem, it follows that the eigenvalue problem admits a principal eigenvalue λ 0 with a positive eigenfunction ψ = (ψ 1 , ψ 2 ).
and denote the spectral bound of A as s(A) := max{Reλ} where λ is any eigenvalue of A. Note that the principal eigenvalue is λ 0 = s(A) for the Neumann boundary conditions because in that case the principal eigenvalue of the Laplacian is 0 and the eigenfunction is constant. We will show that s(A) > 0. Clearly, tr Consequently, A always has two real eigenvalues and s(A) = max{λ}.
For the Dirichlet or Robin boundary condition, the sign of λ 0 is not obvious. Suppose we have Dirichlet boundary conditions. (The case of Robin conditions is similar.) Let (λ d , φ d ) be the principal eigenvalue and associated positive eigenfunction of u = λu, x ∈ , u = 0, x ∈ ∂ . Suppose d 1 λ d + a > 0 and d 2 λ d + d > 0. This guarantees that when α = β = 0, system (2) admits two semi-trivial steady states, and unstable trivial steady state (0, 0). When α and β are small, we still can argue that λ 0 > 0. But what will happen in other cases for (α, β)? Proposition 2.3: If a + d 1 λ d > 0 and d + d 2 λ d > 0, then the trivial steady state (0, 0) is unstable for system (2) under the Dirichlet condition.
Proof: Suppose φ 2 = kφ 1 for some k > 0. Then the eigenvalue problem (5) reduces to and turns out that each of the equations of (6) admits a principal eigenvalue, denoted by Indeed, let k * be the unique positive root of It follows that Therefore, under the Dirichlet condition system (4) admits a principal eigenvalue λ 0 > 0 with a positive eigenfunction (φ d (·), k * φ d (·)).

Remark 2.4:
Indeed, Proposition 2.3 also works for the Robin condition. Moreover, if one of α and β is zero, say, α = 0, then λ = a and λ = a + d 1 λ d are positive eigenvalues for eigenvalue problem (5) under the Neumann and Dirichlet (or Robin) condition, respectively. Therefore, (0, 0) is still unstable.

Minimal patch size under Dirichlet boundary conditions
The sufficient condition for instability of (0, 0) in Proposition 2.3 is simple but not sharp. The eigenvalue λ d depends on , and if λ d is sufficiently negative the equilibrium (0, 0) will be stable. This leads to the phenomenon that under Dirichlet conditions there will be a minimal patch size needed for population growth. In order to study the minimal patch size needed to support a population we account for the size by writing = ˜ 0 with |˜ 0 | = 1 and then rescaling the model on back to˜ 0 (see [ (7) in the proof of Proposition 2.3 in terms of Q 1 and Q 2 , solving, and substituting into the formula for λ 0 , we have which is the larger root of We have In the case max{a − α, d − β} < 0 we have 0 < α − a < α and 0 < β − d < β so the constant term in inequality (8) is negative, so the quadratic must have one positive and one negative root. Denote Then − (α, β) , the smaller root of Equation (9), must be negative. Since λ d < 0 satisfies inequality (8), we obtain the condition − (α, β) ≤ λ d . Therefore, in this case λ 0 ≥ 0 is equivalent to λ d ≥ − (α, β), that is, 2 ≥ −λ 0 / (α, β). This then gives us the minimal patch size n where n is the number of space dimensions and = −λ 0 / (α, β).

Now direct computation indicates
.

Dynamics under Neumann boundary conditions
Define Then L 1 and L 2 divide the first quadrant of α, β plane into at most four regions (see Figure 1). Set We are ready to state an important observation for system (1) with the Neumann condition.
) is the solution of system (1) with Neumann boundary conditions and any bounded nonnegative initial data = (φ 1 , φ 2 ) with ≡ 0. Then the following statements are valid.
Proof: For (a), in the system written in the form (2), we have Adapting the proof in Proposition 2.1, we obtain that (a) is valid.
Employing the same comparison argument as in Proposition 2.1, we can show that lim inf t→∞ u(t, ·, ) > β/c and lim inf t→∞ v(t, ·, ) > α/e. Remark 2.6: Proposition 2.5 (a) is also valid for system (1) under Dirichlet boundary conditions provided a + d 1 λ d > 0 and d + d 2 λ d > 0, and similarly for Robin conditions, based on analogous arguments in those cases.

Existence of positive equilibria with one component small
When α = β = 0, the system (1) becomes a competitive Lotka-Volterra model. By standard results for diffusive logistic equations, as discussed for example in [9], Chapter 3, the system (1) will have semi-trivial equilibria (u 0 , 0) and (0, v 0 ) in the boundary of the positive cone under Neumann boundary conditions if a > 0 and d > 0, and under Dirichlet boundary conditions if a + d 1 λ d > 0 and d + d 2 λ d > 0. We next look at the behaviour of boundary equilibria for α = β = 0 as those parameters become positive. For convenience, we only consider the Dirichlet boundary case. Other boundary conditions can be treated in a similar manner.
) when α and β are sufficiently small.

Proof:
Here, we only consider the case, where a + d 1 λ d > 0, the other case can be treated in a similar manner. Let F : where α and β are smooth, α(0) = β(0) = 0 and α (0) > 0, β (0) > 0. Clearly, the mapping F is continuous and differentiable. Indeed, we have It then follows that F(u 0 , 0, 0) = 0 0 and If 0 is not an eigenvalue of then has a unique solution for any q ∈ C α (¯ , R). Also, since u 0 > 0 in , the principal eigenvalue for is τ = 0, so the principal eigenvalue of will be negative. Therefore, if ψ is determined uniquely by Equation (13) then the equation has a unique solution for any p ∈ C α (¯ , R). This implies D (u,v) F(u 0 , 0, 0) is invertible, and hence, it follows from the implicit function theorem that in some neighbourhood of When (u 0 , 0) is linearly stable for system (1), the principal eigenvalue for Equation (12) is negative, so the principal eigenvalue and hence all other eigenvalues of are positive, so the operator L has a positive resolvent (in other words, Equation (17) has a maximum principle) which implies that v > 0 in . Thus, increasing z slightly from z = 0 will produce a positive steady state (u * 0 , v * 0 ) with both components positive in . In other words, if (α, β) are small enough, there is a positive steady state close to (u 0 , 0). The argument in the case of (0, v 0 ) is analogous to that for (u 0 , 0).

Spatially homogeneous case
To better understand the effect of the switching rates α and β on the dynamics of the system (1), we first study the ODE system In view of Proposition 2.2, it follows that (0, 0) is unstable. Next, we discuss the existence and stability of positive equilibria.

Basic properties of equilibria
Let (u * , v * ) be a positive equilibrium of system (19). Then we have Set k = v * /u * . Substitute k into system (20) and simplify the equations. We get Thus, it suffices to find the positive values of k satisfying the equation Simplifying the above formula, we get the cubic equation: Since P(0) = −αb < 0 and P(∞) = ∞, we immediately obtain that there exists at least one positive root k * for P(k). Therefore, system (19) admits (at least) one positive equilibrium (u * , v * ), and at most three positive equilibria. Proposition 3.1: Suppose (u * , v * ) is any positive equilibrium for system (19).
is hyperbolic, then it is either a saddle or a stable node. Moreover, it is always a stable node provided bf ≥ ce.
is hyperbolic and could be a stable node or stable spiral.
Proof: Consider the Jacobian matrix at (u * , v * ), namely In the case where L 1 (α, β)L 2 (α, β) > 0, it follows from Proposition 2.5 that a 12 a 21 > 0, and hence, (tr(J)) 2 It then follows that det J = a 11 a 22 − a 12 a 21 ≥ a 11 a 22 > 0, and hence, the equilibrium is hyperbolic. However, it is not obvious that what type of equilibrium it is. Indeed, numerical computation indicates that (u * , v * ) could be a stable node or stable spiral.

Conditions for global stability
Next, we use a Lyapunov function to prove the global stability of the positive equilibrium (u * , v * ) for α and β in certain regions of the αβ-plane. This analysis will also extend to the PDE case with Neumann boundary conditions, which we will discuss later. Let The inequality above is valid since if u(t 0 ) > u * > 0 for some Similarly, we also have Note that if AD > BC and A,D are positive, then there exists K > 0 such that Ax 2 + Bxy + K(Cxy + Dy 2 ) > 0, for every x ≥ 0, y ≥ 0 except (0, 0). Clearly, stays positive for every x ≥ 0, y ≥ 0 except (0, 0) provided we can choose K > 0 so that Such is the case if the discriminant for −B 2 /4A In region II a , II b , it is easy to see that BC ≤ 0 < AD, and hence, there exists K 1 > 0 such that The equality holds true if and only if (u, v) ≡ (u * , v * ).
In Region III, B and C are positive, we need AD > BC, that is, bf The equality holds true if and only if (u, v) ≡ (u * , v * ).

Remark 3.3:
If (α, β) ∈ III and bf ≥ ce, e.g. in the weak competition case (b/e > a/d > c/f , a unique stable coexistence), the condition in part (b) is automatically valid because

Bifurcations
Next, we will explore some bifurcations for the system (19). Consider the strong competition case, that is, bf < ce. First, we give the number of positive equilibria when α and β are small enough.

Proof:
The proof is based on using the cubic discriminant and checking the signs of coefficients in Equation (22). Recall that the discriminant of the cubic , and that when it is positive the cubic has three distinct real roots. For the polynomial P(k) all terms in the discriminant except the one corresponding to B 2 C 2 go to zero as α and β go to zero, and when a/d = c/f and a/d = b/e the limit of the remaining term is positive, so when α, β are small enough P(k) always has three distinct real roots k i (i = 1, 2, 3). Now when c/f > a/d > b/e and α = β = 0, system (19) admits two linearly stable boundary equilibria. Adapting the proof in Lemma 2.7, together with the fact k 1 k 2 k 3 = αb/βf > 0, we see that system (19) admits three positive equilibria. When α, β are small enough and a/d > c/f > b/e (or c/f > b/e > a/d) holds, we have k 1 + k 2 + k 3 = −(af − cd − f α + (c + e)β)/βf < 0 (or k 1 k 2 + k 1 k 3 + k 2 k 3 = (ae − bd − (c + e)α + bβ)/βf < 0). Since k 1 k 2 k 3 > 0 is positive, we see system (19) has a unique positive equilibrium. Now we are ready to introduce the main result on the bifurcation of system (19).
Theorem 3.5: Suppose bf < ce, then a (codimension-two) cusp bifurcation occurs at the bifurcation curve in the αβ-plane with parametric form where m = af−cd and n = ae − bd, and k satisfies one of the following cases: Moreover, any non-hyperbolic equilibrium only occurs when (α, β) lies on the bifurcation curve (24) provided one of (A1)-(A3) holds.

The cubic equation is βfk
Rewrite the equation as follows: Note that for any k ≥ 0, fk 2 + (e + c)k + b > 0, thus we have The left-hand-side of Equation (25) represents a curve which is fixed, independent of α and β, and the right-hand side gives a line with two parameters. Moreover, Since m < 0 and n > 0, it follows that fn − m(e + c) > 0. Let k * 1 be the positive root of For k ∈ (0, k * 1 ), V (k) > 0, and for k ∈ (k * 1 , ∞), V (k) < 0, and V(0) = 0, V(∞) = m/f < 0. This will help to plot the curve V(k). Note that as α and β vary, the line on the righthand side of Equation (25) moves but the curve does not. This suggests that analysing the intersections of the curve and the line graphically should be possible. The essential arguments used to find the number of intersections of two curves are similar to those in the Spruce Budworm model [28].
For example, fix a sufficiently small α, let β increase (or fix a relatively large α, let β increase), then the line rotates clockwise around (0, α). Now we have the following pictures.
We need to find the bifurcation turning point. This requires us to check the tangency condition. Taking the derivative (with respect to k) for both sides of Equation (25), we have So for any k > k * 1 , we have β > 0. On the other hand, we have Note that mf < 0 and n(e + c) − mb > 0. Let k * 2 be the positive root of equation Then k * 2 > −n/m > 0. Actually, we have k In a word, we have where k ∈ [k * 1 , k * 2 ]. Checking the derivatives of α(k) and β(k) with respect to k, we have α (k) = −kV (k), Since n > 0 and m < 0, it is easy to see that there exists a unique k c > 0 such that V (k) = 0. Therefore, (α(k c ), β(k c )) is the cusp point with α (k c ) = β (k c ) = 0. Moreover, α(k) and β(k) are increasing on (k * 1 , k c ) and decreasing on (k c , k * 2 ). See Figure 3. In a word, a typical cusp bifurcation occurs.  In the competitive exclusion case (A2) or (A3), see Figure 4 for illustration. We have (A2) m ≥ 0, n > 0 and ce > bf. Now, mimicking our previous analysis, we have Since m ≥ 0, n > 0, a necessary condition for β ≥ 0 in Equation (32) is fn − m(e + c) > 0. If fn − m(e + c) ≤ 0, then V (k) > 0 for k ∈ R. It then follows that Equation (25) has only one root, and hence, system (19) admits a unique positive equilibrium. Since n/m > (e + c)/f > b/(e + c), it follows that mfk 2 + 2nfk + n(e + c) − mb > 0 for k ≥ 0. Now if k > k * 1 , we have β > 0 and α > 0 in Equation (32). (A3) m < 0, n ≤ 0 and ce > bf. The tangency condition will also give the same parametric form as Equation (32). (32) for any k > 0. If n(e + c) − mb ≤ 0, then for any k > 0, α < 0 in Equation (32). Thus it is necessary that n(e + c) − mb > 0; that is, ce − bf > −ne/a ≥ 0. This leads to α > 0 in Equation Since k = v * /u * and Equation (21) holds, it follows from Equation (33) that and hence, It is easy to see that det(M) = (ce − bf )[fk 2 + (c + e)k + b] 2 > 0. Thus, we can directly solve (α, β) in terms of k. A tedious computation also gives the parametric form Based on the previous analysis, it follows that a non-hyperbolic equilibrium can only occur on the bifurcation curve.
By Theorem 3.5 and Proposition 3.4, we have the following result.

Global dynamics
We can now make some general statements about the dynamics of (19). (19) has no periodic orbit in the first quadrant.
Similarly, we can also show that the first part of (b) is valid by Poincaré-Bendixson theorem, and the unique equilibrium is a stable node, too. Remark 3.9: Indeed, we can get a similar result to that in Theorem 3.8, showing that nonzero boundary equilibrium (0, v * ) or (u * , 0) is globally attractive when α > 0, β = 0 or α = 0, β > 0, respectively.

Dynamics and diffusion
In this section, we mainly study the system (1) with Neumann boundary conditions, in the case where (α, β) is quantitatively large, that is, (α, β) ∈ I, see Figure 1. Some of the results are also valid for other boundary conditions.

Spatially constant solutions
Let X = C(¯ , R 2 ) and X + = C(¯ , R 2 + ). Then (X, X + ) is a strongly ordered Banach space. The system (1) with Neumann boundary conditions generates a semiflow on X, see [30]. Note that any solution of system (19) is a solution of system (1) with Neumann boundary conditions. If (α, β) ∈ I or III, Proposition 2.5 shows that system (1) can be simply treated as a competitive or cooperative system. Then we can use comparison arguments as in [38,42] to get the global dynamics of system (1). In particular, if the system (19) has a unique globally attracting equilibrium, then comparison principles in the appropriate ordering imply that it is also globally attracting and therefore the unique equilibrium for the reaction-diffusion model (1). Note that the Lyapunov function in Lemma 3.2 is convex, so it also works for the reaction-diffusion model, for example as in [8]. Thus, in view of Lemma 3.2 and Theorem 3.8, we are ready to state the main result on the global dynamics of system (1).

Theorem 4.1: The following statements hold for the reaction-diffusion model (1) under Neumann boundary conditions:
(a) If (α, β) ∈ I or II, then system (1) admits a unique positive constant steady state (u * , v * ), which is globally attractive for any initial data ∈ X + \ {0}. These results follow immediately from Theorems 3.5 and 4.1. Suppose U(t, ·, ) is a solution of system (1) with homogeneous Neumann boundary conditions and U(0, ·, ) = ∈ X + . Then the following statements hold.
Next, we discuss the existence and nonexistence of a non-constant steady state when bf < ce and one of (A1)-(A3) is valid. Clearly, if one of (A1)-(A3) holds, the bifurcation curve should lie in region III of the αβ-plane.

Nonconstant solutions
In the sequel, we just consider the solution of system (1) in the rectangle (1) is a competition-diffusion system in B, many conclusions follow from what we know about competition-diffusion systems in general.
Note that the u-nullcline is defined by If 1 and 2 have three intersections, system (1) admits three constant steady states, . Suppose that u * 1 < u * 2 < u * 3 , then the expression of 1 and 2 indicates that v * 1 > v * 2 > v * 3 . By Proposition 3.1 and Poincaré-Bendixson theorem, we know it must be the case that there are two stable nodes and one saddle. Suppose that (u * i , v * i ), i = 1, 2 are stable nodes. Clearly, there is no further equilibrium point ] for the ODE system (19). It follows immediately from [22, Proposition 9.1] that there is an entire orbit connecting these two stable nodes, a contradiction. Thus, (u * 2 , v * 2 ) has to be a saddle. Now we have the following estimate for the non-constant positive steady state of system (1).
Likewise, we have Note that Then we see the only possibility for a non-constant solution (u(x), v(x)) is This establishes the proposition.
Consider the elliptic problem We will now give conditions for the nonexistence of non-constant solutions of system (40) under various conditions on α and β related to the number of equilibria for system (19), as discussed in Section 3.3, specifically Theorem 3.5 and Corollary 3.6, and illustrated in Figures 3 and 4.
constant solution (u * , v * ), we always have (a − bu * − cv * )u * + (d − eu * − f v * )v * = 0, so that one of u * , v * must be greater than a positive constant independent of (α, β). Then having one of the diffusion rates large (for u * bounded below, having d 2 large, and for v * bounded below, having d 1 large) excludes the existence of non-constant solutions by the arguments used in the proof of Lemma 4.4.
The following result on the existence and stability of non-constant positive equilibria for (40) is based on well-known results in the theory of competitive reaction-diffusion systems given in [29, Theorem A] and [25]. A typical example of the type of domains in two space dimensions that may admit stable non-constant solutions is a 'dumbbell'-shaped region consisting of two large roughly circular regions connected by a sufficiently narrow strip.

General boundary conditions
Throughout this subsection, we will assume that (α, β) ∈ I, that is, the switching rates α and β are high enough that the system is in effect asymptotically cooperative.
Consider the system (1) with boundary conditions where ⊂ R n (n ≥ 1) is a bounded domain, and if n > 1, suppose that ∂ is of class C 2+θ (0 < θ ≤ 1). Furthermore, assume that either B i w = w or B i w = ∂w/∂n + m i (x)w for some nonnegative function m i ∈ C 1+θ (∂ , R), where ∂/∂n denotes differentiation in the direction of outward normal n to ∂ . Let X = L p ( ) with p ∈ (n, ∞). For each γ ∈ ( 1 2 + n/2p, 1], let X γ i be the fractional power space of L p ( ) with respect to (−d i , B i ), that is, the fractional power space associated with the operator −d i with boundary conditions given by B i (see, e.g. [21,22]). 2 with continuous inclusion for ν ∈ [0, 2γ − 1 − n/p) (see, e.g. [21]). For the case of Neumann or Robin boundary conditions define the positive cone of X + γ to be the set of all functions in X γ that are nonnegative on , and for Dirichlet boundary conditions define the positive cone to be the set of functions that are nonnegative on and whose outward normal derivatives are nonnegative on ∂ . The positive cone on X + γ then has nonempty interior Int(X + γ ). Let · γ be the norm on X γ . It then follows that there exists a constant C γ > 0 such that for all = (φ 1 , φ 2 ) ∈ X γ the inequality ∞ := max max x∈¯ |φ 1 (x)|, max x∈¯ |φ 2 (x)| ≤ C γ γ holds. It is well known that reaction-diffusion systems with smooth coefficients on smooth bounded domains generate semidynamical systems on fractional power spaces such as X γ ; see, for example, [21,38,42]. Indeed, they generate flows on subspaces of [C k (¯ ))] 2 that incorporate the boundary conditions; see [30]. By Proposition 2.1 solutions of system (1) are globally bounded and all of them eventually take values in a particular bounded set. Combining that with standard parabolic regularity implies that the system is dissipative. In fact, since we assume (α, β) ∈ I, it follows from Proposition 2.5 and Remark 2.6 that the set A := [0, β/c] × [0, α/e] is positively invariant and the values of any positive solution lie in A for sufficiently large t and system (1) is cooperative for (u, v) in that invariant set. It follows that solutions of any equilibria of the system (1) must take values in A, and that system (1) generates a semiflow on the subset of X + γ whose elements take values in A. Furthermore, because the system is cooperative, standard parabolic comparison principles can be used to show that it is strongly monotone, and parabolic regularity theory implies that forward semiorbits are precompact. ( (1) and U * = (β/c, α/e). Denote X U * := { ∈ X γ : 0 ≤ (x) ≤ U * }(so that X U * is the subset of X + γ whose elements take values in A. ) We have the following: Theorem 4.7: Suppose that (α, β) ∈ I and that (0, 0) is linearly unstable. For any ϕ ∈ X U * , let U(t, x, ϕ) be the solution of system (1) with boundary conditions (54) and U(0, x, ϕ) = ϕ(x) for x ∈ . Then system (1) has a unique positive steady stateˆ (x) ∈ X U * such that lim t→∞ U(t, ·, ϕ) −ˆ (·) γ = 0 for any ϕ ∈ X U * \{0}. Proof: (sketch) The dissipativity and precompactness of forward orbits for system (1) with boundary conditions 54 imply that the semidynamical system defined by U(t, x, ϕ) has a compact global attractor; see [4]. If (0, 0) is linearly unstable, the principal eigenvalue of the linearization of system (1) at (0, 0) is positive. Let = (ψ 1 , ψ 2 ) be the associated positive eigenfunction. It is easy to see that is a strict subsolution of (1) for > 0 sufficiently small. Also, by the strong comparison principle, U(t, x, ϕ) will lie in the interior of the positive cone for t > 0. (This is where we use the nonpositivity of the normal derivative on ∂ in the definition of the positive cone in the Dirichlet case). It follows that the semidynamical system defined by U(t, x, ϕ) has a nonempty attractor, which is global relative to orbits with nonzero initial data, and which is contained in the order interval defined by < (u, v) < U * . Parabolic regularity implies that the image of X U * under the semiflow is compact for any t > 0. Also, the dynamical terms of system (1) are subhomogeneous, so by Theorem 3.1 of Hirsch [23] the semiflow U(t, x, ϕ) is, also. It then follows from Theorem 5.5 of Hirsch [23] that the attractor in the interior of the positive cone consists of a single fixed point, which proves the desired result. (Related arguments are discussed in [42], Theorems 1.3.6 and 2.3.2.)

Discussion
We have analysed a class of models for the dynamics of populations that can switch between two dispersal modes, specifically between two different rates of diffusion, in the context of populations inhabiting a fixed bounded region as opposed to that of biological invasions, as treated in [7,12,15,16,31]. The idea behind the models is to divide a population into subpopulations dispersing by diffusion but at different rates, and then allow individuals to select their dispersal rate either by behavioural switching (which would typically be reflected by high rates of switching relative to the timescale of population dynamics), or by phenotypic plasticity (where switching rates would be on a timescale comparable to population dynamics) or perhaps by evolution (where 'switching' would reflect the population-level effects of mutation and selection), which might occur at a slower timescale than population dynamics, but might also occur on a timescale as fast as that of dispersal (see, for example, [7,12,34].) Our results show that such models, as formulated in model (1), can behave essentially as either cooperative systems (when the switching rates are high compared to the rates of population dynamical processes) or as competitive systems (when they are low), or display both behaviours depending on the densities of the two subpopulations.
In all cases, the model (1) will predict persistence and have at least one positive equilibrium if (0, 0) is unstable. In the cooperative case, the model (1) typically has either a unique positive equilibrium that is globally stable among positive solutions (if (0, 0) is unstable) or no positive equilibrium, with (0, 0) globally stable. Thus, the model behaves much like a single logistic equation. In the competitive case, the dynamics may be more complicated. Specifically, there may be multiple positive equilibria. However, even in cases where the dynamics are of mixed cooperative and competitive types, the corresponding models without diffusion such as system (19) cannot have periodic orbits. In the case of Neumann boundary conditions, the equilibrium (0, 0) is always unstable if the growth rates at low density of the subpopulations are positive, so the model predicts persistence. For Dirichlet or Robin conditions, the stability of (0, 0) depends on the size of the principal eigenvalue of the Laplacian on the underlying spatial domain, relative to the growth rate for the system at low densities. In the case of Dirichlet boundary conditions, we use this observation to analyze the minimum patch size needed to support a population.
In the competitive case, the models can have dynamics analogous to the cases of coexistence, competitive exclusion, or bistability (founder effects) in ordinary Lotka-Volterra competition models. The effective competition arises in the models because we assume that both subpopulations are subject to crowding effects of the sort that lead to logistic models or Lotka-Volterra competition models. In the case of strong competition, the model in the competitive case may have up to three distinct positive equilibria on general spatial domains if the switching rates are low. A typical configuration in that situation would be where there are stable equilibria with one component small. (These arise if a small amount of switching is introduced to a system with competition but no switching that has stable equilibria with one component zero.) In the case of Neumann boundary conditions, if there are two distinct stable spatially constant equilibria, for non-convex domains in two or more dimensions the model can have stable non-constant solutions which are close to one of the spatially constant equilibria on some subdomains and close to the second on others. This could conceivably suggest a mechanism for allopatric speciation based on dispersal rate, for example, if the movement rate is associated with competing hunting strategies of ambush or pursuit (as described in [14], for example). In the cooperative case, the model (1) typically has either a unique positive equilibrium that is globally stable among positive solutions (if (0, 0) is unstable) or no positive equilibrium, with (0, 0) globally stable. Thus, in that case, the model behaves much like a single logistic equation.
In the present study, we restrict our attention to the case where the coefficients describing population dynamics and switching rates are constant. For the case of Neumann boundary conditions, we show that for some ranges of the switching and population dynamical parameters, the dynamics of the system with diffusion are exactly the same as those of the corresponding model without diffusion. We obtain that result by using a Lyapunov functional. The dynamic behaviours of the nonspatial models are always a subset of those for the spatial models with Neumann boundary conditions in the case of constant coefficients, but in certain situations (for example, the bistable case in non-convex domains, as in Lemma 4.6), the dynamics of the spatial model may be more complex.
It is natural to think that organisms might switch movement modes in response to environmental conditions, or have different modes for searching for resources and exploiting them, and indeed this idea is developed in [13,41]. In a different direction, it is known that for populations that move by diffusion at a fixed rate, spatial heterogeneity favours slower diffusion (see [11]) In future work, we plan to explore these sorts of issues in models analogous to Equation (1) in spatially heterogeneous environments. One specific question we plan to study is competition between a population which can switch dispersal rates and another with a fixed dispersal rate, from the viewpoint of Dockery et al. [11].

Disclosure statement
No potential conflict of interest was reported by the authors.

Funding
Research partially supported by National Science Foundation (grant no. DMS-1514752).