DYNAMICS OF A THREE SPECIES RATIO-DEPENDENT FOOD CHAIN MODEL WITH DIFFUSION AND DOUBLE FREE BOUNDARIES

. This paper focuses on the dynamics of a three species ratio-dependent food chain model with diﬀusion and double free boundaries in one dimensional space, in which the free boundaries represent expanding fronts of top predator species. The existence, uniqueness and estimates of the global solution are discussed ﬁrstly. Then we prove a spreading–vanishing dichotomy, speciﬁcally, the top predator species either successfully spreads to the entire space as time t goes to inﬁnity and survives in the new environment, or fails to establish and dies out in the long run. The long time behavior of the three species and criteria for spreading and vanishing are also obtained. Besides, our simulations illustrate the impacts of initial occupying area and expanding capability on the dynamics of top predator for free boundaries.


Introduction
In ecosystems, Prey-predator interaction is one of the basic interspecies relations, and it is also the basic block of more complicated food chain, food web, and biophysical network structures (e.g. [16]). A general food chain model is described by the system of differential equations      N t = N ϕ( N ) − 1 e1 ρ 1 ( N , P ) N , P t = ρ 1 ( N , P ) P −m 1 P − 1 e2 ρ 2 ( P , Z) P , Z t = ρ 2 ( P , Z) Z −m 2 Z, (1.1) where N , P , Z stand for the population density of prey, predator, and top predator, respectively (cf. [1,4,29]). ρ 1 ( N , P ) and ρ 2 ( P , Z) are the functional response of predators; e i , i = 1, 2 are the yield constants;m i , i = 1, 2 are the normal death rates of predators. The function ϕ( N ) represents the per capita growth rate of N . A Table 1. Parameters table of problem (1.2).

