Dynamic analysis for a class of hydrological model with time delay under fire disturbance

1.   Introduction

Savanna is a mixed ecosystem consists of trees and grasses, and its vegetation dynamic change characteristics are mainly regulated according to the time variability of resource competition, climate, and environmental factors (such as rainfall and fire) ^[1]. The research shows that fire disturbance plays an essential role in the dynamic behavior of vegetation patterns, which is reflected in that fire can maintain the co-existence of trees and grasses and inhibit the extreme ecological phenomenon of single grassland or single forest. Therefore, based on the tree-grass model established by Accatino Francesco and his team, we analyzed the spatial pattern and mechanism of the tree-grass system driven by time delay effect caused by fire and the propagation diffusion effect of vegetation ^[2].

When tree species coexist with grass species, considering some climatic conditions and ecological environment factors, it becomes more difficult to predict the change of vegetation pattern. Based on the data analysis of $854$ observation points in Africa, Sankaran and other scholars concluded that, for the Savanna with annual rainfall less than $650$–$700$ mm, fire disturbance is the main factor affecting the proportion of trees and grasses, and maintaining the coexistence balance of trees and grasses ^[3]. In "balanced competition" models, as we know, the superior competitor limits its own abundance and the inferior competitor can grow ^[4]. Although the water scarcity can limit the abundance of trees in savanna so that grasses can grow. However, in the complex and changeable ecosystem, the competition for limited water resource through niche separation (competitive mechanism) play a limited role in maintaining the co-existence balance of trees and grasses, and is not a necessary condition ^[4]. In addition, for the Savanna with annual rainfall more than $650$–$700$ mm, besides the fire factor, hydrological condition is another major determinant factor affecting on the composition, structure and function of the ecological systems. In terms of vegetation succession stages of Savana, water resources become more critical, especially soil water. It is the connection of hydrological cycle and ecosystem, which has a large impact on the vegetation dynamics ^[5,6]. Therefore, we should clearly take the "soil water" as a state variable, establish a tree-grass model under the coupling of hydrology and fire, and analyze its dynamic characteristics.

The remainder of the article is structured as follows. Section $1$, considering the coupling effect between soil water, tree species and grass species. Also the biological significance of the diffusion coefficient among theses species is discussed, and the influence of the delay effect induced by the fire frequency on the vegetation growth process is clarified. Section $2$, studying the existence and uniqueness of the non-negative solution and the corresponding biological significance for a class of the hydrological model. The distribution of characteristic roots of the system and the Hopf bifurcation caused by time delay are also discussed in detail. Section $3$, takes the time delay $\tau$ as bifurcation parameter, and discusses the direction and stability of Hopf bifurcation by using center manifold and normal form theory ^[7]. Section $4$, the sensitivity analysis of fire frequency and fire intensity is obtained respectively, and the vegetation dynamic simulation are carried out based on the analysis results. Moreover, the variation trend of soil water, trees, and grasses with time-periodic oscillation in the process of ecosystem cycle are obtained. Finally, some conclusions and discussions are given in Section $5$.

2.   The model

In the analysis of the carbon cycle in ecosystems, Accatino Francesco and his team consider not only the coupling effect of tree species and grass species, but also the influence of soil water. Thus, they established a dimensionless hydrological model (2.1) ^[8] based on the tree-grass model, in which the dynamic analysis of the tree-grass system has been discussed in detail in ^[2]. The new part is that the model (2.1) not only describes the differential sensitivity of trees and grasses to fire, but also includes the competition of trees and grasses for water resources. The corresponding vegetation structure of the hydrological model (2.1) is shown in Figure 1, and it shows the vegetation succession and dynamic mechanism under the joint regulation of rainfall and fire in the Savanna ^[8]. The vegetation landscape characteristics it describes are more abundant and complex.

Figure 1 describes the vegetation framework of the ecological cycle system under the coupling of fire and soil water, which shows that the fire intensity depends on the biomass of grass, the moisture content in soil, and relative humidity in the air. It verifies the different feedback of fire disturbance to trees and grasses ^[9]. Among them, rainfall increases soil water content, soil water provides nutrients for trees and grasses, and then promotes the growth of vegetation (inorganic salts, etc.). On the contrary, trees and grasses absorb soil water and then reduce the content of soil water. Trees and grasses compete for survival resources (soil water, inorganic salts, sunlight, and space, etc.) and inhibit each other's growth. Fire utilizes grasses as fuel and directly consumes grasses. Other interference items (grazers, population mortality, etc.) reduce the biomass of grass species, which is equivalent to reducing the accumulation of fire fuel and controlling the spread and intensity of a fire. According to the process of vegetation growth shown in Figure 1, the novelty of the model (2.1) is that rainfall (soil water) is taken into account in the growth function and biomass carrying capacity of trees and grasses.

where $S$ represents the degree of saturation of the profile available water capacity, which is dimensionless and satisfies the condition $S\in\left[0, \; 1\right]$. In particular, $S = 0$ which corresponds to completely dry soil, and $S = 1$ which corresponds to saturated soil. $T$ and $G$ represent the fraction of the area occupied by trees and grasses respectively, they are dimensionless and satisfy the condition $T, \; G\; \in \left[0, \; 1\right]$. In particular, $T = 0\; \left(G = 0\right)$ means there are no trees (grasses) growing in the site, and $T = 1\; \left(G = 1\right)$ means that the area is completely covered by trees (grasses). The first equation of the hydrological model (2.1), $\varepsilon S\left(1-T-G\right)$ represents the evaporation of exposed soil, which depends on the proportion of exposed soil $\left(1-T-G\right)$ and the dirt soil moisture $S$. $\tau_{T}ST$ and $\tau_{G}SG$ represent the water uptake by the trees and grasses, respectively. In addition, $z$ is the root zone depth, $n$ is the porosity (fractional pore volume) and $w_{1} = z\times n$ is the pore space in the volume. $\gamma_{T}$ and $\gamma_{G}$ are the maximum colonization rates, and $\frac{p}{w_{1}}$ is the deep percolation. The ecological significance of the growth term, the removal term and parameters $\delta_{T}$, $\delta_{GO}$, $\delta_{F}$, $f$ in the second equation and the third equation of the hydrological model (2.1) are explained in detail in ^[2].

