Spreading Speeds and Traveling Waves for Monotone Systems of Impulsive Reaction–Diffusion Equations: Application to Tree–Grass Interactions in Fire-prone Savannas

Many systems in life sciences have been modeled by reaction–diffusion equations. However, under some circumstances, these biological systems may experience instantaneous and periodic perturbations (e.g. harvest, birth, release, fire events, etc) such that an appropriate formalism like impulsive reaction–diffusion equations is necessary to analyze them. While several works tackled the issue of traveling waves for monotone reaction–diffusion equations and the computation of spreading speeds, very little has been done in the case of monotone impulsive reaction–diffusion equations. Based on vector-valued recursion equations theory, we aim to present in this paper results that address two main issues of monotone impulsive reaction–diffusion equations. Our first result deals with the existence of traveling waves for monotone systems of impulsive reaction–diffusion equations. Our second result tackles the computation of spreading speeds for monotone systems of impulsive reaction–diffusion equations. We apply our methodology to a planar system of impulsive reaction–diffusion equations that models tree–grass interactions in fire-prone savannas. Numerical simulations, including numerical approximations of spreading speeds, are finally provided in order to illustrate our theoretical results and support the discussion.


Introduction
In nature, all organisms migrate or disperse. This can take diverse forms such as walking, swimming, flying, or being transported by wind or flowing water, see Shigesada and Kawasaki [44]. Such a migration, or dispersion, can be to some extent related to human activities that bring drastic changes in the global environment. According to [44], dispersive movements become noticeable when an offspring or a seed leaves its natal site, or when an organism's habitat deteriorates from overcrowding. The spatially explicit ecological theories of such events have been made possible thanks to successful development of mathematical models that have played a central role in the description of migrations, see Friedman [18], Shigesada and Kawasaki [44], Okubo and Levin [37], Cantrell and Cosner [12], Volpert [52], Logan [33,34], Perthame [39], and the references therein.
Mathematical literature dealing with the species' spread mostly relies on reaction-diffusion equations that assume that the dispersal is governed by random diffusion and that it, along with the growth processes, take place continuously in time and space, (Cantrell and Cosner [12], Lewis and Li [27]). This approach has had a remarkable success in explaining the rates at which species have invaded large open environments, see Shigesada and Kawasaki [44], Okubo and Levin [37], Cantrell and Cosner [12], Lewis and Li [27], Volpert [52], Logan [33,34] and Perthame [39]. However, it is well-known that ecological species may experience several phenomena that, depending on circumstances, can be either time-continuous (growth, death, birth, release, etc.), or time-discrete (harvest, birth, death, release, etc.), see also Ma and Li [35], Dumont and Tchuenche [13], Yatat Djeumen et al. [60], Yatat Djeumen [58] and the references therein. In the case of time-discrete perturbations, whose duration is negligible in comparison with the duration of the process, it is natural to assume that these perturbations act instantaneously; that is, in the form of impulses (Lakshmikantham et al. [26], Bainov and Simeonov [6]). Hence, there is a need to create a meaningful mathematical framework to analyze models leading to, say, impulsive reaction-diffusion equations.
There are several works that considered the impact of impulsive perturbations on the dynamics of a system, both in the space-implicit and space-explicit case (see also Sects. 1.1 and 1.2). Very often, impulsive perturbations in space-implicit mathematical models result in the occurrence of periodic solutions in the model (e.g. Ma and Li [35], Yatat Djeumen et al. [60] and references therein). For space-explicit impulsive mathematical models in bounded domains, in addition to periodic solutions that may occur, the issue of minimal domain has been also addressed (e.g. Lewis and Li [27], Yatat Djeumen and Dumont [61]). On the other hand, in unbounded domains the problems of the existence of traveling wave solutions and the computation of the spreading speeds are rarely addressed (see below).
The study of traveling waves, as well as the computation of spreading speeds for monotone systems of reaction-diffusion equations, have been done by several authors (e.g. Weinberger et al. [56], Weinberger [54,55], Lewis et al. [28], Li et al. [29], Volpert [52], Yatat Djeumen et al. [59] and references therein). However, little is known about the case of impulsive reaction-diffusion equations (see e.g. Yatat Djeumen and Dumont [61], Fazly et al. [15]). For a scalar impulsive reaction-diffusion equation, Yatat Djeumen and Dumont [61] considered the Fisher-Kolmogorov-Petrowsky-Piscounov (FKPP) equation and obtained conditions under which an invasive traveling wave, connecting the extinction steady state and the positive steady state may exist (see also Lewis and Li [27]). To the best of our knowledge, the existence of traveling waves for system of impulsive reaction-diffusion equation was studied only in Huang et al. [23], and only in a particular case. Precisely, the authors considered a stage-structured population model, were only one stage (or state variable) experiences a spatial diffusion, while the others are stationary. This assumption leads to a partially degenerate system of impulsive reaction-diffusion equations. Moreover, they also assumed the reaction term for the diffusing state variable to be linear.
In the following, we provide a brief literature review of both space-implicit and spaceexplicit impulsive mathematical models.