Parameters
Biological meanings r The prey intrinsic growth rate K The prey carrying capacity a 1 (resp., a 2 ) The maximum growth rate of predator (resp., top predator) e 1 The yield conversion factor for predator feeding on the prey e 2 The yield conversion factor for top predator feeding on the predator b 1 (resp., b 2 ) The half saturation constant for predator (resp., top predator) m 1 (resp.,m 2 ) The death rate of predator (resp., top predator) d 1 (resp., d 2 , d 3 ) The diffusivity of the prey (resp., predator, top predator) distinct feature of the simple food chain is secondary extinction because of trophic cascading effect [13,14]. That means if one species dies out, all the species at higher trophic levels die out as well. Especially, there has been considerable interest in ratio-dependent food chain model for last two decades [9,10,17,18,20,29]. Peng et al. [18] investigated the following coupled reaction diffusion food chain model with ratio-dependent functional response: N (0, x) =φ 1 (x), P (0, x) =φ 2 (x), Z(0, x) =φ 3 (x), x ∈ Ω. (1.2) Here, Ω is a fixed bounded open domain in R N with smooth boundary at any given time, and ν is the outward unit normal vector on ∂Ω. For i = 1, 2, 3, the initial data are continuous functions and satisfyφ i (x) ≥ 0, ≡ 0 for x ∈ Ω. Parameters values are assigned according to biological meanings in Table 1.
Taking into account the spatial spreading, we suppose that the prey and predator are native and the top predator is invasive. At initial time (t = 0), the prey and predator live in the whole line R, the top predator exists in a bounded interval and spreads from double boundaries to expand its habitats. That is, it will move outward along the unknown curves (free boundaries) as time t increases. To understand the dynamics of these three species and free boundaries, we formulate the following free boundary problem Z(t, x) = 0, t ≥ 0, x / ∈ (g(t),h(t)), g (t) = −μ Z x (t,g(t)),h (t) = −μ Z x (t,h(t)), t ≥ 0, whereg(t) andh(t) represent the left and right moving boundaries, respectively; r, K, d i , a i , b i , e i ,m i , µ, andh 0 are given positive constants for i = 1, 2. The free boundary conditionsg (t) = −μ Z x (t,g(t)) and h (t) = −μ Z x (t,h(t)) are special cases of the well-known Stefan type conditions, and their ecological background and derivation can be referred to [2,15,28]. It is worthwhile to mention free boundary problems have been introduced to describe the population spreading process. This paper can be seen as further research of [33], in which the authors studied a ratio-dependent preypredator model with a free boundary. More precisely, Du and Lin [6] firstly discussed a free boundary problem for the diffusive one-species Logistic model. They investigated the existence and uniqueness, regularity and uniform estimates, and long-time behavior of global solution. Some sufficient conditions for spreading and vanishing were also established. When spreading occurs, they provided some estimates of asymptotic spreading speeds of the free boundary. For additional relevant works in this direction can be found in [2,3,5,7,12,15,[21][22][23]. Free boundary problems for the diffusive two-species competition models and prey-predator models have been systematically investigated in literatures, refer to [8,24,25,27,28,30,33,34]. Recently, a free boundary problem for the diffusive three-species model has been studied in [31].
For simplicity, we non-dimensionalize problem (1.3) with the following scaling: Dropping the caps of t and x, then problem (1.3) is rewritten in the following form: Furthermore, we suppose that φ 1 (x), φ 2 (x), and φ 3 (x) satisfy is the space of continuous and bounded functions in R. Inspired by [10,29], we may define that E 0 = (0, 0, 0) is an trivial steady state of (1.4). Besides, problem (1.4) also admits the semi-trivial steady states E 1 = (1, 0, 0) and if and only if 0 < c 1 (k 1 − m 1 ) < k 1 . It is easily shown that (1.4) has a unique positive constant steady state if and only if the following are satisfied: Moreover, the unique positive constant steady state can be expressed as In addition, we note that k 2 > m 2 and A > 1 imply k 1 > m 1 .
In the absence of diffusion, Hsu et al. [10] dealt with problem (1.2). Obviously, the steady states E 0 and E 1 always exist. For some suitable conditions, they concluded that E 0 is globally asymptotically stable when c 1 > 1 and E 1 is globally asymptotically stable when c 1 < 1, please refer to Theorems 2.1 and 2.3 in [10]. Then, Peng et al. [18] investigated the stationary pattern of problem (1.2). They obtained that E 1 is globally asymptotically stable if c 1 ≤ 1 and k 1 ≤ m 1 , see Proposition 3.1 in [18]. As shown in Proposition 3.1 in [18], if k 1 ≤ m 1 or k 2 ≤ m 2 , the two predators or the top predator will face extinction, respectively. In this work, we mainly analyze that the influence of free boundaries on the spreading process of the top predator. Therefore, we only consider the assumptions that k 1 > m 1 , k 2 > m 2 , and c 1 < 1.
The organization of this paper is organized as follows. In Section 2, we first give the local existence and uniqueness of the solution to (1.4) and then show that the solution exists for all time t ∈ (0, ∞). Section 3 is devoted to the long time behavior of (N, P, Z). In Section 4, we shall provide the criteria governing spreading and vanishing. Furthermore, two numerical examples to illustrate the theoretical analysis are carried out in Section 5. The last section is a brief conclusion and discussion.

1)
here positive constants C and T only depend on Σ, α, and p.
Proof. The proof is motivated by [23,28,32]. Firstly, we straighten the free boundaries. Let so that the transformation (t, x) → (t, y) helps us straighten the left free boundary x = g(t) to the line y = −1 and the right free boundary x = h(t) to the line y = 1. Set Therefore, the free boundary problem (1.4) is equivalent to where , The rest of the proof is similar to that of Theorem 1.1 in [23] and Theorem 2.1 in [28]. By using the upper and lower solutions method, we may get the local existence of the solution (Q, R, W ) to the initial-boundary value problems (2.2) and (2.3) with fixed boundary for given g(t) and h(t). Then, the local existence and uniqueness of the solution (Q, R, W, g, h) to (2.2), (2.3), and (2.4) is obtained by the contraction mapping theorem, L p theory, Schauder theory and the Sobolev imbedding theorem [11]. See also Du and Lin [6], Wang and Zhang [24]. We omit the details.
To show the global existence of the solution, we need the following priori estimates. .
Proof. Firstly, as long as the solution exists, one can see that N, , h(t)) by the strong maximum principle. Note that N , P and Z satisfy In view of the comparison principle, and yet straighten the free boundary x = h(t) to the line y = 1. Hence, Q, R and W obey To verify this assertion, compare Z with the following auxiliary function. Define and construct an auxiliary function Direct calculations give that, for (t, x) ∈ Γ, Therefore, we observe that Consequently, we expect to find some M 4 independent of T such that Thus, from (2.5), (2.6), (2.7), and the comparison principle, we may deduce V ≥ Z in Γ. Then, It remains to show that (2.7) holds. To this aim, choose Theorem 2.3. Assume that k 1 > m 1 and k 2 > m 2 . Then the solution (N, P, Z, g, h) of (1.4) exists and is unique for all t ∈ (0, ∞).
Proof. Let T max be the maximum existence time of (N, P, Z, g, h). Taking advantage of Theorem 2.1, we observe T max > 0. We now wish to show T max = ∞. Suppose that T max < ∞. In terms of Lemma 2.2, we have Choose σ ∈ (0, T max ). Arguing as the same as the proof of Theorem 2.1, one can see that . Then, following the proof of Theorem 2.1 again, there exists τ > 0 independent of T max such that the solution of (1.4) with initial time T max − τ 2 can be extended uniquely to the time T max − τ 2 + τ . This means that the solution of (1.4) can be extended as long as N , P and Z remain bounded, which is in contradiction with our hypothesis. Thus, T max = ∞. The proof is complete.
3. Long time behavior of (N, P, Z) Then, the subsequent analysis is divided into two parts: vanishing case (h ∞ − g ∞ < ∞) and spreading case We now present the following results that will be used frequently.
This implies The proof of Lemma 3.2 is similar to that in Lemmas A.5 and A.4 in [33] with some modifications and we omit it here.
To discuss the vanishing case, we first exhibit the following estimate and remember a basic result.
Proof. By similar arguments in the proof of Theorem 2.1, we adopt the following transformation to straighten the free boundaries x = h(t) and x = g(t) to the lines y = 1 and y = −1. The direct calculations show that R, W obey (2.3) and g, h satisfy (2.4). Consequently, in light of Proposition 7.1 from [28], there exists a constant ξ 0 > 0 such that holds. Then, we may find a constant M > 0 such that −M ≤ g (t) < 0 as in Lemma 2.2 and yet (2.4) yields Combining all above, there exists a constant ξ 1 > 0 such that We have g ∞ > −∞ and h ∞ < ∞ according to the condition h ∞ − g ∞ < ∞. Assume that lim t→∞ g (t) = 0 does not hold. Then, for some δ > 0, there exists a sequence {t n } such that t n → ∞ and g (t n ) → −δ as n → ∞. In terms of (3.4), we choose > 0 so small that In addition, p ((0, T ) × (0, ρ(t))) for some p > 1 and any T > 0, and for some constant ξ > 0. Then Proof. First, from Theorem 3.3 and Lemma 3.4, we deduce (3.8). Fix 0 < η 1. In view of Z(t, x) = 0 for t > 0 and x / ∈ (g(t), h(t)), there exists T 0 > 0 such that Step 1: The constructions of N 1 and P 1 .
Step 5: Continuing the above routine, we obtain four sequences {N i }, {P i }, {N i } and {P i }, whose monotonicity is straightforward; namely, for all i = 2, 3, · · · , uniformly on the compact subset of R. Furthermore, these sequences can be determined by the following iterative formulas: Because of the constant sequences {N i }, {P i } are monotone non-increasing and bounded, {N i }, {P i } are monotone non-decreasing and bounded. Hence, all of those sequences limits of exist. Let us denote their limits by N , P , N and P , respectively. Then, This, together with the conditions k 1 > m 1 and c 1 < 1, yields N = N = N * 2 and P = P = P * 2 . Thus, from (3.18) and (3.19), we get (3.6) and (3.7).

