Dynamic properties of a discrete population model with diffusion

We study the dynamical properties of a discrete population model with diffusion. We survey the transcritical, pitchfork, and flip bifurcations of nonhyperbolic fixed points by using the center manifold theorem. For the degenerate fixed point with eigenvalues ±1 of the model, we obtain the normal form of the mapping by using the coordinate transformation. Then we give an approximating system of the normal form via an approximation by a flow. We give the local behavior near a degenerate equilibrium of the vector field by the blowup technique. By the conjugacy between the reflection of time-one mapping of a vector field and the model we obtain the stability and qualitative structures near the degenerate fixed point of the model. Finally, we carry out a numerical simulation to illustrate the analytical results of the model.


Introduction
In recent years, discrete dynamic systems developed rapidly and achieved fruitful results (see [1][2][3][4][5][6][7][8][9]). Particularly, the dynamic properties of discrete population models received extensive attention by many researchers [10][11][12][13][14][15]. In [10] the authors considered the following discrete population model with diffusion: with Dirichlet boundary condition μ t 0 = 0 = μ t n+1 , where d > 0 is the diffusion coefficient, t ∈ Z + = {0, 1, 2, . . .}, p > 1, q > 0, 2 is second-order difference operator, and 2 μ t i-1 = μ t i+1 -2μ t i + μ t i-1 , i ∈ {1, 2, . . . , n}. When n = 2, p = q = b, and μ t 1 , μ t 2 are written as μ t , ν t , respectively, and model (1) can be transformed into a particular two-patch discrete-time metapopulation model of the form (2) © The Author(s) 2020. This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/.
In this paper, we study the following problems. First, we focus on the stability and bifurcations of the nonhyperbolic fixed point E 0 (0, 0) and on the case (i.e., b = (d+1) 2 d-1 , d > 1) of flip bifurcation at nontrivial fixed point E 1 . Second, we investigate the stability and qualitative structures near a degenerate fixed pointẼ 0 (0, 0) formed by fixed points E 0 and E 1 , which collide and bond together when (d, b) = (1,2).
The stability and qualitative structures near a degenerate fixed point of a discrete dynamic system are paid attention by many researchers; for example, Elaydi and Luís [7] proposed some open problems about stability of degenerate fixed point of discrete-time Guzowska-Luís-Elaydi competition model and discrete-time Ricker competition model. Due to the complexity of the degenerate fixed point, it is necessary to combine a variety of mathematical tools to study, for example, coordinate transformation [3], normal form theory (see [3,16]), Picard iteration (see [3,8]), approximation by a flow (see [3,16]), qualitative theory of ordinary differential equations (see [17]), and blowup technique (see [17][18][19]). In this paper, we survey the local behavior near a degenerate equilibrium of a vector field by homogeneous polar blowup (see [18,19]), which can avoid tedious calculations. This paper is organized as follows. In Sect. 2, we give topological types of fixed point of model (2). In Sect. 3, we study transcritical bifurcation, pitchfork bifurcation, flip bifurcation, and the stability of the fixed point E 0 of model (2). In Sect. 4, a particular case (b = (d+1) 2 d-1 , d > 1) of flip bifurcation and stability of fixed point E 1 of model (2) are surveyed. In Sect. 5, we apply the coordinate transformation to obtain the normal form of the mapping of model (2) near the degenerate fixed pointẼ 0 . The approximate system (also called a vector field) of the normal form of the mapping is given by Takens's theorem (see [16, pp. 142-148]) and Picard iteration. The local behavior near the degenerate equilibrium (0, 0) of the vector field is obtained by using the blowup technique. The stability and qualitative structures near the degenerate fixed pointẼ 0 of model (2) are obtained by the conjugacy between reflection of time-one mapping of the vector field and model (2). Finally, we perform numerical simulation analysis.

Topological types of the fixed point
In this section, we give topological types of fixed points E 0 and E 1 on the parameter (d, b)plane. We write model (2) as the following planar mapping T : R 2 → R 2 : .
The mapping T expanded by Taylor series at the fixed point E 0 is where O(3) represents the terms of orders greater than or equal to 3. The Jacobian matrix of T at the point E 0 is Translating the fixed point E 1 into the origin O by coordinate translation y 1 = μμ * , y 2 = νν * , the Taylor expansion at the fixed point E 1 is T : The Jacobian matrix of T at the fixed point E 1 (i.e., mappingT at the fixed point O) is as follows: -3d. For convenience, we give some notations: D i (i = 1, 2, . . . , 6) (see Fig. 1), and ℵ i (i = 1, 2, . . . , 12) (see Fig. 2) are the regions of the parameter (d, b)-plane of E 0 and E 1 , respectively, divided by the curves β i (i = 1, 2, 3, 4) and γ i (i = 1, 2, . . . , 6):