Modelling Releases as Pulse Events
In the framework of biological control of pests or vectors of infectious diseases, the sterile insect technique (SIT) is one of the promising ones. SIT control generally consists in massive releases of sterile insects in the targeted area in order to eliminate, or at least to lower the pest population under a certain threshold (Anguelov et al. [4]). Generally, SIT releases are done periodically. That is why several authors modelled the release process as a periodic impulsive event, while keeping the continuous-time differential equation framework for the birth, growth, death, or mating process (e.g. White et al. [57], Dumont and Tchuenche [13], Strugarek et al. [45], Bliman et al. [9] and the references therein). Based on the qualitative analysis of their impulsive models, the authors were able to derive meaningful relations between the period and the size of the releases in order to achieve the elimination of the vector or pest population in the long term.
In the context of interacting species such as prey-predator interactions, there exist mathematical models that tackled periodic releases of one of the interacting species (e.g. prey only, predator only). They considered periodic impulses as to model periodic release events. The authors found relations involving the pulse time period and the amount of released species that precluded extinction, in the long term dynamics, of interacting populations (see for instance Zhang et al. [64], Zhao et al. [65]).

Modelling Fires as Pulse Events in Tree-Grass Interactions in Fire-prone Savannas
Maintaining the balance between the grass and the trees in savanna is of utmost importance for both human and animal populations living in such areas. The problem is that in typical circumstances the trees encroach on the grassland making the environment inhabitable for many species. It turns out that periodic fires, either natural or man-made, are one of the way to maintain an acceptable equilibrium. Thus, several mathematical models have been developed to study tree-grass interactions in fire-prone savannas (see the review of 1 3 Yatat Djeumen et al. [63]). It is well-known that in such savannas, grass biomass that dries up during dry seasons is used as fuel for fires (see op. cit. and references therein). Hence, low grass biomass will result in fires of lower intensity. We emphasize that here we are not dealing with 'forest fires', like that usually occurring in Europe, North America, or, sometimes, in Australia, but only with grass-induced fires. Thus, in particular, a dense tree population results in a lower fire intensity and consequently, in a lower tree destruction. Some of the models take into account fire as a time-continuous forcing in tree-grass interactions. However, as pointed out in Yatat Djeumen et al. [63], it is questionable whether it makes sense to model fire as a permanent forcing that continuously removes a fraction of the fire sensitive biomass. Indeed, since several months and even years can pass between two successive fires, they can be rather considered as instantaneous perturbations of the savanna ecosystem (see also Yatat Djeumen [58], Yatat Djeumen et al. [60], Tchuinté et al. [47,48]). Several recent papers have proposed to model fires either as stochastic events, while keeping the continuous-time differential equation framework (Baudena et al. [7], Beckkage et al. [8], Synodinos et al. [46]), or by using a time-discrete model (Higgins et al. [22], Accatino et al. [1,2], Klimasara and Tyran-Kamińska [24]). However, a drawback of many of the aforementioned recent time-discrete stochastic models (Higgins et al. [22], Baudena et al. [7], Beckage et al. [8]) is that they hardly lend themselves to analytical treatment. Thus, on the basis of recent publications (see for instance Yatat Djeumen et al. [63] and the references therein), we consider fires as impulsive time-periodic events. While certainly an approximation, such an approach results in impulsive differential equation models which are a good compromise combining the impact of time-discrete fires with a time-continuous process of the vegetation growth. As such, they are analytically tractable, while at the same time remain reasonably realistic. Now we recall the minimalistic tree-grass interactions model with pulse fires that will be used later in the paper (see Sect. 3). We assume that the trees and grass form an amensalistic system in which grass is harmed by trees (which, for instance, block the sunlight) but itself does not affect them. The fires occur periodically every ̃ units of time. We denote by T n (resp. G n ) the tree (resp. grass) biomass during the inter-fire season number n ∈ ℕ * , where ℕ * = {1, 2, …} . Following the formalism of the recursion equations (see Weinberger et al. [56], Yatat Djeumen and Dumont [61]), the resulting minimalistic system of equations governing the trees-grass interactions with periodic fires is given by (see also Yatat Djeumen et al. [63] and references therein) with non negative initial conditions (T 1 (0), G 1 (0)).
Here, between two successive fires; that is, in the inter-fire season n, the dynamics of both the tree and grass biomasses is modelled by the first two equations of system (1), where G ( ) and T ( ) denote the unrestricted rates of growth of the grass and the tree (1)