For convenient of the mathematical derivation, let $\frac{p}{w_{1}} = w$, $\tau_{T} = a_{1}$, $\tau_{G} = a_{2}$, $\delta_{T} = b_{1}$, $\delta_{GO} = b_{2}$, $\gamma_{T} = c_{1}$, $\gamma_{G} = c_{2}$, $\delta_{F} = \delta$, $S = u_{1}(t)$, $T = u_{2}(t)$ and $G = u_{3}(t)$. Considering the ecological significance, the environmental parameters $w$-$f$ (including $w$, $a_{2}$, $b_{1}$, $b_{2}$, $c_{1}$, $c_{2}$, $\delta$ and $f$) defined by rainfall and fire are positive. The hydrological model (2.1) becomes

According to the dynamic changes of grass biomass in the model (2.2), for the main fuel of fire is grass biomass, then the frequent occurrence of the fire directly leads to the reduction of grass biomass $u_{3}(t)$, which is reflected in item $u_{3}(t)f$. According to the dynamic changes of tree biomass in the model (2.2), for the combined effect of fire frequency and fire intensity, the reduction of tree biomass $u_{2}(t)$ caused by fire disturbance is closely related to grass biomass $u_{3}(t)$, which is reflected in the item $u_{2}(t)u_{3}(t)\delta f$. It further shows that the frequency and intensity of fire are positively correlated with grass biomass, that is, grass biomass determines the frequency and intensity of fire, which is reflected in $u_{3}(t)\delta f$.

For the hydrological model (2.2) under the coupling of the fire and hydrology, tree species and grass species may diffuse in the direction of regions rich in resources, or may diffuse in areas with few competitors because of the asymmetric competition between species, or may diffuse with the influence of external environmental factors such as water flow. Based on the above complex ecological process, it isn't easy to study the influence of diffusion on the evolution trend of vegetation pattern by using the bifurcation theory of reaction-diffusion equation. It is worth noting that the delay effect caused by fire frequency is one of the key parameters affecting the dynamic change of vegetation landscape pattern in Savanna. The vegetation density distribution before the occurring time of fire directly affects the next state of vegetation growth, that is, a series of consequences caused by the time delay "$\tau$" of fire frequency can not be ignored. Based on the above considerations, this paper introduces the effect of time delay of fire frequency into the model (2.2). It discusses the potential dynamic mechanism in soil water-tree-grass growth without considering the effects of propagation and diffusion of vegetation. It is found that the response of grass species to fire feedback is more significant than that of tree species ^[1], so it is more representative to discuss the effect of time delay $\tau$ on grass species. Next, $\tau$ is introduced into the reduction item $u_{3}(t)f$ caused by fire to grass species $u_{3}(t)$ in the model (2.2), then we have established a class of hydrological model with time delay as follows

As we all know, the study of dynamical systems not only includes the discussion of stability, attractiveness and persistence, but also includes many dynamic behaviors, such as periodic phenomena, bifurcation and chaos. Because the bifurcation phenomenon corresponding to high codimension bifurcation is beneficial to the study of the dynamic characteristics of system, that is to say, it is very meaningful to study the direction and stability of the Hopf bifurcation periodic orbit of the model, but at the same time, this research is also a difficult and challenging task. In order to obtain abundant bifurcation phenomena of the hydrological model (2.3) with time delay, the nature (bifurcation direction and stability) of Hopf bifurcation of the system under the framework of Faria and Magalhaes' norms-based method ^[7] will be considered in the following analysis.

3.   Dynamic analysis of the hydrological model with time delay

In this section, we analyze the dynamic behavior near the non-negative solution of the hydrological model (2.3) from two aspects: one is not considers the delay effect of fire frequency ($\tau = 0$), and the other is considers the delay effect of fire frequency $\tau\neq0$. For convenience, let

When $\tau = 0$, the internal equilibrium points of the model (2.3) are obtained by calculating equations $g_{1}(t) = 0$, $g_{2}(t) = 0$, and $g_{3}(t) = 0$, including $E_{1} = (s_{1}, 0, 0)$, $E_{2} = (s_{2}, t_{2}, 0)$, $E_{3} = (s_{3}, 0, g_{3})$ and $E_{0} = (s_{0}, t_{0}, g_{0})$, where $c = {c_1} + {c_2}$ and

There are always three semi-trivial solutions $E_{1} = (s_{1}, 0, 0)$, $E_{2} = (s_{2}, t_{2}, 0)$, $E_{3} = (s_{3}, 0, g_{3})$ and one nontrivial solution $E_{0} = (s_{0}, t_{0}, g_{0})$ in the system, and their biological significance are shown in Table 1.

Obviously, from the perspective of environmental protection and ecological sustainable development, people do not want to see the extreme natural phenomena represented by the equilibrium points of $E_{1}$, $E_{2}$ and $E_{3}$, but hope that species can coexist and grow together. Clearly $E_{0} = (s_0, t_0, g_0)$ is feasible and has biological significance if and only if $s_0 > 0$, $t_0 > 0$ and $g_0 > 0$ are all true together. To investigate the stability of the equilibrium $E_{0}$, we compute the eigenvalues of the Jacobian matrix at this equilibrium based on the dynamic theory. In what follows, we first consider the sufficient conditions of existence of positive equilibrium point $E_{0}$, that is

When the conditions $\left(H_1\right)$ holds, the model (2.3) has a unique positive equilibrium solution $E_{0}$.

3.1.   Stability analysis

This section will study the dynamic properties of the hydrological model (2.3) near the coexistence equilibrium point $E_{0}$. The key of theoretical research is to calculate the distribution of eigenvalues $\lambda$ corresponding to the Jacobian matrix of the system in the equilibrium state. Let ${\tilde u_1}(t) = {u_1}(t) - {s_0}, \; {\tilde u_2}(t) = {u_2}(t) - {t_0}, \; {\tilde u_3}(t) = {u_3}(t) - {g_0}$, drop out $"\backsim"$, and the linearization equation corresponding to the model (2.3) is

where

Next, calculate the eigenvalue of Jacobian matrix at $E_{0}$. The eigenvalue satisfies the characteristic equation

where

When $\tau = 0$, the delay effect of fire frequency on the soli water-tree-grass dynamic balance mechanism of the model (2.3) is not considered. With the condition $\left({\bf{H}}_1\right)$, the model satisfies the existence, uniqueness and continuous dependence on initial condition of solutions. Moreover, there exists a positively invariant box ${\bf{B}}$, and such that all solutions $\{ {u_1}(t), \; {u_2}(t), \; {u_3}(t)\}$ with non-negative initial conditions approach ${\bf{B}}$ as $\tau\rightarrow 1$ ($\tau$ is the non-negative constant, and $t\in [-\tau, 0]$). Namely