Bifurcation at fixed point E 0
In this section, the parameter b is regarded as a bifurcation parameter. We survey stability and bifurcation properties of the mapping T at the fixed point E 0 . Proof Mapping (3) can be diagonalized by the linear transformation μ = x + y, ν = xy when b = d + 1. We obtain the following mapping: We choose b as a bifurcation parameter. According to the center manifold theorem [20], a center manifold of mapping (5) is where b 0 := d + 1, and ε and δ are sufficiently small positives. Therefore assume that the center manifold form is which must satisfy By comparing the coefficients we obtain d 1 = λ 2 d 1 , e 1 = λ 2 e 1 , f 1 = λ 2 f 1 . We havee λ 1 = 1 and λ 2 = 1 -2d when (d, b) lies on the curve β 1 . Therefore d 1 = e 1 = f 1 = 0. The center manifold is By substituting the center manifold y = h 1 (x, b) into mapping (5) we obtain a onedimensional mapping reduced to the center manifold as follows: We can verify that the transversality and nondegeneracy conditions of transcritical bifurcation (see [4, pp. 504-507]) are established. Therefore the mapping T undergoes a transcritical bifurcation at the fixed point E 0 as (d, b) passes through the curve β 1 . We obtain by Theorem 1.15 in [5, p. 29] and λ 1 = 1, λ 2 = 1-2d = 1 that the fixed point E 0 is stable and nonhyperbolic for 0 < d < 1 when (d, b) lies on the curve β 1 . The fixed point E 0 is unstable and nonhyperbolic for d > 1 when (d, b) lies on the curve β 1 . The proof of the theorem is completed.
The mapping can be seen as Therefore Proof Similarly, the center manifold has the form can be obtained by the same calculation method as in Theorem 1. By substituting the center manifold into mapping (5) we obtain the following one-dimensional mapping reduced to the center manifold: We can verify that the transversality and nondegeneracy conditions of flip bifurcation (see Theorem 4.3 in [3]) are established. Note that the mapping T undergoes a flip bifurcation as (d, b) passes through the curve β 2 . By Theorem 3.5.1 in [6] and F 2 > 0 the mapping T appears in a stable 2-period orbit as (d, b) passes through the curve β 2 . We have λ 1 = -1 and λ 2 < -1 when (d, b) lies on the curve β 2 . So the nonhyperbolic fixed point E 0 is unstable. The proof of the theorem is completed. (2) The fixed point E 0 is nonhyperbolic and unstable when (d, b) lies on the curve β 3 .
Proof Applying the same calculation method as in Theorem 1, we obtain the center manifold (2) yields a one-dimensional mapping reduced to the center manifold We can check that the transversality and nondegeneracy conditions of pitchfork bifurcation (see [4, p. 511]) are established. Therefore the mapping T undergoes a pitchfork bifurcation at fixed point E 0 as (d, b) passes through the curve β 3 . We have λ 1 = 3d + 1 > 1 and λ 2 = 1 when (d, b) lies on the curve β 3 . Therefore the fixed point E 0 is nonhyperbolic and stable when (d, b) lies on the curve β 3 . The proof of the theorem is completed. Proof Similarly, when b = b 3 := 3d -1, we obtain the center manifold as follows: Substituting the center manifold x = h 4 (y, b) into mapping (5), we obtain a one-dimensional mapping reduced to the center manifold, We can check that the transversality and nondegeneracy conditions of flip bifurcation The mapping T undergoes a flip bifurcation as (d, b) passes through the curve β 4 at fixed point E 0 . By Theorem 3.5.1 in [6] the mapping T appears in a stable 2periodic orbit for d ∈ ( 1 3 , 1) and in a repellent 2-periodic orbit for d ∈ (1, +∞). We have λ 1 = 2d -1 and λ 2 = -1 when (d, b) lies on the curve β 4 . By Theorem 1.61 in [5] the fixed point E 0 is unstable and nonhyperbolic for d ∈ (1, +∞) and stable and nonhyperbolic for d ∈ ( 1 3 , 1). The proof of the theorem is completed.