3
biomass, respectively, while K G ( ) and K T ( ) are the carrying capacities for grass and the trees, respectively. All these functions are assumed to be increasing and bounded functions of the water availability which is supposed to be known: see [62] for further details about the definition of the growth rates and the carrying capacities with respect to . Further, G and T denote, respectively, the rates of the grass and the tree biomass loss due to natural causes, herbivores (grazing and/or browsing) or human actions, while TG denotes rate of the loss of the grass biomass due to the existence of trees per units of the biomasses.
At the end of the inter-fire season a fire occurs and impacts both the tree and grass biomasses. Thus there is an update of the biomasses at the beginning of the next inter-fire season. This event is modelled by the last two equations of system (1). We model the fire intensity, denoted by w G , using the Holling III functional response (see (66)) ; that is, we assume that w G is an increasing function tending to 1 as the argument increases and satisfying w G (0) = w � G (0) = 0 . This reflects the fact that the fire becomes significant only when the grass biomass becomes sufficiently large. Impulsive fire-induced tree/shrub mortality, denoted by , is assumed to be a positive, decreasing, and nonlinear function of the tree biomass (see (67)). Indeed, fires affect differently large and small trees since fires with high intensity (flames more than ca. 2 metres high) cause greater mortality of shrubs and topkill of trees while fires of lower intensity (flames less than ca. 2 metres high) kill only shrubs and subshrubs (see also [62,63] and references therein). Further, is the specific loss of the grass biomass due to the fire. To avoid the extinction of either T n , or G n , we assume that (see also Yatat Djeumen et al. [63]) System (1) will be studied in Sect. 3.
Readers are referred to Yatat Djeumen et al. [63] for the derivation and analysis of system (1) following the formalism of Lakshmikantham et al. [26].

Space-Explicit Impulsive Models
The formulation of space-explicit impulsive models generally consists in the addition of local or non-local spatial operators to a temporal impulsive model. In Akhmet et al. [3], Li et al. [30] and Liu et al. [32], the authors considered impulsive reaction-diffusion equations to model spatio-temporal dynamics of ecological species with prey-predator interactions and experiencing pulse and periodic perturbations like harvest, release, etc. The spatial movement of species was modelled by the Laplace operator with a constant diffusion rate. Qualitative analysis of these models was done by using the theory of sectorial operators (Henry [21], Rogovchenko [41,42], Li et al. [30]) and comparison arguments (Rogovchenko [40], Walter [53], Liu et al. [32], Akhmet et al. [3]). More precisely, the authors obtained some conditions involving the pulse time period that ensured the permanence of the predator-prey system and the existence of a unique globally stable periodic solution (Akhmet et al. [3], Li et al. [30], Liu et al. [32]). Vasilyeva et al. [50] dealt with the question of persistence versus extinction in a single population model featuring a non-local impulsive reaction-advection-diffusion model for an insect population. The non-local term was used to describe the dispersal of the adult insects by flight. The authors employed a dispersal kernel that gave the probability density function of the signed dispersal distances.
We note that the study of the species spread and their wave speeds, when they experience impulsive and periodic perturbations, in the case of scalar equations was done in Vasilyeva et al. [50], Lewis and Li [27], or Yatat Djeumen and Dumont [61]. However, until recently, systems of impulsive reaction-diffusion have not received much attention (Huang et al. [23] and Fazly et al. [15]). We aim to address this question here by extending the minimalistic trees-grass interactions model (1). To this end, we assume that both the woody and herbaceous plants can propagate in space through diffusion; see Yatat Djeumen et al. [59] for a discussion of the construction of the trees-grass interactions partial differential equations models. The resulting minimalistic system of impulsive reaction-diffusion equations is then given by with given non negative initial conditions In system (3), d T ( ) and d G ( ) denote the woody, respectively, herbaceous biomass spatial vegetative diffusion coefficient, while the remaining coefficients and assumptions on them are as in (1). The aim of this paper is to give some insights into the existence of traveling wave solutions for monotone systems of impulsive reaction-diffusion equations as well as the computation of their spreading speeds. Precisely, we use the vector-valued recursion theory proposed by Weinberger et al. [56] (see also Weinberger [54,55], Lewis et al. [28], Li et al. [29], Lewis and Li [27], and references therein) to develop a framework in which we are able to deal with the existence of traveling wave solutions for monotone systems of impulsive reaction-diffusion equations and the computation of spreading speeds. After this long introduction, Sect. 2 presents the framework for the computation of the spreading speeds and the existence of traveling wave solutions for monotone systems of impulsive reaction-diffusion equations. Section 3 deals with the application of this framework to a system of two impulsive reaction-diffusion equations that models tree-grass interactions in fire-prone savannas. We also provide some numerical illustrations of our theoretical results and, in particular, we show an approximation of the spreading speeds. The paper ends with a conclusion, where further extensions are discussed.