The characteristic equation of hydrological model (2.3) with $\tau = 0$ at $E_{0}$ is

where ${p_1} = {q_1} + {q_3}, {p_2} = {q_2} + {q_4}, {p_3} = {q_5} + {q_6}$. According to the Routh-Hurwitz stability criterion, when ${p_1} > 0, \; {p_2} > 0, \; {p_3} > 0$ and ${p_1}{p_2} - {p_3} > 0$, all eigenvalues of the equation $DET\left({\lambda, 0} \right) = 0$ have negative real parts, that is, the system is locally asymptotically stable at $E_{0}$. Let

Combined with the above analysis, the following theorem exists.

Lemma 3.1. when $\tau = 0$, if conditions $\left({\bf{H}}_1\right)$ and $\left({\bf{H}}_2\right)$ hold, that the hydrological model (2.3) is locally asymptotically stable at uniform state $E_{0}$.

From the perspective of ecology, when conditions $\left({\bf{H}}_1\right)$ and $\left({\bf{H}}_2\right)$ hold, the value range of system parameter and fire frequency $w$-$f$ are more constrained, in which, the determined parameter range is called threshold interval. $E_{0}$ is asymptotically stable inside the threshold range, which means that the population number of trees and grasses gradually increases to $E_{0}$ over time, tending to a constant and fixed value. $E_{0}$ is unstable outside the threshold range, which means that there is a large amplitude disturbance in the dynamic changing of population density of trees and grasses.

3.2.   Hopf bifurcation

When $\tau\neq0$, the delay effect of fire frequency on the soli water-tree-grass dynamic balance mechanism of the model (2.3) is considered. According to the stability theorem ^[10], the necessary and sufficient conditions for $E_{0}$ to be asymptotically stable for all $\tau \geq0$ are:

$1)$ The real parts of all the roots of the characteristic equation $DET(\lambda, \; 0) = 0$ are negative.

$2)$ For all real $\kappa$ and any $\tau > 0$, $DET(i \kappa, \tau) \ne 0$, where $i^2 = -1$.

Lemma 3.1 obviously verifies the case $1)$. In the following, the condition of the case $2)$ will be analyzed. When $\kappa = 0$, we have $DET(0, \tau) = {p_3} > 0$. When $\kappa > 0$, let $\lambda = \gamma + i \kappa$ be the characteristic root of Eq (3.2), then after separating the real part and the imaginary part, we get a transcendental equations as follows

among them, the parameters $\gamma$ and $\kappa$ in the eigenvalue $\lambda$ are functions of time delay $\tau$. When $\gamma = 0$ and $\kappa\neq0$, the model (3.4) becomes

After simplification, the characteristic equation becomes

where ${\omega ^2} = s,$ ${m_1} = q_1^2 - 2{q_2} - q_3^2$, ${m_2} = q_2^2 + 2{q_3}{q_5} - q_4^2 - 2{q_1}{q_6}$, ${m_3} = q_6^2 - q_5^2$. By using the theory of Veda, the roots $s_{1}$, $s_{2}$, $s_{3}$ of the Eq (3.6) should satisfy the following requirements.

combined with Routh-Hurwitz, the sufficient conditions for the non-existence of a real numbers satisfying $D\left(i\kappa, \tau\right) = 0$ can be expressed by the sufficient conditions for non-existence of a real $s$ as

then the Eq (3.6) has a unique positive real root $s_0 = \kappa_0^2$. Thus, we have the following conclusion:

Lemma 3.2. If the condition $({\bf{H}}_{3})$ holds, then $E_{0}$ of the hydrologic model (2.3) is asymptotically stable for all $\tau > {\tau _0}$.

Let the time delay $\tau$ be the bifurcation parameter, and $\tilde \tau$ be the critical value for the occurrence of Hopf bifurcation. Next, the transversality condition of Hopf bifurcation will be verified, such as $d \gamma /d\tau {{{\rm{|}}}_{\tau = \tilde \tau }} \ne 0$. First, we are going to take derivative with respect to $\tau$ of the model (3.4), and let $\tau = \tilde \tau$, $\gamma = 0$, then we have

in which

By solving the Eq (3.7), then have

in which, the signs of $d\gamma /d\tau$ is as same as $LU + MV$ ^[11], and $LU + MV$ is equivalent to the following equation

If the conditions $({\bf{H}}_{1})$–$({\bf{H}}_{3})$ hold, obviously there have $LU + MV\neq 0$, then

therefore, the transversality condition of Hopf bifurcation holds. Suppose $\pm {i}{\kappa _0}$ are a pair of complex conjugate numbers of the Eq (3.7), then the critical value $\tilde \tau$ satisfies the following equation

The minimum value of time delay is recorded as ${\tau _0}$, as follows

Remark 3.1. For convenient to get Hopf bifurcation analysis of model (2.3), the precondition of the main parameter selection is to guarantee that the time delay $\tau _0 > 0$.

According to the Hopf bifurcation theory ^[12], we have the following result on the stability of the hydrologic model (2.3).

Theorem 3.1. If the conditions $({\bf{H}}_{1})$–$({\bf{H}}_{3})$ hold, the expressions for ${\kappa _0}$ and ${\tau _0}$ are (3.6) and (3.8), respectively. Then

Case i: When $\tau \in [0, \; {\tau _0})$, then the soil water-tree-grass coexistence equilibrium point $E_{0}$ of the hydrologic model (2.3) is asymptotically stable, while when $\tau > {\tau _0}$, $E_{0}$ is unstable.

Case ii: When $\tau = \tilde \tau,$ $j = 0, 1, 2\cdots$ the hydrologic model (2.3) undergoes a sequence of Hopf bifurcations at the soil water-tree-grass coexistence equilibrium point $E_{0}$.

4.   Stability of Hopf bifurcation periodic solutions

In this section, the stability of Hopf bifurcation periodic solutions of the hydrologic model (2.3) will be analyzed by using the central manifold theorem and the normal form theory proposed by Faria and Magalhaes ^[7]. In addition, when $\tau = {\tau _0}$, the characteristic equation (3.2) has a pair of pure imaginary roots $\pm$ $i{\kappa _0}$, and the other roots have negative real parts. By time-scaling $t \to t/\tau$, then the hydrologic model (2.3) is rewritten as a functional differential equation in space ${\bf{C}} = \left(\left[ - 1, \; 0\right], \; \mathbb{C}^2\right)$.