Spreading case
For convenience of later analysis, we now present the following comparison principle which can be proved by the same way as that of Lemma 3.5 from [6] and the details will be omitted.

The criteria governing spreading and vanishing
In this section, we focus on the criteria governing spreading and vanishing. We first give a necessary condition for vanishing.
Proof. In view of Theorem 3.5, the condition h ∞ − g ∞ < ∞ implies that lim t→∞ Z(t, x) C([g(t),h(t)]) = 0, lim t→∞ N (t, x) = N * 2 and lim t→∞ P (t, x) = P * 2 uniformly in any compact subset of R. Suppose that h ∞ − g ∞ > Λ. For any given 0 < δ 1, there exists T 1 such that Let ω(t, x) be the solution of the following initial-boundary value problem Then, by using the comparison principle, we deduce . Due to h(T ) − g(T ) > Λ, one can easily see that ω(t, x) → ϑ(x) as t → ∞ uniformly in the compact subset (g(T ), h(T )), where ϑ(x) is the unique positive solution of which is in contradiction with (3.8). So, (4.1) holds. We now pay close attention to the case h 0 < Λ 2 . Lemma 4.3. Assume that k 1 > m 1 and k 2 > m 2 . Let (N, P, Z, g, h) be the unique global solution of (1.4). If h 0 < Λ 2 , then there exists µ * > 0 depending on (φ 1 , φ 2 , φ 3 , h 0 ) such that g ∞ = −∞ and h ∞ = ∞ if µ ≥ µ * . Proof. The proof is inspired by Lemma 3.6 from [19]. Take M = max{M 1 , M 2 , M 3 } is defining as in the proof of Theorem 3.5. Observing the third equation of (1.4), for all P, Z ∈ [0, M ], we may find a positive constant δ * such that Consider the following auxiliary free boundary problem Arguing as in the proof of the existence and uniqueness of the global positive solution to (1.4), we deduce that (4.2) also admits a unique positive solution (ω, ρ 1 , ρ 2 ) defined for all t > 0. In addition, ω is bounded above by some positive constants and ρ 1 (t) < 0 < ρ 2 (t) for all t > 0. Now we write (Z µ , g µ , h µ ) and (ω µ , ρ µ 1 , ρ µ 2 ) in place of (Z, g, h) and (ω, ρ 1 , ρ 2 ) to clarify the dependence of the solutions on the parameter µ.
Applying the comparison principle, we derive for all 0 < µ ≤ µ * . Therefore, On the other hand, Z(t, x) satisfies Thus, we conclude from Lemma 3.6 that This implies that The proof is finished.