Preliminaries and Hypotheses
In this section we introduce the basic notation, definition and results that allow us to deal with the issues of the existence of traveling wave solutions and/or computation of the spreading speeds for impulsive reaction-diffusion (IRD) systems. Let us denote: In the case when there are no impulsive perturbations, the reaction-diffusion system is written as together with sufficiently smooth and nonnegative initial condition For -periodic impulsive perturbations, we consider the whole time interval, [0, +∞) , as a succession of inter-perturbation seasons of length . Let us denote the state variables at time t ∈ [0, ] and location x during the inter-perturbation season n ∈ ℕ * as n (t, x) = (u n,1 , u n,2 , ..., u n,N ) . Following the recursion formalism (Weinberger [54,55], Weinberger et al. [56], Lewis and Li [27], Vasilyeva et al. [50], Fazly et al. [16], Huang et al. [23], Yatat Djeumen and Dumont [61]), the impulsive reaction-diffusion system is written as together with the updating condition and with sufficiently smooth and nonnegative initial data 1 (x, 0) = 0 (x) . We note that (3) is a special case of (8).
Our work is based on the results of Li et al. [29] (see also Weinberger et al. [56]) concerning the existence of monostable traveling wave solutions and the computation of spreading speeds for systems of reaction-diffusion equations. We note that they focused on the case, when was just the time--map operator solution, Q , of system (7). Since, however, the reduction of an IRD system to the recursion form does not depend on the updating condition, the results of op. cit. on the existence of traveling waves for recursions and the spreading speeds determined by them can be used verbatim to the recursion obtained from (7)-(8). Thus we recall the relevant results from [29].
We first assume that , and the initial data are such that the IRD system (7)-(8) admits a unique nonnegative classical solution for each n ∈ ℕ * (Zheng [66], Volpert [52], Logan [33,34], Perthame [39]). We begin with some notation (see Weinberger et al. [56], Li et al. [29]). For two vectorvalued functions (x) and for all i and x. We use the notation 0 for the constant vector whose all components are 0. If ≫ is a constant vector, we define the set of functions Let Q be the time--map operator solution of system (7). Then (8) can be written as with the initial condition 0 (x).
In the sequel, we recall the key assumptions of Li et al. [29], related to the operator Q defined in (9). For a fixed y ∈ ℝ, the translation operator by y is defined by

