Spatiotemporal dynamics for a diffusive mussel–algae model near a Hopf bifurcation point

In this paper, the Hopf bifurcation and Turing instability for a mussel–algae model are investigated. Through analysis of the corresponding kinetic system, the existence and stability conditions of the equilibrium and the type of Hopf bifurcation are studied. Via the center manifold and Hopf bifurcation theorem, sufficient conditions for Turing instability in equilibrium and limit cycles are obtained, respectively. In addition, we find that the strip patterns are mainly induced by Turing instability in equilibrium and spot patterns are mainly induced by Turing instability in limit cycles by numerical simulations. These provide a comprehension on the complex pattern formation of a mussel–algae system.


Introduction
As early as 1952s, Turing [1] found that the pattern formation in biological systems can be expressed by the coupled reaction-diffusion equations. Turing theory shows that pattern formation is possible if two substances interact differently in their diffusion range.
Since then, the dynamics of pattern formation has attracted wide attention in biology [2][3][4], chemistry [5], aquatic ecosystem [6,7], embryogenesis [8][9][10], and so on. With the application of spatial patterns in multifarious ecosystems systems, the stability of pattern solutions [11,12], the stability of the positive equilibrium solutions, the Hopf bifurcation at the critical point, the Turing bifurcation in different reaction-diffusion systems, and the dynamical properties in Turing-Hopf bifurcation point have been getting a significant progress on the stability of pattern formation [13]. As well as implication of the mussel-algae bed patterns, Van De Koppel [14] proposed a simple mathematical model to illustrate the pattern formation process, and they investigated the influences of selforganization on the emergent properties of nature systems. The reaction-diffusion equations for the original mussel-algae model with weakly nonlinear diffusive instability have Although algae are considered advection of tidal currents on a large scale, they are actually dispersed in the fluid in the form of Brownian particles on a small scale [13]. Numerous laboratory results show that mussels can actively move both within and between clusters, the influence of advection with tidal flow at a small-scale space on the mussel bed is very small. These imply that advection and diffusion are two different ecological processes; in real mussel bed ecosystems, these processes normally coexist and share the same activator-inhibitor mechanism [13,16]. Advection and diffusion are equivalent due to the emergence of spatial self-organizing patterns. Therefore, we study system (2) with self-diffusion and Neumann boundary conditions in one-dimensional space x, which has the following form: Models (1)- (3) have been investigated by scholars. Wang et al. [17] derived the conditions for differential-flow instability that causes the formation of spatial patterns, and then systematically investigated the influence of parameters on pattern formation. Liu et al. [16] verified that the dimensionless model is an extension of the original model and in exact accordance with their laboratory experiment. Ghazaryan and Manukian [15] used the geometric singular perturbation theory to study the nonlinear mechanisms of pattern and wave formation. Sherratt and Mackenzie [18] considered the influence of advection of algae for pattern formation. Song et al. [13] applied the normal form method to study the Turing-Hopf bifurcation with Neumann boundary conditions. Based on the previous research, it has been observed that (3) possesses a spatially homogeneous periodic solution for some parameter ranges, and the interaction between Hopf bifurcation and Turing bifurcation may be the driving force for more complex spatiotemporal phenomena for (3) [13].
For the Hopf bifurcation of a reaction-diffusion system, Yi et al. [19] determined the direction of Hopf bifurcation and the stability of bifurcating periodic solutions for a reaction-diffusion Lengyel-Epstein model. Particularly, they showed the existence of multiple spatially non-homogeneous periodic orbits while the system parameters are all spatially homogeneous. Early, Yi et al. [20,21] considered a reaction-diffusion Lengyel-Epstein system, they investigated the global asymptotic behavior and gave Turing instability conditions for the spatially homogeneous equilibrium solution and the spatially homogeneous periodic solution. Similarly, Liu and Yi et al. worked along this line and have done the Gierer-Meinhardt model and the bimolecular model [22]. Recently, Wang and Hou et al. [23,24] developed this method and have done the Gierer-Meinhardt model and the tussock wedge model. They showed sufficient conditions for Turing instability in the spatially homogeneous equilibrium solution and the spatially homogeneous periodic solution, respectively.
The purpose of this paper is to study the stability of the stable positive equilibrium as a spatial homogeneous steady solution, and the stability of the periodic orbit as a spatial homogeneous periodic solution of the reaction-diffusion mussel-algae model. In this article, we analyze the stability and Hopf bifurcation of the positive equilibrium in both ODE (the corresponding ordinary differential equation model of system (3)) and PDE (the partial differential equation system (3)) models and determine the conditions for the direction of bifurcation and the stability of the bifurcating periodic solution. For PDE system (3), we also derive the precise conditions on the parameters so that the spatially homogeneous steady solution and the spatially homogeneous periodic solution become Turing unstable. By both theoretical analysis and numerical simulations, we show the coexistence of an unstable equilibrium point and a stable limit cycle for ODE, and a Turing unstable periodic solution is attracted by a stable non-constant steady state for PDE. The discussion of the pattern formation mechanism in equation (3) mainly concerns the dynamics of the corresponding ODE system and the spatial diffusion system. The paper is organized as follows: In Sect. 2, the existence and stability of the positive equilibrium and Hopf bifurcation are analyzed. In Sect. 3, sufficient conditions for Turing instability in equilibrium and spatial homogeneous periodic solutions are given theoretically, respectively. In Sect. 4, some simulations are given to illustrate the theorems in Sects. 2 and 3.

