Dynamic analysis of a fractional-order single-species model with diffusion ∗

Abstract. In this paper, we consider a fractional-order single-species model, which is composed of several patches connected by diffusion. We first prove the existence, uniqueness, non-negativity, and boundedness of solutions for the model, as desired in any population dynamics. Moreover, we also obtain some sufficient conditions ensuring the existence and uniform asymptotic stability of the positive equilibrium point for the investigated system. Finally, numerical simulations are presented to demonstrate the validity and feasibility of the theoretical results.


Introduction
In the natural world, dispersal often occurs among patches in ecological environments.In recent years, due to the ecological effects of local construction and development, such as the locations of manufacturing industries, the construction of dams and highways, as well as the development of tourism, habitats have been increasingly separated into isolated patches.Thus, realistic models should include dispersal processes that take into consideration the effects of spatial heterogeneity [11,16,18,19,[29][30][31]33].In [30], Takeuchi considered the following integer-order single-species system: where x i denotes population density in the ith patch, g i (x i ) denotes the density dependent growth rate in the ith patch, constant D ij 0 is the dispersal rate from patch j to patch i, and D ii = 0.By applying cooperative system theory, Takeuchi showed that species can survive in all patches at a globally asymptotically stable equilibrium point for any dispersal rate whenever it can survive in each isolated patch.In [16], Li and Shuai proved the global asymptotic stability of the positive equilibrium point for system (1) by using the graph theoretic approach.Many authors pointed out that fractional-order calculus is very suitable for the description of memory and hereditary properties of various materials and processes, which are neglected in classical integer-order models [9,10,13,26].Although many interesting and important results have been done in modelling the dynamics of population models, they have been restricted to integer-order differential equations.Recently, fractional-order differential equations have garnered a lot of attention and appreciation since they are naturally related to systems with memory, which exists in most biological systems.Also they are closely related to fractals, which are abundant in biological systems.For these reasons, many phenomena in mathematical biology [1,2,7,12,24,27,28] and some other interdisciplinary fields [5,6,8,[20][21][22][23] can be described very successfully by the models using fractional-order differential equations.Due to the limited theories for analyzing the dynamics of fractional-order systems, stability study of fractional-order population models is only the beginning.In [7], El-Sayed et al. considered the following fractionalorder logistic system: where 0 < q < 1 and c 0 D q t is in the sense of the Caputo fractional derivative, x(t) denotes population density at time t, and r is the intrinsic rate of population or the strength of intraspecific competition of population.They showed existence and uniqueness of solution and studied the asymptotical stability of the positive equilibrium point.However, to the best of the authors knowledge, to this day, still no scholar has investigated the dynamic behavior of fractional-order single-species model in a patchy environment.Motivated by the above considerations, in this paper, we consider the following fractional-order single-species model, which is composed of several patches connected by diffusion: where 0 < q < 1, t 0 0 is the initial time, and c t0 D q t is in the sense of the Caputo fractional derivative, x i (t) denotes population density in patch i at time t, a i is the intrinsic growth rate of population in patch i, b i is the intraspecific competition coefficient of population in patch i, the dispersal rate of population from patch j to patch i is given by d ij .Once an individual enters into a patch, it is counted as an individual in the patch and homogeneously mixed with other individuals in the patch, The contributions of this paper are as follows.Firstly, we prove the non-negativity and boundedness of solutions for a fractional-order population model, as desired in any population dynamics.Secondly, combining Lyapunov method and graph theoretic approach, we obtain a uniform asymptotic stability principle, which has a close relation to the dispersal matrix.Finally, by numerical simulations, we show the effects of dispersal rates on the solutions of the considered system.This paper is organized as follows.In Section 2, some definitions and lemmas are reminded.In Section 3, we prove the existence, uniqueness, non-negativity and boundedness of solutions for system (3).In addition, we consider the stability of the positive equilibrium point for the considered system, and some criteria are derived to ensure the uniform asymptotic stability of the positive equilibrium point.In Section 4, numerical simulations are presented to illustrate the effectiveness of our theoretical results.Finally, conclusions are given in Section 5.

Preliminaries
In this section, we will introduce some definitions and some useful lemmas.It is well known that the initial conditions of fractional differential equations with Caputo derivatives take on the same form as for integer-order ones, which have well-understood practical meanings and more applications in modeling and analysis.Throughout the paper, we use the Caputo fractional-order derivative.Definition 1. (See [26].)The fractional integral with non-integer order q > 0 of function f (t) is defined as follows: [26].)The Caputo derivative of fractional order q of function f (t) is defined as follows: and n is a positive integer such that n − 1 < q < n. (See [17].)The constant x * is an equilibrium point of Caputo fractional dynamic system c t0 D q t x(t) = f (t, x), 0 < q < 1, with the initial condition x(t 0 ) = x t0 if and only if f (t, x * ) = 0.
Remark 1.When 0 < q < 1, it follows from Definition 2 that the Caputo fractional-order system (3) has the same equilibrium points as the the integer-order system To prove the non-negativity and boundedness of solutions of system (3), we need the following lemmas.
Lemma 2. (See [15].)Let V (t) be a continuous function on [t 0 , +∞) and satisfying where 0 < q < 1, λ ∈ R, and t 0 is the initial time.Then For stability analysis of fractional-order single-species model in a patchy environment, we need the following lemmas.Lemma 3. (See [14].)Suppose that the following assumptions are satisfied: ) is a Lyapunov function, and the trivial solution of coupled system on network is uniformly asymptotically stable.
https://www.mii.vu.lt/NARemark 2. If the weighted digraph (G, A) is strongly connected, then c i > 0 for 1 i n, where c i is the cofactor of the ith diagonal element in the Laplacian matrix of (G, A), the basic concepts and notations on graph theory can be seen in [16].
Lemma 4. (See [32].)Let x(t) ∈ R + be a continuous and derivable function.Then, for any time instant t t 0 ,