Hypothesis 2.1
i. The operator Q is order preserving in the sense that if and are any two functions In biological terms, the dynamics is cooperative. ii. Q[ ] = , there is a constant vector ≫ such that Q[ ] = , and if 0 is any constant vector with 0 ≫ , then the constant vector n , obtained from the recursion (9), converges to as n approaches infinity. This hypothesis, together with (i), imply that Q takes C into itself, and that the equilibrium attracts all initial functions in C with uniformly positive components. There may also be other equilibria lying between and the extinction equilibrium 0, in each of which at least one of the species is extinct. iii. Q is translation invariant. In biological terms this means that the habitat is homogeneous, so that the growth and migration properties are independent of location. iv. For any , ∈ C and any fixed y, | is sufficiently small on a sufficiently long interval centered at y.
v. Every sequence n in C has a subsequence n l such that Q[ n l ] converges uniformly on every bounded set.
We are now in position to recall results of Li et al. [29] that deal with spreading speeds as well as traveling wave solutions for the IRD system (7)-(8), rewritten following the recursion formalism (9). In the sequel, we assume that Hypothesis 2.1. holds for the recursion operator Q of system (9).
Following Li et al. [29] (see also Weinberger et al. [56]), we consider a continuous ℝ N -valued function (x) with the properties We let, for all fixed c ∈ ℝ , 0 (c;s) = (s) , and define the sequence n (c;s) by the recursion Li et al. [29] showed that the sequence n converges to a limit function (c;s) such that (c; ± ∞) are equilibria of Q and (c;∞) is independent of the initial function . Following this, they defined the slowest spreading speed c * ≤ ∞ by the equation The following result holds.
Theorem 1 [29, Theorem 2.1] There is an index j for which the following statement is true: Suppose that the initial function 0 (x) is 0 for all sufficiently large x, and that there are for all x and 0 ≥ for all sufficiently negative x. Then for any positive the solution n of recursion (9) has the properties and That is, the j th component spreads at a speed no higher than c * , and no component spreads at a lower speed.
In order to define the fastest speed c * f , we choose with the properties (10), and let n (x) be the solution of the recursion (9) with 0 (x) = (x) . Following Li et al. [29], we define the function Li et al. [29] showed that (c;∞) is independent of the choice of the initial function as long as has the properties (10). We therefore can define the fastest spreading speed c * f ≥ c * by the formula The following result holds.

Theorem 2 [29, Theorem 2.2]
There is an index i for which the following statement is true: Suppose that the initial function 0 (x) is 0 for all sufficiently large x, and that there are positive constants 0 < ≤ < 1 such that ≤ 0 ≤ for all x and 0 ≥ for all sufficiently negative x. Then for any positive the solution n of the recursion (9) has the properties and That is, the ith component spreads at a speed no less than c * f , and no component spreads at a higher speed.
Let ̂ n , n ≥ 0, be the solution to the recursion with ̂ 0 (x) = 0 (x) . Recall that a traveling wave of speed c is a solution of the recursion (9) which has the form n (x, 0) = (x − nc) with (s) a function in C ; that is, the solution at time n + 1 is simply the translate by c of its value at n. Then such a traveling wave defines a traveling wave solution for (7)-(8) in the following sense. By (9) we have and thus, by Hypothesis 2.1. (iii)., for n ≤ t < (n + 1) we have n+1 (x, t) = 1 (x − nc, t) . We observe that since the model is translation invariant, we obtain a traveling wave for the system (7) without the updating conditions.
Using the definitions of c * and c * f , we have the following result that deals with the existence of traveling wave solutions for the IRD systems (7)- (8). In practice, assumptions (iii), (iv). and (v) are typically satisfied for biologically reasonable (impulsive) models. The most challenging assumptions for IRD systems (7)-(8) are (i) and (ii).

Application to a Minimalistic Trees-Grass Interactions IRD System
In this section, we consider the minimalistic tree-grass interactions IRD system (3)-(4). Using a similar normalization procedure to Yatat Djeumen et al. [59] (see also Appendix A), system (3)-(4) becomes together with the updating conditions and sufficiently smooth and nonnegative initial data U 1 (x, 0), V 1 (x, 0) . We are now looking for traveling wave solutions as well as the spreading speeds involving semi-trivial steady states; that is, steady states where either U n = 0, or V n = 0 but not simultaneously U n = 0 = V n = 0.

Basic Properties of (19)-(20)
Let C ub (ℝ) be the Banach space of bounded, uniformly continuous function on ℝ and C ub (ℝ) and C 2 b (ℝ) are endowed with the following (sup) norms and is a Banach space. We recall that we assumed that w V was an increasing C 1 (ℝ) function such that for all V ∈ ℝ, Similarly, is a decreasing C 1 (ℝ) function such that for all U ∈ C ub (ℝ), For simplicity we note (U n (x, 0), V n (x, 0)) = (U n,0 (x), V n,0 (x)). In the sequel, we first address the question of the existence and uniqueness of solutions of the reaction-diffusion (RD) system (19) in unbounded domains.
For fixed n ∈ ℕ * , we set = (w 1 , w 2 ) ∶= (U n , V n ) . System (19) can be written as the abstract Cauchy problem where in the Banach space B = C ub (ℝ) × C ub (ℝ) we have We shall consider (25) as a nonlinear perturbation of the linear part that, in this case, consists of two uncoupled diffusion equations. Thus, the corresponding semigroup is the diagonal semigroup consisting of the Gauss semigroups where, for d = d u , d v and f ∈ X, Our problem is posed on the whole line and thus comparison theorems for the solutions are a little more delicate. Though in various forms they appear in many papers, see e.g. [17,38,49,51] and references therein, a comprehensive proof of them, starting from the first principles, is difficult to find. Therefore, for pedagogical reasons as well as because we shall need some intermediate estimates later on, we decided to provide a simple proof for the problem at hand that uses the positivity of the semigroup {S d (t)} t≥0 and the triangular structure of the nonlinearity in (19). In fact, the semigroup for the scalar problem, Since the updating conditions (20) are non-increasing, if the initial data U 1 (⋅, 0), V 1 (⋅, 0) satisfy then, for each n ∈ ℕ * , the solutions (U n , V n ) of system (25) satisfy