Let $\mu$ be the bifurcation parameter and $\mu = \tau-\tau_{0}$, then the linearization equation of the model (4.1) near the origin is

Let $\eta (\theta) = \mathbb{A} \delta (\theta) + \mathbb{B} \delta (\theta + 1)$, where

the higher order terms of model (4.1) are

where

Suppose $U(t) = {(\begin{array}{*{20}{c}} {{u_1}(t)}, & {{u_2}(t)}, & {{u_3}(t)} \end{array})^T}$, then model (4.2) can be written as the vector form in phase space ${\bf{C}}$:

where $\mathfrak{L}$ is linear bounded operator. By the Riesz representation theorem and Riemann-Stieltjes integral, then

In which $\eta (\theta)\; (\theta \in [ - {{\rm{1}}}, \; 0])$ is $2\times 2$ matrix function of bounded variation, and

According to the semigroup theory of linear operators ^[13], the infinitesimal generator ${A^0}$ of the solution semigroup ${C^0}$ of model (4.2) is obtained, which satisfies equation ${A^0}\varphi {{\rm{(}}}\theta {{\rm{)}}} = \dot \varphi {{\rm{(}}}\theta {{\rm{)}}} + {U_{{\rm{0}}}}[\mathfrak{L}\varphi {{\rm{(}}}\theta {{\rm{)}}} - \dot \varphi (0)]$, and the domain of ${A^0}$ is defined as follows

In addition, the bilinear form on ${{{\rm{C}}}^*} \times {{\rm{C}}}$ ($*$ represents the adjoint matrix) is

From the discussions in Section $2.2$, we suppose that $\mathfrak{L}$ has a pair of simple purely imaginary eigenvalues $\pm$ $i\kappa \; (\kappa > 0)$, and all other roots have negative real parts. Let $\Lambda = \{ i\kappa, \; - i\kappa \}$, and the two-dimensional matrix ${{\rm{P}}}$ is the generalized eigenspace associated with $\Lambda$, and note that ${{{\rm{P}}}^*}$ is the space adjoint with ${{{\rm{P}}}}$. Then ${{\rm{C}}}$ can be decomposed as ${{\rm{C}}} = {{\rm{P}}} \oplus {{\rm{Q}}}$, where ${{\rm{Q}}} = {{\rm{\{ }}}\varphi \in {{\rm{C}}}:\; \; < \psi, \; \varphi > \; = 0, \; \forall \psi \in {{{\rm{P}}}^*}{{\rm{\}} }}$. Suppose $\Phi$ and $\Psi$ as the bases of ${{{\rm{P}}}^*}$ and ${{\rm{P}}}$ respectively, and such that $< \Psi, \; \Phi > \; = {{\rm{I}}}$, $\dot \Phi = \Phi {{\rm{J}}}$, $\dot \Psi = - {{\rm{J}}}\Psi$, where ${\rm{I}}$ is $2\times 2$ identity matrix, ${\rm{J}}$ is a diagonal matrix, written as $diag(i{\kappa _0}{\tau _0}, \; - i{\kappa _0}{\tau _0})$. Next, we should enlarge the phase space ${\rm{C}}$ to the space ${\rm{BC}}$

On the phase space ${\rm{BC}}$, the model (4.3) is equivalent to the following abstract ordinary differential system

where

in the form of bilinear integral, the infinitesimal generator ${A^0}$ and its adjoint operator ${A^0_*}$ satisfy

Considering the complex coordinates, and suppose that $\left({\varphi, \; \bar \varphi } \right)$ and ${\left({\bar \psi, \; \psi } \right)^T}$ are the eigenvectors of the operator ${A^0}$ and ${A^0_*}$ corresponding to eigenvalues $i{\kappa _0}{\tau _0}, - i{\kappa _0}{\tau _0}$, respectively.

Lemma 4.1. Let $\Phi$ = ${\left({\varphi, \; \bar \varphi } \right)}$ is a set of bases for ${\rm{P}}$, $\Psi$ = ${\left({\bar \psi, \; \psi } \right)^T}$ is a set of bases for ${{{\rm{P}}}^{{\rm{*}}}}$, and such that $\dot \Phi = \Phi {{\rm{J}}}, \; \dot \Psi = - {{\rm{J}}}\Psi,$ where

Proof. Suppose $\varphi$ = ${\left({1, \; {\kern 1pt} {\eta _1}, \; {\eta _2}} \right)^T}{e^{i{\kappa _0}{\tau _0}\theta }}$ is an eigenvector of ${A^0}$ corresponding to $i{\kappa _0}{\tau _0}$, and satisfy ${A^0}\varphi \left(\theta \right) = i{\kappa _0}{\tau _0}\varphi \left(\theta \right)$. From the definition of operator ${A^0}$, there are

after solving the above equation, we get

let $\Phi$ = ${\left({\varphi, \; \bar \varphi } \right)}$ such that $\dot \Phi = \Phi {\rm{J}}$.

Similarly, suppose $\psi$ = ${D}\left({1, \; {\kern 1pt} {\sigma _1}, \; {\sigma _2}} \right){e^{i{\kappa _0}{\tau _0}s}}$ is an eigenvector of $A_{{\rm{*}}}^0$ corresponding to $- i{\kappa _0}{\tau _0}$, and satisfy $A_{{\rm{*}}}^0\psi \left(s \right) = - i{\kappa _0}{\tau _0}\varphi \left(\theta \right)$. From the definition of operator $A_{{\rm{*}}}^0$, there are

after solving the above equation, we get

let $\Psi$ = ${\left({\bar \psi, \; \psi } \right)^T}$ such that $\dot \Psi = - {{\rm{J}}}\Psi$.

From $\left\langle { \Psi \Phi } \right\rangle = {{{\rm{I}}}_{2 \times 2}}$, it is easy to check $\left\langle { {\bar \psi, } \varphi } \right\rangle = \left\langle { {\psi, } {\bar \varphi } } \right\rangle = 1, \; \left\langle { {\psi, } \varphi } \right\rangle = \left\langle { {\bar \psi, } {\bar \varphi } } \right\rangle = 0$. To ensure $\left\langle { {\bar \psi, } \varphi } \right\rangle = 1$ holds, we first calculate the coefficient $D, \bar D$. Since