Bifurcation at fixed point E 1
In this section, we consider b as a bifurcation parameter and apply the center manifold theorem to survey the stability and flip bifurcation of the fixed point E 1 of the mapping T. Proof The mappingT can be diagonalized into the following form by linear transformation y 1 Similarly to the previous proof method, we get that the center manifold of mapping (6) is Substituting the center manifold into mapping (6), we obtain the one-dimensional mapping reduced to the center manifold We can verify that the transversality and nondegeneracy conditions of flip bifurcation

Qualitative structures and stability of degenerate fixed point
In this section, we apply Picard iteration, approximation by a flow, and near identity transformation (see [3,8,16]) to investigate the degenerate fixed point of the mapping T when (d, b) = (1, 2). We obtain qualitative structures and the stability of the degenerate fixed point of model (2).

Theorem 6
The smooth mapping F can be transformed into the following normal form by reversible near identity transformation: to mapping F to eliminate quadratic and cubic items as much as possible. Then the resonance term is left, and the normal form is obtained. The proof is complete.

Approximation by a flow
For dealing with stability and qualitative structures near the degenerate fixed points of a planar mapping with eigenvalues ±1, the general idea is embedding the reflection of the normal form of a mapping into the flow of a vector field. Then we study the dynamic behavior of the mapping by the properties of the vector field. In [3, p. 424], it is stated that "any mapping sufficiently close to the identity mapping can be approximated up to any order by a flow shift. " That is, an approximating system (see [3]) can be constructed, in which the unit time shift ϕ 1 along orbits coincides with (or, at least, approximates) the near identity mapping up to order k. In this paper, we use the reflection transformation (see [3,8]) to transform the truncated normal form of the mappingF into a near identity mapping and embed it into the flow of a vector field. We use the vector field to investigate the stability and qualitative structures near the degenerate fixed point of the mapping T when (d, b) = (1, 2). The Taken theorem and the following lemma come from [21, p. 38].
Lemma 3 (Taken's) Let F : R n → R n be the C r -diffeomorphism (r ≥ 2) defined by where A = S(I + N), S is semisimple, N is nilpotent, SN = NS, and I is the identity. Let 1 < l < r be any integer, and let F k ∈ H k n , where H k n is the vector space of homogeneous polynomials of order k in n variables with values in C n . Then there exist a diffeomorphism ψ l : ⊆ R n → R n , where is a neighborhood of the origin in R n , and a vector field X(x) on R n such that Sx)), where j l is the truncation operator up to order l, and X (t, x) is the flow of X(x). Furthermore, for such a vector field X(x), j l X(x) is uniquely determined by j l F(x).
If there exist constants k, σ , η such that then the C ∞ vector field is of Łojasiewicz type, where · is the Euclidean norm on R 2 . (2) with associated formal normal formG = R •X 1 (X → X 1 , i.e., time-one mapping), whereX has in 0 a singularity of Łojasiewicz type with a characteristic orbit and R n = I for some n > 0. Suppose, moreover, that all singularities in some nice decomposition ofX are of type I and, restricted to a fundamental conic domain , we have α(g n | ) ∈ [Diff 0 rot (2)] k (this especially is the case whereX has no elliptic sectors). Then there exists a C ∞ -coordinate change H with j ∞ (H -I)(0) = 0 and an R-invariant representative X ofX such that

Lemma 4 Let g ∈ Diff
Remark 2 Lemmas 3 and 4 deal with many concepts, the details of which can be found in [16] and [21], respectively.
For the truncated normal form of the mappingF, we can see that and RN(x) = N(Rx). For R • N embedded in the flow of a vector field, we have the following results.

Theorem 7 The mapping R • N satisfies
in a small neighborhood of the origin O, where ϕ t is the flow generated by the planar vector field The ϕ 1 is the unit-time shift along orbits of planar vector field (10), and R = diag(1, -1), Proof By Lemma 3 the reflection of the truncated normal form N can be embedded in the flow of a vector field, and ϕ t is constructed as follows: where J = 0 and Y k (X) = (Z k 1 (X), Z k 2 (X)) T , where Z k 1,2 (X) are homogeneous polynomial functions with unknown coefficients. To facilitate the expression of the Picard iteration process, define M(X) = N(η), which is introduced to solve the explicit expression of the vector field Y such that RM(X) = φ t (X) + O( X 4 ). Because the highest order of N(η) expressions is 3, it is only necessary to perform three Picard iterations on system (11) to achieve the approximate purpose. We start the iterative from X (1) (t) = e Jt X. The calculation shows that the linear part of X (1) (1) coincides with RM(X). The expression of Y 2 can be assumed of the following form Therefore by Picard iteration we obtain By comparing the quadratic terms of RM(X) and X (2) (1) we obtain A 20 = 1, A 02 = 32, and B 11 = -1. Then we put By iterating and letting t = 1 we get By comparing the cubic terms of MR(X) and X (3) (1) we obtain By substituting the expressions of Y 2 , Y 3 into system (11) we obtain a planar vector field (10) and the following C ∞ vector field: We can see from this expression that the C ∞ vector field X is of Łojasiewicz type. The proof of the theorem is completed. Proof Vector field (10) is η 2 → -η 2 invariant, and equilibrium (0, 0) is degenerate (i.e., equilibrium (0, 0) without linear part), so the phase portraits of vector field (10) are symmetric about the η 1 -axis, so we consider only η 2 ≥ 0. To study the local behavior near the degenerate equilibrium (0, 0) of the vector field X, we will desingularize it by using the homogeneous polar blowup technique. Consider the mapping This mapping transforms {r = 0} into (0, 0), and inverse mapping φ -1 blows up the degenerate equilibrium (0, 0) of the vector field X to a circle (see [18, pp. 91-92]).
By using coordinate transformation η 1 = r cos θ , η 2 = r sin θ of vector field X, we obtain the system .
To investigate the phase portrait of the vector field X in a neighborhood of the origin O, it clearly suffices to investigate the phase portrait of system (13) in the neighborhood φ -1 ( ) of the circle S 1 × {0}. The phase portrait of the vector field X near the degenerate equilibrium (0, 0) is easily obtained by shrinking the circle S 1 × {0} to a point. Since dθ dt θ=0,π > 0, dr dt θ=0 > 0, dr dt θ=π < 0, the trajectory of the vector field X on the η 1 -axis tends to +∞ when θ = 0 and η 1 > 0 and tends to the origin O when θ = π and η 1 < 0. The vector field X has an orbit connecting with origin O in a definite direction. If r is sufficiently small, then dr dt > 0, θ ∈ 0, π 2 and dr dt < 0, θ ∈ π 2 , π .
Therefore we can obtain the phase portrait of system (13) on the circle S 1 × {0}. By shrinking the circle S 1 × {0} to a point we can obtain the local behavior of the vector field X near the degenerate equilibrium (0, 0). Figure 3 also can be seen as the phase portrait of the vector field X near the degenerate equilibrium (0, 0). From the previous proof process we can see that the vector field X has no elliptic sectors and that the vector field X obtained by desingularization has only hyperbolic equilibrium. Therefore the vector field X satisfies the assumption of Lemma 4. From Lemma 4 we see that there is a C ∞ diffeomorphism satisfying j ∞ ( -I) = 0 such that at the origin O. Thus the degenerate fixed pointẼ 0 of model (2) is unstable. We see from [3] that the orbit unit-time shift ϕ 1 along vector field (10) coincides with the near-identity mapping at least up to order 3. For a sufficiently small initial value P(μ 0 , ν 0 ) ∈ R 2 + , by means of C ∞ diffeomorphism , linear transformation 1 , and the local behavior of the vector field X near the origin O we can see that the sequence {T n (P)} +∞ 0 enters along both sides of straight line μ = ν "repel" and then "attract". The proof is completed.
Remark 3 For dynamic properties of the normal form N , the general description is given in [8] when a 0 b 0 = 0. In this paper, the normal form N satisfies a 0 = 1 2 , b 0 = 16, which is one of the situations given in [8].
Remark 4 In [17] the author investigates qualitative properties of the vector field near degenerate equilibrium by desingularization. This desingularization is obtained by homogeneous directional blowup. In this paper, we realize desingularization by homogeneous polar blowup. As can be seen from Remark 3.1 in [19], the homogeneous directional blowup and the homogeneous polar blowup are equivalent.

Numerical simulation
In this section, we select appropriate parameters and initial values for the numerical simulation to verify our results. The details are as follows. is constantly "oscillating" along the straight line μ = ν, which also confirms Theorem 8.