Numerical simulations
In this section, two numerical examples are carried out by Matlab to support our theoretical results. To be computable, we fix some coefficients and initial functions as follows: D 1 = 0.5, D 2 = 0.8, c 1 = 0.838, c 2 = 0.52, k 1 = 1.5, k 2 = 2, m 1 = 0.61, m 2 = 1.05; It clearly follows from the above parameters that k 1 > m 1 , k 2 > m 2 , k 1 > m 1 + c 2 and c 1 < 1. This means that the conditions are maintained in Theorems 3.5 and 3.8. Then, problem (1.4) admits the semi-trivial steady state Now, we assume that the initial habitat h 0 and expanding capability µ are both changeable and design the following experiments:  k2−m2 ≈ 1.61. Furthermore, Theorem 4.5 shows that if h 0 < Λ 2 then vanishing happens with small expanding capability µ and spreading happens with large expanding capability. From Figure 1a, we observe that the smaller the expanding capability µ is, the slower the free boundaries x = g(t) and x = h(t) increase, and Z(t, x) decays to zero quickly. The numeric result shows that the top predator Z(t, x) will eventually go to extinction, agreeing with Theorem 4.5. However, in Figure 1b, the larger the expanding capability µ is, the faster the free boundaries x = g(t) and x = h(t) increase, and Z(t, x) stabilizes to a positive solution Z * 3 . This means that the top predator Z(t, x) will eventually survive in the new environment. The numerical simulations illustrate the corresponding ecological explanation as follows. If the initial habitat h 0 is less than expanding barrier Λ 2 , then the expanding capability µ dramatically affect the top predator Z(t, x) whether spreading or vanishing occurs.
Example 5.2. Take h 0 = 2 > Λ 2 . Comparing with Figure 1a, we can see that in Figure 2 the free boundaries x = g(t) and x = h(t) increase fast, and Z(t, x) stabilizes to a positive solution Z * 3 . It clearly shows that the top predator Z(t, x) successfully spreads to the entire space as time goes to infinity. The corresponding ecological explanation depict that when the initial habitat h 0 is more than expanding barrier Λ 2 , the top predator Z(t, x) spreading always happens regardless of the expanding capability µ. It turns out that Remark 4.2 is valid.

Conclusion and discussion
In this work, we have studied the dynamics of a three species ratio-dependent food chain model with diffusion and double free boundaries in one dimensional space, in which the free boundaries represent expanding fronts of top predator species. We are mainly concerned with the case of a weak predation rate for the prey species (i.e., c 1 < 1). It is assumed that the top predator initially occupies a finite region [−h 0 , h 0 ] and has a tendency to expand its territories. The dynamic behavior of (1.4) exhibits a spreading-vanishing dichotomy, i.e., as t → ∞, either (i) the top predator Z successfully establishes itself in the new environment in the sense that g(t) → −∞ and h(t) → ∞. Moreover, N → N * 3 , P → P * 3 and Z → Z * 3 when k 2 > m 2 , k 1 > m 1 + c 2 and c 1 < 1 hold (see Thm. 3.8); or (ii) the top predator Z spreads within a bounded area and vanishes eventually, i.e., h(t) − g(t) → h ∞ − g ∞ ≤ π 1 k2−m2 and Z(t, x) C([g(t),h(t)]) → 0 when k 1 > m 1 , k 2 > m 2 and c 1 < 1 hold. Furthermore, N → N * 2 and P → P * 2 uniformly in any compact subset of R (see Thm. 3.5). The phenomena look more realistic and may play an important role in the understanding of ecological complexity.
To the best of our knowledge, Hsu et al. [10] formulated the ratio-dependent food chain model and thought that N as a plant species, P as pest species and Z as a species used to control the pest. The purpose of the work was to study its applications to biological control from a mathematical point of view. Subsequently, a growing number of papers on the model have been published (e.g. [9] and [18]). A notable problem is under what conditions, all three species coexist. That is, the introduction of a natural enemy Z for pest P helps the population level of prey N . The authors gave some suitable coefficients conditions to ensure that all three species coexist. However, for problem (1.4) which is closer to reality than (1.2), main conclusions tell us that in order to control the pest species we further should put natural enemy at the initial state at least in one of the following three ways: (i) expand the initial habitat of natural enemy; (ii) increase the expanding capability; (iii) augment the initial density of the natural enemy.
The results in the present paper may enrich biological control in ecosystems.