then get $\bar D{{\rm{ = }}}{1 / {\left({1 + {\eta _1}{{\bar \sigma }_1} + {\eta _2}{{\bar \sigma }_2} + {\tau _0}{\eta _2}{{\bar \sigma }_2}{b_{33}}{e^{ - i{\kappa _0}{\tau _0}}}} \right)}},$ the Lemma 3 has been proved.

Define $U\left(t \right) = \Phi x + y$, $x \in {\mathbb{C}^2}, \; y \in {Q^1}{{\rm{ = }}}\left\{ {\varphi \in Q:\; \; \dot \varphi \in C} \right\}$, then there is the following decomposition formula

Let

then the model (4.4) is decomposed into

it can also be written as follows

where

and

The model (4.6) can be transformed as the following normal form

in order to obtain the specific expression of normal form Eq (4.7), $g_2^1\left({x, 0, \mu } \right)$ and $g_3^1\left({x, 0, \mu } \right)$ need to be calculated, in which, $g_j^1\left({x, 0, \mu } \right)$ is a $j$-order $\left({j = 2, 3 \cdots } \right)$ homogeneous polynomial about $\left({x, \mu } \right)$.

Define operator $M_j^1$

in addition, define

and $e_{1} = \left(1, \; 0\right)^{T}$, $e_{2} = \left(0, \; 1\right)^{T}$. Then

First, we should calculate $g_2^1\left({x, 0, \mu } \right)$. From the Eq (4.5), we have

where

Since the second-order terms of $\left(x, \mu\right)$ on the central manifold can be obtained by

Second, we should calculate $g_3^1\left({x, 0, \mu } \right)$, and $g_3^1\left({x, 0, \mu } \right) \in Ker\left({M_3^1} \right)$. The higher-order terms $o\left({\left| x \right|, \; \mu } \right)$ of parameter $\mu$ is independent of the general Hopf bifurcation. Therefore, only the coefficients of the term without parameter $\mu$ in space $Ker\left({M_3^1} \right)$ need to be calculated, and simplified as $Sp = \left\{ {x_1^2{x_2}{e_1}, \; {x_1}x_2^2{e_1}} \right\}$, then

where $\tilde f_3^1\left({x, 0, 0} \right)$ represents the third-order term of the system, which is obtained by calculating the second-order term of the normal form. The results are as follows

in which, $U_2^1\left({x, \mu } \right)$ is a solution of $M_2^1U_2^1\left({x, \mu } \right) = f_2^1\left({x, 0, 0} \right)$, and $h = h\left(x \right)\left(\theta \right) = U_2^2\left({x, 0} \right)$ is the coefficient of the second-order homogeneous polynomial with respect to $\left({{x_1}, {x_2}, \mu } \right)$.

According to (4.1)–(4.5), then

then

where ${A_3}{{\rm{ = }}} - {\tau _0}\bar D\left[ {{c_1}\eta _1^2{{\bar \sigma }_1} + {c_2}\eta _2^2{{\bar \sigma }_2} + \left({c{{\bar \eta }_1} + 2{c_2}{{\bar \eta }_2}} \right){\eta _2}{{\bar \sigma }_2} + 2{c_1}{{\left| {{\eta _1}} \right|}^2}{{\bar \sigma }_1} + 2c{\mathop{\rm Re}\nolimits} \left[ {{\eta _2}} \right]{\eta _1}{{\bar \sigma }_2}} \right]$. It is easy to get $g_2^1\left({x, 0, 0} \right){{\rm{ = }}}0$, then $\left({{D_x}U_2^1} \right)g_2^1\left({x, 0, 0} \right){{\rm{ = }}}0$.

Third, we should calculate $U_2^1\left({x, 0} \right)$ and $U_2^2\left({x, 0} \right)$ to get the normal form. From (4.8), we have

according to the definition of $M_2^1$, then $M_2^1U_2^1\left({x, 0} \right)$ is written as the following partial differential equation

From (4.11), we get

then

where ${A_4} = \frac{i}{{2{\kappa _0}}}\left({{b_{11}}{b_{20}} - \frac{1}{3}{{\left| {{b_{02}}} \right|}^2} - 2{{\left| {{b_{11}}} \right|}^2}} \right)$.

Finally, we should calculate $\left[ {{D_y}f_2^1\left({x, 0, 0} \right)} \right]U_2^2\left({x, 0} \right)$. Define $h = h\left(x \right)\left(\theta \right) = U_2^2\left({x, 0} \right)$, and $h$ is a second-order homogeneous polynomial of $\left({{x_1}, {x_2}, \mu } \right)$, let

where $\left\{ {{h_{11}}\left(\theta \right), {h_{20}}\left(\theta \right), {h_{02}}\left(\theta \right)} \right\} \in {Q^1}$. From (4.5), we have

then

where

Since

where

Combine with the Eqs (4.10), (4.12) and (4.13), we have

where ${B_1} = {A_3} + {A_4} + {A_5}$.

In the following, we should calculate the unknown term ${h_{11}}\left(\theta \right), {h_{20}}\left(\theta \right)$ ($\theta \in \left[ { - 1, 0} \right]$) of the coefficient ${A_5}$. According to ^[7], $h$ is uniquely determined by $\left({M_2^2h} \right)\left({x, \mu } \right) = f_2^2\left({x, 0, 0} \right)$, and is equivalent to

from the definition of $A_{Q^1}$ and $\pi$, we get

in which, $\dot h$ is the derivative of function $h\left(\theta \right)$ with respect to $\theta$. Let $A_{11} = \left({\begin{array}{ccc} A_{11}^{1}, &A_{11}^{2}, &A_{11}^{3}\end{array}}\right)^{T}$ and $A_{20} = \left({\begin{array}{ccc} A_{20}^{1}, &A_{20}^{2}, &A_{20}^{3}\end{array}}\right)^{T}$, then

and ${A_{ij}} \in {\mathbb{C}^3}, 0 \le i, j \le 2$. Comparing the coefficients of $x_1^2, x_2^2$ and ${x_1}{x_2}$, there is ${\bar h_{20}} = {h_{02}}$, and ${h_{20}}, \; {h_{11}}$ satisfy the following differential equations, respectively.

solve Eqs (4.16) and (4.17) respectively, and finally get

then

where ${c^1}, \; {c^2} \in {\mathbb{C}^3}$, and satisfy the following equations

then we have

for the detailed derivation of Eqs (4.18)–(4.20), refer to ^[14], and coefficients are obtained after easy but long computation as follows

in which