The First Coordinate Change
System (19) is monotone competitive and system (20) is not monotone. Therefore, the full system (19)-(20) is not monotone. Hence, in order to be able to apply results of Li et al. [29], we first apply a coordinate change in order to obtain a monotone cooperative system. We set so that system (19)-(20) is transformed to together with the updating conditions Properties (34), (33) and (23) imply that w v is a decreasing C 1 (ℝ) function such that We also deduce that system (36) is monotone cooperative and the sequence defined in (37) is monotone increasing. Hence system (36)-(37) is monotone cooperative as long as the initial conditions belong to [0, 1].

Space Homogeneous Steady States of System (36)-(37)
In this section we compute the space homogeneous steady states of system (36)-(37) by solving the fixed point problem associated to system (43).

and only if
In this study, we are mainly concerned with the existence of traveling wave solutions of system (36)-(37) involving steady states u and v computed in Lemma 1. In order to use the results of Li et al. [29], we translate the steady state v to 0 through another coordinates change.

The Second Coordinate Change
Recall that and Recall also that 0 ≤v < 1 is equivalent to R 0 > 1. Therefore, in this section we assume that We set Hence, system (36)-(37) becomes together with the updating conditions As previously, we deduce from properties (38) that w is a decreasing C 1 (ℝ) function such that For simplicity, we set (u n,0 , q n,0 ) = (u n (x, 0), q n (x, 0)) . We also set n,0 = (u n,0 , q n,0 ) and let Q denote the time--map operator solution of system (47). Then where (u, q) = (H 1 (u, q), H 2 (u, q)).
Using the coordinate change (46), system (43) becomes and, using Lemma 1, it is straightforward to deduce that system (47)-(48) admits as space homogeneous steady states
(1 − w(q n ( )) (u n ( )))u n ( ) e − + (1 − e − )(1 − w(q n ( )) (u n ( )))u n ( ) =∶ F 1 (u n ( ), q n ( )), For every u, q ∈ ℝ + , the function s ↦ Z(u, q, s) is continuous on the interval [0, ] . In addition, exists and is continuous for all (u, q, s) Similarly, For convenience, we set (see (51)) Computing the partial derivatives of A 1 , A 2 , B 1 and B 2 defined in (51) gives denote the Jacobian matrix of (51). Using (52), the quotient rule and the properties of w (see (49)) and (see (24)), we obtain the following results: -Local stability of 0 . At the steady state 0 = (0, 1 −v) the matrix J has the following entries: Indeed, for example, (49)). Since e is an eigenvalue of J at 0 and e > 1, the steady state 0 is unstable. -Local stability of u . At the steady state u = (1, 1 −v) the matrix J has the following entries: Eigenvalues of the Jacobian matrix at u are e − and R 1 with e − < 1 . Therefore, u is locally asymptotically stable (LAS) whenever R 1 < 1. -Local stability of v . J at the steady state v = (0, 0) has the following entries: Knowing the explicit value of J 21 is not necessary since J 12 = 0 . Eigenvalues of the Jacobian matrix at v are J 11 and J 22 . Recall that we assumed R 0 > 1 . Hence, J 22 < 1 . Therefore, the steady state v is LAS whenever