Main results
In this section, we will prove the existence, uniqueness, non-negativity and boundedness of solutions for system (3).In addition, some sufficient conditions are derived to ensure the uniform asymptotic stability of system (3).
To prove the the existence and uniqueness of the solution for system (3), we need the following lemma.
Proof.We denote x = (x 1 , x 2 , . . ., x n ) and y = (y 1 , y 2 , . . ., y n ) and consider a mapping , where For any x, y ∈ R n , we have where L = max i∈L {a i + 2b i c i + 2 n j=1 d ji }, L = {1, 2, . . ., n}.Thus, H(x) satisfies the Lipschitz condition with respect to x.It follows from Lemma 5 that there exists a unique solution x(t) of system (3) with initial condition x t0 .This completes the proof of Theorem 1.
Considering the biological significance, we are only interested in solutions that are nonnegative and bounded.The following result ensures the non-negativity and boundedness of solutions for system (3).Let R + denote the set of all non-negative real numbers and Theorem 2. All the solutions of system (3), which start in Ω + , are non-negative and uniformly bounded.
Proof.Firstly, we prove that all the solutions of system (3), which start in Ω + , are nonnegative.That is, x i (t) 0 for all i ∈ I = {1, 2, . . ., n} and t t 0 .Suppose that is not true, then there exist a positive t 1 > t 0 and some i ∈ I 1 ⊆ I such that x i (t) > 0, i ∈ I, t 0 t < t 1 , x j (t 1 ) > 0, j ∈ I/I 1 , x i (t 1 ) = 0, i ∈ I 1 , Based on (8) and system (3), we have According to Lemma 1, we have x i (t + 1 ) 0, which contradicts the fact x i (t + 1 ) < 0. Therefore, we have x i (t) 0 for all i ∈ I and t t 0 .
Next, we show that all the solutions of system (3), which start in Ω + , are uniformly bounded.Denote W (t) = n i=1 x i (t), then we have x i https://www.mii.vu.lt/NA Let By Lemma 2, we have Therefore, all the solutions of system (3), which start in Ω + , are confined to the region Γ , where This completes the proof of Theorem 2.
For system (4), we consider a mapping F (x) = (F 1 (x), F 2 (x), . . ., F n (x)), where x = (x 1 , x 2 , . . ., x n ) and The Fréchet derivative of Clearly, F is continuously differential, F (0) = 0, F i (x) 0 for all x ∈ R n + with x i = 0, i = 1, 2, . . ., n.For any F (x) = (F 1 (x), F 2 (x), . . ., F n (x)) and x = (x 1 , x 2 , . . ., x n ), ∂F i /∂x j = d ij 0, i, j = 1, 2, . . ., n, and i = j, so F is cooperative.Since matrix D = (d ij ) n×n is irreducible, it follows that the matrix DF (0) is also irreducible.In addition, it follows from (4) that which shows the boundedness of solutions of system (4).Let s(M 0 ) denote the stability modulus of an n × n matrix M 0 , which is defined by [34,Cor. 3.1], system (4) admits a positive equilibrium point x * = (x * 1 , x * 2 , . . ., x * n ).Since the Caputo fractional-order system (3) has the same equilibrium points as the integer-order system (4), we know that system (3) has a positive equilibrium point x * = (x * 1 , x * 2 , . . ., x * n ).We define a function It can be vertified that V i (x i ) > 0 for all x i > 0 and V i (x i ) = 0 if and only if x i = x * i .Calculating the fractional-order derivatives of V i (x i ) along the solutions of system (3), we have where In fact, 1 − a + ln a 0 for a > 0 and 1 − a + ln a = 0 for a = 1.If D = (d ij ) n×n is irreducible, we know that A = (a ij ) n×n is irreducible.That is to say, the weighted digraph (G, A) is strongly connected, it follows from Remark 2 that c i > 0 for 1 i n.
https://www.mii.vu.lt/NAClearly, we have shown that V i , F ij , G i , and a ij satisfy the conditions of Lemma 3. Therefore, V (x 1 , x 2 , . . ., x n ) = n i=1 c i V i (x i ) is a Lyapunov function for system (3), and c t0 D q t V (x) = 0 implies that x i = x * i for all 1 i n.It follows from Lemma 3 that the positive equilibrium point x * = (x * 1 , x * 2 , . . ., x * n ) of system (3) is uniformly asymptotically stable.Therefore, we thus have the following result.
Remark 3. The dispersal matrix D = (d ij ) n×n is not required to be symmetric, namely, the dispersal rate from patch i to patch j may not be the same as that from patch j to patch i.A typical assumption we impose on dispersal matrix D is that it is irreducible.From biological point of view, this means that population individuals in each patch can disperse between any two patches directly or indirectly.Remark 4. Evidently, the derived result in this paper is still true for q = 1 although the above discussion is based on 0 < q < 1.However, the case q > 1 is not considered in this paper since Lemma 3 is not solved for q > 1, this deserves our further consideration.