Combine (4.9) and (4.14), the model (4.7) can be reduced to the following form

since ${\bar x_1} = {x_2}$. Through the change of ${x_1} = {w_1} - i{w_2}, {x_2} = {w_1} + i{w_2}$, and let ${w_1} = \nu {\mathop{\rm Cos}\nolimits} \varsigma$, ${w_2} = \nu {\mathop{\rm Sin}\nolimits} \varsigma$, then (4.21) becomes

where ${\chi _1}$ = ${\mathop{\rm Re}\nolimits} [{A_1}], {\chi _2} = {\mathop{\rm Re}\nolimits} [{B_1}], {\chi _3} = {\mathop{\rm Im}\nolimits} [{A_1}]$. Following the conclusion of Hale ^[15], the sign of $\chi_{1}\chi_{2}$ determines the direction of the Hopf bifurcation. And the sign of $\chi_{2}$ determines the stability of the nontrivial periodic solution generated by Hopf bifurcation. When $\chi_{1}\chi_{2} < 0$, Hopf bifurcation is supercritical, and when $\chi_{1}\chi_{2} > 0$, Hopf bifurcation is subcritical. In addition, the bifurcation periodic solution on the central manifold is stable when $\chi_{2} < 0$ and unstable when $\chi_{2} > 0$. The central manifold theory shows that the stability of the periodic solution of the model (2.3) in the whole phase space is consistent with the normal form projected on the central manifold. Therefore, the following theorem exists.

Theorem 4.1. When $\mu = 0$, the stability of the bifurcation periodic solution of the model (2.3) at $E_{0}$ is as follows

Case i: If $\chi_{1} < 0, \; \chi_{2} < 0$, the Hopf bifurcation of the model (2.3) is subcritical, and the bifurcation periodic orbit is asymptotically stable.

Case ii: If $\chi_{1} < 0, \; \chi_{2} > 0$, the Hopf bifurcation of the model (2.3) is supercritical, and the bifurcation periodic orbit is unstable.

Case iii: If $\chi_{1} > 0, \; \chi_{2} < 0$ the Hopf bifurcation of the model (2.3) is supercritical, and the bifurcation periodic orbit is asymptotically stable.

Case iv: If $\chi_{1} > 0, \; \chi_{2} > 0$ the Hopf bifurcation of the model (2.3) is subcritical, and the bifurcation periodic orbit is unstable.

5.   Dynamic experimental and simulation

This section mainly analyzes two problems. One is to discuss the effect of fire frequency and fire intensity on each variable in a hydrological model (2.3) without time delay, including soil water, tree species, and grass species. After that, we conclude different feedback effects of fire frequency and fire intensity on the ecological balance of vegetation through sensitivity analysis. The other is to analyze the corresponding numerical simulation experiments on the potential ecological processes in the hydrological model (2.3) with time delay, and verify the accuracy of the theoretical analysis.

5.1.   Sensitivity analysis

The model (2.3) contains many parameters, among which the fire frequency $f$ and fire intensity $\delta$ have a great effect on the soil water-tree-grass co-existence mechanism. In this section, without considering the time delay, the sensitivity of fire frequency and fire intensity to soil water, tree species, and grass species is carried out by using the Runge-Kutta method (RK-4) ^[16]. Assuming ${u_i}\left(t \right) = {u_i}$, ${S_j}\left(t \right) = {S_j}$ ($i = 1, 2, 3,$ $j = 1, 2, \cdots, 6$), and the sensitivity equation of the model (2.3) includes $9$ equations, as follows

where

Let the time step $\Delta t = 2000$, fire frequency $f = 1.35$, the values of other parameters are shown in Table 2, and the co-existence equilibrium point ${E_0} = \left(0.6114, 0.6467, 0.0067\right)$. Carry out simulation experiment on the model (5.1) by using the RK-4 method, then the sensitivity and relative sensitivity analysis for fire frequency $f$ and fire intensity $\delta$ to soil water $u_{1}$, trees $u_{2}$ and grasses $u_{3}$ are obtained. After that, the conclusion can reflect the effect of fire frequency and fire intensity on each variable (soil water, trees and grasses) in the hydrological model (2.3) without time delay, as shown in Figures 2–5.

In Figures 2a, 3a, 4a and 5a, the sensitivity ${S_1}, \; {S_4}$ and relative sensitivity $f\; {S_1}/{u_1}$, $\delta\; {S_4}/{u_1}$ of $f$ and $\delta$ to soil water $u_{1}$ are all greater than zero, which indicates that the soil water storage would increase with the increasing of $f$ and $\delta$. When $t = 5$, we have ${S_1} = 0.04808 > 0$, $f\; {S_1}/{u_1} = 0.09661 > 0$ and ${S_4} = 0.07455 > 0$, $\delta \; {S_4} /{u_1} = 0.03858 > 0$. In other words, when $f$ increases by $10{{\rm{\% }}}$, the soil water storage $u_{1}$ will increase by $0.9661\%$, and when $\delta$ increases by $10{{\rm{\% }}}$, the soil water storage $u_{1}$ will increase by $0.3858\%$. At the same time, the increase of soil water storage promoted by fire frequency $f$ is greater than that promoted by fire intensity $\delta$ ($0.9661\% > 0.3858\%$), indicating that in the hydrological model (2.3), $f$ has more significant effect on the dynamic change of soil water than $\delta$. From a biological viewpoint, compared with grass species, tree species have a greater demand for soil water, and the competition for nutrients is in an advantageous position. In addition, frequent fires will eventually lead to the extinction of trees. Therefore, soil water will increase steadily due to the decline of supply-demand, and the continuous supply of rainfall.

In Figures 2b, 3b, 4b and 5b, the sensitivity $S_{2}, \; S_{5}$ and relative sensitivity $f{{\rm{ }}}{S_2}/{u_2}$, $\delta {{\rm{ }}}{S_5}/{u_2}$ of $f$ and $\delta$ to trees $u_{2}$ are all less than zero, which indicates that the tree species would decrease with the increasing of $f$ and $\delta$. When $t = 5$, we have ${S_2} = -0.1778 < 0$, $f\; {S_2}/{u_2} = -0.66 < 0$ and ${S_5} = -0.09075 < 0$, $\delta \; {S_5}/{u_2} = -0.08677 < 0$. In other words, when $f$ increases by $10{{\rm{\% }}}$, the trees $u_{2}$ will decrease by $6.6\%$, and when $\delta$ increases by $10{{\rm{\% }}}$, the trees $u_{2}$ will decrease by $0.8677\%$. At the same time, the decrease of tree species caused by fire frequency $f$ is more significant than that caused by fire intensity $\delta$ ($6.6\% > 0.8677\%$), indicating that in the hydrological model (2.3), $f$ has a greater effect on the dynamic change of trees than $\delta$. From a biological viewpoint, the increase of fire frequency and fire intensity are all bad for the growth of trees and grasses. But the growth cycle and population recovery process of tree species take longer than that of grasses, resulting in the slow succession of tree species. Therefore, the increase of fire frequency will destroy the growth of trees and then lead to the extinction of trees finally.