Application of the Results of [29]
In the sequel, we study the recursion operator Q defined in Eq. (50) and we check if it satisfies Hypotheses 2.1 of Li et al. [29]. Recall that n,0 = (u n,0 , v n,0 ) and Q is the time--map solution operator of reaction-diffusion system (47). We consider the order inter- is the positive coexistence steady state defined in Lemma 2. where as well as the convolution ⋆ , were defined in (28) and the spatio-temporal convolution is given by 3. Let (u, q) and (v, p) be two solutions of system (47) initiated at (u 0 , q 0 ) and (v 0 , p 0 ) respectively. We assume that (u 0 , q 0 ) , (v 0 , p 0 ) ∈ C u , hence (u, q) and (v, p) are also u(y, s), q(y, s))dyds.
uniformly bounded. Let {S c d (t)} t≥0 denotes the (positive) semigroup solving (29) for some function c satisfying |c(x, t)| ≤ L . Using |S d (t)u 0 | ≤ S d (t)|u 0 |, u 0 ∈ C ub (ℝ), where {S d (t)} t≥0 is the diffusion semigroup (28), and the Phillips-Dyson expansion to (29) we ascertain that for any u 0 As before, we begin with solutions u and v to (47) 1 . Repeating the argument leading to (32), we see that z(x, t) = u(x, t) − v(x, t) can be estimated as Let, for > 0 , r > 0 be such that Then let us fix x and let |u 0 (x) − v 0 (x)| ≤ ≤ ∕2e L on (x − r, x + r) so that we obtain for 0 < t ≤ By choosing appropriate r, we see that the estimate is valid for x in any given bounded subset of ℝ.
Considering now (47) 2 , we see that Z(x, t) = q(x, t) − p(x, t) is a solution to and considerations as above show that |Z(x, t)| ≤ e L 1 t (x, t) , where In the above, L 1 is a constant bounding |1 − (p + q) − u|, 0 ≤ p, q, u ≤ 1 , and we used In order to obtain properties of the operator Q defined in (50), we first formulate results for the nonlinear operator ∶ B → B . Recall that the Banach space considered here is B = C ub (ℝ) × C ub (ℝ) endowed with the sup-norm. For = (u, q)  Proof Let = (u, q) and In the sequel, we assume that R 0 = (1 − )e > 1 and R 1 = (1 − )e (1− ) < 1 ; that is, the coexistence steady state u = (1, 1 −v) exists and is stable. We also assume that Taking into account Lemma 5, we deduce that the recursion operator Q defined in (50) verifies all conditions of Hypothesis 2.1. Consequently, we can apply the results of Li et al. [29] that deal with the spreading speeds and existence of traveling wave solutions for systems (47)- (48). Recall that a traveling wave of speed c is a solution of the recursion (50) which has the form n (x, 0) = (x − nc) with (s) being a function in C u . That is, the solution at time n + 1 is simply the translate by c of its value at n. Using the definition of c * (see (12)) and c * f (see (15)), the following result holds true.

Numerical Simulations
In this section we provide numerical simulations of the impulsive tree-grass reaction-diffusion model (3)-(4). We note that the parameters with ̃ below refer to this system and can be derived from the corresponding parameters related to the normalized system (47)-(48), see Appendix A. Thus, we consider fire events as periodic and pulse perturbations with the time period ̃ . The form of the functions G ( ) , T ( ) , K G ( ) , K T ( ) , and w G is considered following Yatat Djeumen et al. [63]. The readers are referred to Appendix A for their definition and parametrization. The parameter values used in the following simulations are also given in Appendix A: see Tables 1 and 2.
Using the parameters values given in Table 1, Fig. 1 depicts the spreading of tree and grass biomasses toward the stable forest homogeneous steady state T = (273.3955, 0) . In this case, R 0 = 1.3667 , R 1 = 0.0058 and R 2 = 3.3801 . Recall that T is LAS whenever R 1 < 1 , while the grassland homogeneous steady state G = (0, 3.3888) exists when R 0 > 1 and is LAS whenever R 2 < 1 . In terms of tree-grass interactions, Fig. 1 illustrates the spreading of forest or the so-called 'forest encroachment' phenomenon (Yatat Djeumen et al. [59]).
In the setting of the forest encroachment phenomenon, we carry out numerical simulations to compute the spreading speed of forest biomass. We investigate the relationship between the tree biomass diffusion coefficient and its spreading speed. To estimate the spreading speed of the tree biomass that undergoes a forest encroachment, grass biomass diffusion coefficient d G is kept constant and equal to 0.002 while tree biomass diffusion coefficient d T varies in the range [0.001, 0.9]. In the diffusive logistic equation, a linear relationship is obtained between the wave speed and the square root of the diffusion coefficient (e.g. Volpert [52], Yatat Djeumen et al. [59], Yatat Djeumen and Dumont [61]). Hence, for the tree biomass, we consider an equation of the form c T (d T ) = a 1 d a 2 T to be fitted for the data shown in Fig. 2(b), where c T ∈ [0.0865, 1.6899] . We found that a 1 ∈ (1.7624, 1.7861) and a 2 ∈ (0.4816, 0.4911) with 95% confidence. In fact, a 1 = 1.7743 and a 2 = 0.48634 , with r 2 = 1 , indicating that 100% of the variance of the data is explained by the equation.
With the parameter values given in Table 2, Fig. 3 illustrates the spreading of both tree and grass biomasses toward the grassland homogeneous steady state G = (0, 2.1096) . In this case, R 0 = 2.0425 , R 1 = 1.2541 and R 2 = 0.9932 . Recall that G exists when R 0 > 1 and is LAS whenever R 2 < 1 . We further investigate the relationship between the diffusion coefficient and the spreading speed of the grass biomass. We assume that d T  G was fitted to the data shown in Fig. 4b, where c G ∈ [0.0465, 1.0877] . We found that a 1 ∈ (1.0829, 1.0935) and a 2 ∈ (0.4839, 0.4915) with 95% confidence. In Remaining parameters are in Table 1 1 3 Remaining parameters are in Table 2 fact a 1 = 1.0882 and a 2 = 0.48769 , with r 2 = 1 , indicating that 100% of the variance of the data is explained by the equation.