Local stability analysis and Hopf bifurcation for the ODE system
For the reaction-diffusion mussel-algae system (3), the local system is an ordinary differential equation, which has the following form: The dynamics of the ODE system mainly includes the stability of the equilibrium and Hopf bifurcation. In the following, we analyze these dynamics in detail.

Local stability analysis of the equilibria
The equilibria of system (4) satisfy the following equation: About the analysis of the equilibria of system (5), we have the following proposition.

Hopf bifurcation analysis
From Theorem 2.1, we can find that the stable equilibrium (m * , a * ) becomes unstable when γ increases from the left to the right of γ H . Hence, let γ be the main bifurcation parameter, the direction and the stability of Hopf bifurcation will be investigated in the following. Mathematically, this parameter shows an equal role to other parameters in system (4), which can also be chosen to make similar demonstration. Firstly, we transform equilibrium from (m * , a * ) to (0, 0), namely let w = mm * and z = aa * , then equation (4) becomes where Based on [25,26], the occurrence of Hopf bifurcation needs the nonzero transversality > 0. Next, we study the stability and direction of the Hopf bifurcation near the critical value γ = γ H . When γ = γ H , we obtain that λ 1,2 = ±iω 0 , and ω 0 = (r-1) T , then system (9) will transform to the following form: where In order to determine the type of Hopf bifurcation at (m * , a * ), we should determine the sign of discriminatory quantity σ . According to [26,27], σ has the following form: By further computation, we have σ = αr 2 8(1-rα) [(1r) + r(rα -1)] < 0. Through the above analysis, we can get the following theorem.
The Hopf bifurcation is supercritical due to σ < 0 and the stable limit cycle will bifurcate from (m * , a * ) at the right of γ H .