In Figures 2c and 4c, the sensitivity $S_{3}$ and relative sensitivity $f{{\rm{ }}}{S_3}/{u_3}$ of $f$ to grasses $u_{3}$ show a change trend of less than zero at first and then greater than zero. The above parametric change interval can be divided into two regions: there were $S_{3}\leq 0$ and $f{{\rm{ }}}{S_3}/{u_3}\leq 0$ for $t\in\left[0, \; 2.12\right]$, and there were $S_{3} > 0$ and $f{{\rm{ }}}{S_3}/{u_3} > 0$ for $t > 2.12$. Taken $t = 2$ in the first region, then the relative sensitivity $f{{\rm{ }}}{S_3}/{u_3} = -0.0156$, which shows that when $f$ increases by $10\%$, $u_{3}$ decreases by $0.156$. And taken $t = 2.5$ in the second region, then the relative sensitivity $f{{\rm{ }}}{S_3}/{u_3} = 0.05142$, which shows that when $f$ increases by $10\%$, $u_{3}$ increases by $0.5142$. From a biological viewpoint, under mild and moderate interference of fire disturbance, the density distribution of dead trees is less than $70\%$, the density distribution of living trees is more than $30\%$, and the consumption of fuel composed of shrub layer, herb layer and litter layer is less than $50\%$ ^[19]. The live trees have a wide distribution in this time and is in a dominant position to compete with grass species for nutrient resources. Therefore, when grass species cannot obtain sufficient nutrition, the density distribution of grass biomass will show a decreasing trend, which corresponds to the theoretical analysis in the case where the values of sensitivity and relative sensitivity of fire frequency $f$ to grasses $u_{3}$ is less than zero. In addition, under severe interference of fire disturbance, the density distribution of dead trees is greater than $70\%$ or more, the density distribution of live trees is less than $30\%$, and the consumption of fuel composed of shrub layer, herb layer and litter layer is almost burned and depleted ^[19]. The distribution of live trees is very scattered, so the tree species no longer have an advantage in competing with grass species for nutrient resources. And because of the slow growth cycle of trees, the frequent occurrence of fire will bring the growth of new trees to a near standstill, which provides vitality for grass recovery quickly and stably. Then the above conclusions are consistent with the theoretical analysis in the case where the values of sensitivity and relative sensitivity of fire frequency $f$ to grasses $u_{3}$ are greater than zero.

In Figures 3c and 5c, the sensitivity $S_{6}$ and relative sensitivity $\delta {{\rm{ }}}{S_6}/{u_3}$ of $\delta$ to grasses $u_{3}$ show a change trend of greater than zero, which indicates that the tree species would increase with the increasing of $\delta$. When $t = 5$, we have $\delta {{\rm{ }}}{S_6}/{u_3} = 0.08665$, that is, when $\delta$ increases by $10\%$, the grasses $u_{3}$ will increase by $0.8665$, which seems to be contrary to our common sense. However, from a biological viewpoint, the herbaceous plants are mostly $1$–$3$ years old and complete their full life span in nearly one growing season. Under the disturbance of proper fire intensity, burning weeds and dead leaves will give vegetation better space to grow. And then, herbs contain inorganic salts and other nutrients absorbed from the soil. Burned herbs turn into ash and return to the soil, which can bring sufficient fertilizer to the soil, replenish the inorganic content of the soil, and loosen the soil. In addition, the proper intensity of fire disturbance can also burn the overwintering eggs in dead grasses and increase fertilizer for new grasses' growth. Therefore, the theoretical analysis of the sensitivity and relative sensitivity of fire intensity to grasses are consistent with the growth pattern of vegetation under the actual background.

In summary, the effect of fire frequency on vegetation ecological balance in the hydrological model (2.3) is more extensive relative to fire intensity. When the time delay effect of fire frequency is not considered, model (2.3) can reasonably reflect the feedback effects of fire frequency and intensity on soil water-tree-grass growth trends.

5.2.   Numerical simulation

In this section, the dynamic properties of the hydrological model under fire disturbance are verified on the basis of sensitivity analysis, and the corresponding numerical simulations are carried out by using MATLAB for the model (2.3) with $\tau\neq 0$. The values of the parameters see Table 2, take $f = 1.35$ and then we have $w = p/w_{1} = 35$, the unique positive equilibrium is $E_{0} = (0.6114, \; 0.6467, \; 0.0067)$. In addition, we have $p_{1} = 160.3177$, $p_{2} = 2220.2434$, $p_{1}p_{2}-p_{3} = 355944$ and $m_{1}m_{2} = -1968073$, that is, the conditions $({\bf{H}}_{2})$ and $({\bf{H}}_{3})$ are satisfied, and then the Theorems 1, 2 hold. Moreover, the cross section condition $(LU+MV)/(M^{2}+L^{2}) = 0.9737 > 0$, the time delay $\tau_{0} = 0.8064$, and then the unique positive root $\kappa_{0} = \sqrt{s_{0}} = 0.4178$ is obtained by solving the characteristic Eq (3.6), and the other two roots are $-3165.2$ and $-0.0067$.

According to Theorem 4.1, main conclusions are described as follows

Case i: When $\tau\in[0, \; 0.8064]$, the hydrological model (2.3) is asymptotically stable near $E_{0}$, as shown in Figure 6. When $\tau > 0.8064$, the model (2.3) is unstable near $E_{0}$.

Case ii: When $\tau = 0.8064$, the Hopf bifurcation of the hydrological model (2.3) occurs near $E_{0}$, as shown in Figure 9.