Conclusion
In this paper, we used the vector-valued recursion equation theory (e.g. Weinberger et al. [56], Weinberger [54,55], Lewis et al. [28], Li et al. [29]) to propose a framework in which we can systematically address the questions of the existence of traveling waves for monotone systems of impulsive reaction-diffusion equations and the computation of spreading speeds. This study extends the previous one that dealt with the impulsive Fisher-Kolmogorov-Petrowsky-Piscounov (FKPP) equation (Yatat Djeumen and Dumont [61]) but is restricted to monostable situations; that is, when only one of the steady states is stable. However, the bistable case (i.e. when two steady states are simultaneously stable) is also meaningful and needs to be studied for monotone systems of impulsive reaction-diffusion equations. Traveling waves in bistable reaction-diffusion systems without impulsive perturbations are treated in Volpert [52] (see also Yatat Djeumen et al. [60] for application in the context of bistable tree-grass reaction-diffusion model).
The computation of spreading speeds and the existence of traveling waves for bistable monotone systems of impulsive reaction-diffusion equations will be the aim of future studies. It first requires to elaborate a recursion equations theory that includes bistable cases and thus extending the results of Li et al. [29] that also deal only with monostable cases. Last but not least, in a recent work (Yatat Djeumen et al. [62]) we studied a space-implicit model that gives quite satisfactory results, in term of the model's predictions, across the rainfall gradient and improving already published results obtained with tree-grass interactions models of similar complexity or less. Hence, it is desirable to also study its spaceexplicit counterpart acknowledging impulsive fire events. The aim will be to analyse the impact of impulsive fires on the dynamics of the vegetation mosaics like forest-grassland and forest-savanna in humid tropical savannas (see also Yatat Djeumen et al. [59]).

3
where G, in tons per hectare (t.ha −1 ), is the grass biomass and G is the value taken by G, when the fire intensity is half of its maximum. However, the proofs are valid for general w G . -the function is defined by where T, in tons per hectare (t.ha −1 ), stands for the tree biomass, min fT (in year −1 ) is the minimal loss of tree biomass due to fire in systems with a very large tree biomass, max fT (in year −1 ) is the maximal loss of tree/shrub biomass due to fire in open vegetation (e.g. for an isolated woody individual having its crown within the flame zone), p T (in t −1 .ha) is proportional to the inverse of biomass suffering an intermediate level of mortality.
(69)  Recall that R 0 , R 1 , v and R 2 are related to the normalized system (47)- (48), while R 0 , R 1 , Ḡ and R 2 are related to the original system (3)-(4). In Tables 1 and 2, we summarize the parameter values that are used for numerical simulations. They are chosen according to [58,60,63].

Another Monotone Increasing Impulsive System
If, instead of the first coordinates change (35), one considers then the normalized systems (19)-(20) becomes together with the updating conditions As we mentioned in Sect. 3.3, this is a monotone increasing system and hence, reasoning as before, one can study the case where the stability of the steady states is reversed. (76) u n+1 (x, 0) = w v (v n (x, )) (1 − u n (x, ))(1 − u n (x, )) + u n (x, ), v n+1 (x, 0) = (1 − )v n (x, ).