Turing instability analysis
In this section, we study the impact of diffusion term on the stability of (m * , a * ) and the limit cycle of system (3), respectively. It should be noted that the solution to system (4) is also a solution of the reaction diffusion system (3) [28,29].
Let u 1 = mm * and u 2 = aa * , then system (3) is transformed into with the boundary conditions The solution of (11)- (12) in H 2 ((0, π)) × H 2 ((0, π)) has the following form [3,4]: where λ k ∈ C is the temporal spectrum, k is the spatial spectrum, M k and A k are complex numbers for k = 0, 1, 2, . . . . Substituting (13) into (11), then Collecting items of the same powers about k, then where The characteristic equation of J k is where If Tr(J k ) < 0 or Det(J k ) > 0 does not hold for some k, Turing instability will occur in equilibrium [8,[31][32][33]. In the following, we give conditions for the occurrence of Turing instability in (m * , a * ). If Det(J k ) > 0 for all k = 0, 1, 2, . . . , (m * , a * ) is stable for system (3), and there are no Turing patterns occurring. If Det(J k ) < 0 for at least one k ∈ N , (m * , a * ) is unstable for system (3), and Turing patterns occur. Through the above analysis, we can get the following theorem.
, m ∈ N + , then (m * , a * ) is asymptotically stable for system

Turing instability of limit cycle
In this subsection, we investigate the stability of limit cycle with spatially heterogenous perturbations. The reaction-diffusion system (3) possesses any periodic solution of system (4) as a spatially homogeneous periodic solution, including the ones bifurcating from Hopf bifurcation. We can also perform a Hopf bifurcation for system (3) at the same bifurcation point in (4), and bifurcating spatially homogeneous periodic solutions exist near γ = γ H [28,29]. But the stability of these periodic solutions with respect to (3) could be different from that for (4). That is, if φ(t) is an unstable periodic solution of (4), then it is clearly also unstable for (3). While if φ(t) is a stable periodic solution of (4), it could be unstable for (3) because of diffusion. For example, if (m * , a * ) is an unstable equilibrium solution of (4) but stable for (3), then the nearby bifurcating periodic solutions through Hopf bifurcation are also unstable, which illustrates the interaction of Hopf instability and Turing instability [20].
Suppose that (H2) and (H3) are satisfied, so the supercritical Hopf bifurcation appears at γ = γ H . It can be confirmed that the limit cycle is stable under spatially homogeneous perturbation.
The linearized system of (17) at the equilibrium (0, 0) is with the boundary conditions The solutions of (18) in H 2 ((0, π)) × H 2 ((0, π)) have the following form: where λ(k) ∈ C is the temporal spectrum, k is the spatial spectrum, m k and a k are complex numbers for k = 0, 1, 2, . . . . Plugging (μ 1 , μ 2 ) into (18), then Collecting terms with the same powers about k, we get where The corresponding characteristic equation of L k has the following form: where If conditions (H2) and γ = γ H are established, we have Tr(L 0 ) = 0, Det(L 0 ) = ω 2 0 > 0 and Tr(L k ) < 0 for k = 1, 2, . . . . Then it follows that, for k = 0, L 0 has eigenvalues with zero real parts. We need to proceed to the center manifold reduction.

Theorem 3.2 When conditions
then the spatially homogeneous periodic solution for system (3) is stable if condition (H8) or (H9) holds. The spatially homogeneous periodic solution for system ≤ (m + 1) 2 , and μ < D m hold, which illustrates that the spot patterns will occur. We show these results in the following section.
Remark 3.2 In [28,29], authors considered the existence and stability of positive steady state solutions and investigated asymptotic behavior of time-dependent solutions by utilizing degree theory. They mainly investigated the direction of Hopf bifurcation and the stability of the bifurcation periodic solution for a reaction-diffusion system. In this subsection, except for the research on direction and stability of Hopf bifurcation as stated in [28,29], we also derive the precise conditions on the parameters so that the spatially homogeneous equilibrium solution and the spatially homogeneous periodic solution become Turing unstable. Turing instability in a stable coexistence equilibrium is determined by Tr(J k ) > 0 or Det(J k ) < 0, and Turing instability in a stable limit cycle is determined by Tr(L k ) > 0 or Det(L k ) < 0.
Compared with the former research works, we find that system (3) with the Neumann boundary condition can also induce very interesting dynamical behaviors. For example, (i) system (3) can undergo supercritical Hopf bifurcation under certain parameter conditions. There exists a stable limit cycle when γ is greater than γ H slightly. (ii) Interestingly, we find that Turing instability can occur in the equilibrium and stable limit cycle, respectively.
Compared with the research for reaction-diffusion in two-dimensional (x, y), the theoretical analyses for the stability of equilibria, the type of Hopf bifurcation, and the homogeneous periodic solutions for the corresponding ODE system are valid. However, the theoretical analysis methods used in this paper are not enough for Turing instability in equilibria or homogenous periodic solutions in two-dimensional (x, y). For example, Song and Liu et al. [3,6,7,[31][32][33][34] used the amplitude equation and multi-scale analysis method to analyze Turing instability in equilibrium. In the two-dimensional space (x, y), Turing instability in equilibrium will induce more complex and rich dynamics than that in a onedimensional space x. For example, there may exist a hexagon pattern, a mixed pattern (spot-stripe pattern), and a stripe pattern [3,6,7,33]. Note that, in the one-dimensional space x, we rigidly show how the stable homogeneous periodic solutions can become un-stable. But in the two-dimensional space (x, y), as far as we know, it needs more knowledge and methods to analyze this dynamic behaviors.

Numerical simulations
In this section, we use numerical simulations to illustrate the results in Theorems 2.1-3.2.
Firstly, we show the critical surface for the occurrence of Hopf bifurcation in a parameter space (α, r, γ ), as shown in Fig. 1(a). Particularly, taking α = 2 3 , the critical curve for the occurrence of Hopf bifurcation is shown in Fig. 1(b). From Theorem 2.1, we know that the equilibrium (m * , a * ) is stable if (r, γ ) is above the blue curve, while unstable if (r, γ ) is below the curve.
In the following, let α = 2 3 and r = 1.0901, we have γ H = 6.00074 and (m * , a * ) = (0.2198, 0.7520). In order to show the dynamic behavior near the Hopf bifurcation points γ H , the bifurcation diagram in a three-dimensional space (γ , m, a) is shown in Fig. 2(a). The bifurcation diagram shows the stable equilibrium, the unstable equilibrium, and the limit cycles. The equilibrium is stable if γ < γ H , which is represented by the black curve. The stable equilibrium becomes unstable and stable limit cycles occur if γ > γ H slightly. Here the unstable equilibrium is represented by the red curve, and stable limit cycles are represented by the blue surface. Particularly, let γ = 5.95 < γ H , so condition (H3) holds, we have that the equilibrium (m * , a * ) is asymptotically stable. Set initial val- ues (m 0 , a 0 ) = (0.3198, 0.8520), we find that the orbit starting from (m 0 , a 0 ) converges to (m * , a * ), which shows that (m * , a * ) is asymptotically stable, as shown in Fig. 2(b).
The supercritical Hopf bifurcation occurs for γ = γ H , and a stable limit cycle bifurcates from (m * , a * ) at the right side of γ H . Let γ = 6.1 > γ H , then condition (H4) holds, we have that the equilibrium (m * , a * ) is unstable. Set initial values (m 0 , a 0 ) = (0.2298, 0.7620), we find that the green orbit starting from (m 0 , a 0 ) goes clockwise away from (m * , a * ), which shows that (m * , a * ) is an unstable one, as shown in Fig. 3(a). Set initial value (m 0 , a 0 ) = (0.2498, 0.7820), the blue orbit starting from (m 0 , a 0 ) goes clockwise inside, as shown in Fig. 3(a). Since the blue orbit disjoints the green orbit, there exists a stable limit cycle. We show the stable limit cycle in Fig. 3(b). The above analyses support Theorem 2.1.

Discussion and conclusion
In this paper, we investigate the spatiotemporal dynamics for a diffusive mussel-algae model. Spot and stripe patterns can be generated by a small disturbance on a homogeneous steady state. The results of pattern formation agree with many biological phenomena [11,12,17]. In order to determine the sufficient conditions for Turing instability, we firstly analyze the stability conditions for equilibrium and the critical parameter conditions for Hopf bifurcation. Then, using the center manifold theorem and normal form theory, we show sufficient conditions for the occurrence of Turing instability in equilibrium and limit cycle, respectively. Moreover, we illustrate the theoretic results by a series of numerical simulations. We show abundant dynamical behavior, such as bifurcations, phase trajectory, and patterns. Finally, we find that Turing instability in equilibrium will occur in a stripe pattern and that in limit cycles will occur in a spot pattern.