In addition, it is also obtained that $Re[A_{1}] = 1.0727 > 0$ and $Re[B_{1}] = -866336 < 0$, the model (2.3) has a supercritical nontrivial periodic solution in $0 < \tau-0.8064\ll1$, and the periodic solution is asymptotically stable on the central manifold. In the savanna, the biological community under the joint action of fire disturbance and groundwater circulation is a transition zone composed of tree layer and grass layer. According to the analysis results of dynamic bifurcation theory, the soil water, the frequency, and the intensity of fire, as well as the resulting distribution trend and growth law of tree and grass species density are simulated with the help of software MATLAB.

Figure 6 shows that in $t\in\left[0, \; 150\right]$, the time delay $\tau$ caused by fire frequency has a significant impact on the ecological coexistence balance of vegetation. We simulated the growth trend of trees and grasses of the model (2.3) in $\tau\in\left[0, \; 0.8064\right]$, and eventually tends to be steady, see Figure 6. In which, Figure 6b is the enlarged area of the Figure 6a in $t\in\left[40, \; 120\right]$. In order to further observe the steady state, the density distribution solution curves of trees $u_2(t)$ and grasses $u_3(t)$ are described with the time delay $\tau = 0.1, \; 0.3, \; 0.5, \; 0.7, \; 0.76$ respectively, as shown in Figure 7. It is found that as the value of $\tau$ approaches the critical value $\tau_0$, the stable state of the solution gradually changes from the stationary state to periodic state, and finally convergence to the equilibrium point $E_{0}$. In addition, Figure 8 shows the phase diagram of $u_2(t)$ and grasses $u_3(t)$ changing with time. The closer to the critical value $\tau_0$, the more obvious the periodic orbit of vegetation ecosystem, which is consistent with the theoretical analysis. Its biological significance is that the delay $\tau$ of fire frequency disturbs in stable region, which will not affect the stable growth of the standing vegetation density in both trees and grasses. As $\tau$ approaches the critical value and disturbs near it, the ecological cycle of trees and grasses under the coupling of fire and hydrology will show an obvious periodic oscillation.

Figure 9 shows that in $t\in\left[0, \; 200\right]$, the change trend of the co-existence equilibrium of soil water, trees and grasses when $\tau\approx0.806$. The figure shows that the solution curve of soil water $u_{1}(t)$ shows an increasing tendency, and eventually tends to be stable. This is because water resources are mainly supplemented by precipitation such as rainstorms and floods, which are introduced into the system as a gradient function ^[1]. Moreover, fire has fewer effects on soil water, so that it can be ignored. After a period of stable growth, the density distribution rate of trees $u_{2}(t)$ (grasses $u_{3}(t)$) will eventually remain close to the initial value $E_{0}$ and coexist in perfect harmony after a short period of fall (rise), as shown in Figure 9a. From the perspective of hydrological effect, in the process of simulation, the annual precipitation $p = 560$ $mm$ is taken, which corresponds to the ecosystem in arid and semi-arid areas. At this time, the limited water resources restrict the density of dominant competitors (tree species) and are not conducive to the stable growth of trees. From the perspective of fire disturbance, the impact of fire on trees is fatal, and the delay effect $\tau$ of fire frequency still oscillates around the critical value $\tau_{0}$ of Hopf bifurcation, that is, the potential dynamic balance among tree species under fire disturbance is about to lose stability, which is not conducive to the stable growth of trees. To sum up, the biomass of tree will gradually decrease, while the biomass of grass will increase with the decrease of competitors (tree species). In addition, compared with grass species, tree species have a greater demand for water resources. Therefore, when the biomass of tree decreases, the storage of soil water will increase gradiently along with the decreasing demand for water resources. Under the background of fire disturbance and hydrological effect, the co-existence mode of trees and grasses of the model (2.3) is simulated, which shows the form of periodic oscillation, as shown in Figure 9b. However, the periodic states is a "transient state". If the time delay $\tau$ caused by fire frequency exceeds the critical value, the distribution, existence and succession of trees and grasses will be unstable and the coexistence balance of biological communities will be destroyed.

In Figures 10–12, the solution curve of soil water $u_{1}(t)$, trees $u_{2}(t)$ and grasses $u_{3}(t)$ along with time are simulated respectively. When $\tau$ is disturbed near the critical value $\tau_{0}$, the structure of the solution shows the obvious periodical change characteristics. In particular, the periodic characteristic of $u_{1}(t)$ is more significant than $u_{2}(t)$ and $u_{3}(t)$ at the same time. This is because water resources are recorded in the form of precipitation and are always "increasing". Trees and grasses not only compete for water resources, but also are affected by fire, so the periodic amplitude of density distribution rate is slightly different. When bifurcation conditions are satisfied, and the selection of environmental parameters has certain practical significance, Hopf bifurcation will lead to the ecological mode of the hydrologic model (2.3) exhibiting periodic oscillation features. That is to say, as the delay $\tau$ approaches the threshold, the soil water, trees and grasses present an obvious form of periodic oscillation, and the numerical simulation results are consistent with the theoretical analysis.

6.   Conclusions

Considering a series of conditions caused by the delay effect of fire frequency, we have established a hydrological model (2.3) with delay effect (time delay) in this paper. We obtain sufficient conditions for local stability of positive equilibrium point $E_{0}$, and the existence condition for Hopf bifurcation. The study found that, when $\tau > {\tau _0}$, a coexistence mechanism of tree species and grass species will lose stability at $E_{0}$, while on one side of the $\tau _0$, the soil water, trees and grasses would change to the periodic oscillation ecological-pattern by Hopf bifurcation. Then the stability of the Hopf bifurcation periodic solution is given by using the center manifold theorem and the normal form theorem. Finally, the paper provides sensitivity analysis for uncertain factors, such as the frequency and intensity of fire in the model (2.3) without $\tau$. The analysis showed that the effect of fire frequency on hydrological effect, and the density distribution of trees and grasses was more significant relative to fire intensity under the same conditions. All of the above further verified that fire is a key control parameter of vegetation coexistence equilibrium mechanism.

To analyze the issue of the time delay driven by fire frequency, we also investigate the potential dynamic properties of soil water, trees and grasses in the model (2.3) under the coupling effect of fire disturbance and hydrological effects are studied in the time domain by using bifurcation theory and sensitivity analysis methods. Dynamic theoretical analysis demonstrated the existence of an ecological phenomenon periodic oscillations of soil water-tree-grass under the fire disturbance and hydrological effect, which is not directly observable in nature.

Acknowledgments

This research was supported by Doctoral Scientific Research Foundation of Inner Mongolia Minzu University (No.BS648).

Conflict of interest

The authors declare there is no conflicts of interest.