Numerical simulations
An effcient method for solving fractional-order differential equations is the predictorcorrectors scheme or more precisely, PECE (Predict, Evaluate, Correct, Evaluate) [3,4], which represents a generalization of Adams-Bashforth-Moulton algorithm.In this section, we use the PECE method for the numerical simulations of the considered model.
In the following, we consider a fractional-order single-species model, which is composed of two patches connected by diffusion as follows: where the parameters are chosen in Table 1.
We keep all parameter values in Table 1 unchanged.By virtue of Matlab, we can compute that x * = (10.48,6.47).Since By simple calculation, c 1 = 2.096, c 2 = 2.588.Obviously, (G, A) is strongly connected, it is very easy to verify that condition s(M 0 ) = 0.8899 > 0 holds.Then it follows from Theorem 3 that system (12) is uniformly asymptotically stable, numerical simulations can be seen in Fig. 1, where the initial values are taken as (8 + 2k, 6 + 2k) (k = 1, 2, 3).From Fig. 1 we find that all solutions of system (12), which through these initial values will converge to the positive equilibrium point Next, we will show the effect of fractional order on the solutions of system (12).We vary the parameter q, while keeping the other parameters fixed as mentioned in Table 1. Figure 2 shows the state trajectories of system (12) with initial value (8,6) when q = 0.4, q = 0.7, and q = 0.98.We can see from Fig. 2 that q is a significant factor, which affects the convergence speed of solutions.

n
j=1 d ij x j denotes https://www.mii.vu.lt/NA the total immigration of population individuals from every other patch to patch i, and n j=1 d ji x i denotes the total emigration of population individuals from patch i to every other patch.Since n i=1 n j=1 d ij x j = n i=1 n j=1 d ji x i always holds, the exchange of population individuals among all patches is always in balance.All parameter values are assumed to be nonnegative and b i > 0 for all i.

Figure 3 . 2 .Figure 4 .
Figure 3. Influence of dispersal rate d 12 on system (12).Here we vary the parameter d 12 , while keeping the other parameters fixed as mentioned in Table 1.(a) Population density in patch 1.(b) Population density in patch 2.

Table 1 .
(12)meter values used in the simulations of system(12).Figure 1.Uniform asymptotic stability of system (12) with parameter values as given in Table 1.(a)Time series of x 1 and x 2 .(b)Phase portrait of x 1 and x 2 .Figure 2. State trajectories of system (12) with parameter values as given in Table 1.(a) Time series of x 1 .(b) Time series of x 2 .