SOME REMARKS ON AN ENVIRONMENTAL DEFENSIVE EXPENDITURES MODEL

. In this paper, we consider the environmental defensive expenditures model with delay proposed by Russu in [16] and obtain diﬀerent results about stability of equilibria in the case of absence of delay. Moreover we provide a more detailed analysis of the stability for equilibria and Hopf bifurcation in the case with delay. Finally, we discuss possible modiﬁcations of the model in order to make it more accurate and realistic.


1.
Introduction. We consider a model introduced in [16] that deals with management of tourism in protected areas (PAs). The model is based on the interactions among visitors V , quality of ecosystem goods E, and capital K, intended as accommodation and entertainment facilities in PAs, which is described by the system: where V d = V (t − τ ), τ ≥ 0 is the time delay andP is the maximum tolerable pollution stock P ; m 1 , m 2 , a, b, c, η, δ,P are strictly positive constants, 0 ≤ η < 1 and 0 < r < 1. A detailed discussion of the model can be found in [16] but it is worth mentioning that a part, η, of total revenues is invested to defend the environmental resources while the other part,(1 − η), is used to increase the capital stock K. Also, we point out that the quality of ecosystem goods is measured as E =P − P . The idea is to study the stability of the model, that is the sustainability of the tourism in PAs, depending on the parameters involved (see also [1], [2], [6]- [10], for more general discussions on this the topic). For this reason the existence and stability of positive equilibria play a crucial role in this context. Our analysis of the equilibria is different from that in [16], in particular we find a necessary and sufficient condition for the existence of a positive equilibrium (V ∞ , E ∞ , K ∞ ) inside R + × (0,P ) × R + that is (see next section for details): which is in contrast with (a) of Proposition 1 in [16]. We also consider interesting to revise (in section 3 below) the bifurcation analysis for τ > 0 providing more detailed results with respect to those in [16]. Moreover we prove the existence of stability switches for the fixed point in the model with delay. We did not prove that the hypotheses are compatible with the positiveness of E ∞ , but still we consider the result interesting enough. In fact, in an attempt of using the model in real cases, the existence of stability switches makes the strategy of the managers of PAs more complicate in order to stabilize the system (see section 3 below). In section 5 of this work we discuss some critical issues about the model and propose some changes in order to improve it.
2. Fixed point and stability analysis with τ = 0. The steady states of system (1) in absence of delay are obtained as non-negative solutions of the algebraic system As 0 ≤ P ≤ P and E = P − P , one has that P − E ≥ 0. Hence, when b − cη ≤ 0, it is immediately seen that system (2) has no solution (we look for interior steady states, so we exclude the solution E = P obtained for b − cη = 0). Let b − cη > 0. Then, the existence of fixed points for (1) is related to existence of solutions for the following second order equation in V A direct calculation shows that Eq. (3) possesses only one positive solution, independently of the sign of the coefficient of V.
If we want E ∞ > 0 and b − cη > 0 we obtain the following necessary condition: The study of local stability of equilibrium solution is based on the localization on the complex plane of the eigenvalues of the Jacobian matrix of (1) evaluated at (1) is locally asymptotically stable in absence of delay if and only if the following condition holds Proof. By Routh-Hurwitz's condition, all roots of Eq. (4) have negative real parts if and only if A > 0, B > 0, C > 0 and AB > C. Since A > 0 and C > 0, it is clear that AB > C implies B > 0. Expliciting the terms A, B and C in AB > C gives The statement follows by rewriting the previous inequality.
Consequently, if r ≥ δ, then condition (S) holds true and the fixed point is always stable.
In [16] the claim (without proof) is that AB > C and as a consequence the fixed point is always stable. However, we only were able to prove the case r ≥ δ leaving the case r < δ unsolved.
3. Stability analysis and Hopf bifurcation with τ > 0. The appearance of constant time delay does not affect equilibria. Hence, under the above conditions (V ∞ , E ∞ , K ∞ ) is still the unique positive equilibrium of system (1). In this case the associated characteristic equation of the linearization of system (1) at the equilibrium is given by i.e. where The stability of the equilibrium point (V ∞ , E ∞ , K ∞ ) will change when the system under consideration has zero or a pair of purely imaginary eigenvalues. The former occurs when λ = 0 in (6), i.e. if a 0 + b 0 = 2aδrV ∞ − (1 − η) rm 2 + (b − cη) δm 1 = 0, but we know this value to be positive from (5). Thus, λ = 0 is not a root of the characteristic equation (6). The latter deals with assuming that (6) has a purely imaginary root λ = iω, with ω > 0. Then, it follows from (6) that By separating real and imaginary parts in (7), we obtain Adding up the squares of Eqs. (8), we obtain Setting z = ω 2 , Eq. (9) becomes where Lemma 4. 1) Let q ≥ 0 and s ≥ 0. Then Eq. (10) has no positive roots. 2) Let q ≥ 0 and s < 0 or q < 0 and s ≤ 0. Then Eq. (10) has one positive root, say z 0 . 3) Let q < 0 and s > 0. Set Then Eq. (10) has no positive roots if h(z * ) > 0, Eq. (10) has one positive root z * if h(z * ) = 0, and it has two positive roots, say z − and z + , Proof. 1) is straightforward. For 2) and 3), we apply Descartes' rule of signs that states that the possible number of positive roots (10) is equal to the number of sign changes in the sequence formed by the polynomial's coefficients, or less than the sign changes by a multiple of 2. Finally, notice that h(0) = s, h(+∞) = +∞, h (z) = 3z 2 + 2pz + q, h (z) = 6z + 2p > 0; h(z) is a convex function, and it has a unique minimum z * = (−p + p 2 − 3q)/3 if q < 0 in the positive semiaxes.

TOMÁS CARABALLO, RENATO COLUCCI AND LUCA GUERRINI
Using (6), a direct calculation yields This completes the proof. In case Eq. (6) has two positive roots ω − , ω + , with ω − < ω + , then crossing from left to right with increasing τ occurs whenever τ assumes a value corresponding to ω + , and crossing from right to left occurs for values of τ corresponding to ω − . This implies that there exists a finite number of delayed intervals in which the equilibrium point is locally asymptotically stable, while unstable for the outside of the delayed ranges. Hence, the system dynamics switches from stable to unstable, and then back to stable when delay increases and crosses the critical delayed values.
Based on the transversality conditions (16) and the Hopf bifurcation theorem (see [12]), we have the following results on the stability of the equilibrium of system (1).
Proof. The proof of the theorem follows from Lemma 4 and Proposition 5. The expressions of the bifurcating values for the delays is given by (14) and (15).

Remark 8.
We note that the bifurcation analysis (together with the study of roots of eq. (10)) we propose here is more detailed and different with respect to that contained in [16].
4. Numerical simulations. In this section we will provide some numerical experiments to illustrate the results of the previous sections.
The fixed point is: In order to check the stability we need to compute the following quantities: Then, the stability condition is not fulfilled because The solution is represented in Fig. 1 below: the fixed point is not stable, the solution diverges to −∞.
The fixed point is: In order to check the stability we need to compute the following quantities: Then, the stability condition is not fulfilled either because The solution is represented in figures 2 below: the fixed point is not stable, there exists an attractive limit cycle.  thus, the fixed point is The initial data has been chosen as E, V, K = (1, 1, 1) for t ∈ [−τ, 0], while the critical value of the delay parameter for stability is τ 0 0 = 9.7592. The hypotheses of 2) in theorem 7 are fulfilled, then we expect that the fixed point (V ∞ , E ∞ , K ∞ ) is stable for τ < τ 0 0 and unstable otherwise. In Fig. 4 below we represent the solutions for τ = 8.6 and τ = 9.8 respectively: the fixed point looses stability and the solution diverges to −∞.  5. Comments on the model. In this section we would like to point out some features of the model that could be improved. In particular, the quantities V, E, K are positive and the system should preserve their sign. The system is well posed if for any positive initial conditions the solution remains positive, moreover by the modelisation strategy it is required that E ∈ [0,P ]. Suppose that we start with a positive initial condition V 0 , K 0 , E 0 > 0. Then if V or K become zero thenV > 0 (respK > 0). A contradiction of the requirement that E ∈ [0,P ] may occur when E reaches the value 0 orP .
In fact if there exists t 1 such that E(t 1 ) =P , thenĖ(t 1 ) > 0 if b − cη < 0 and E takes values bigger thanP . Then we need: Conversely, if there exists t 2 such that E(t 2 ) = 0 theṅ we have thatĖ(t 2 ) < 0 and E takes values less than 0. A sufficient condition (not necessary but that can be considered in the modelisation process) to avoid this should be b − cη < 0.
Of course these two conditions are not compatible. In order to preserve the positiveness of the system (and of its solutions) we consider a modified version of the second equation: It is worth highlighting that, for the original model, even if the condition of existence of the stable positive fixed point is fulfilled, the solutions can become negative. We provide a simulation of this problem in figure 6  This example suggests that a modification of the model should be necessary, in particular our suggested modification solves this problem.
Further research concerning the analysis of system (19) and a comparison with the results about the original system will be carried out somewhere else. 6. Further development: Non autonomous dynamics and some open directions. We have seen that the introduction of a delay may affect the stability of the interior fixed point. Moreover, Theorem 7 provides a very complex and varied picture. We have already observed that the stability switches introduce a serious difficulty for the managers of PAs: if they act on the delay and let it change, then it can bring both stabilizing and destabilizing effects. We consider that is very important to deal with this problem. However, in many situations it seems more realistic to consider a non-constant delay, in other words, instead of a constant delay τ > 0, we can consider a time dependent delay τ (t). More precisely, one could consider the following variation of (1), and study its asymptotic dynamics: where τ (·) is a continuous function (although in many real applications it uses to be continuously differentiable) from R to [0, h] where h > 0 is the maximum admissible delay.
The variable delay case in (20) possesses the same equilibria as the constant delay case, however the problem becomes automatically non-autonomous and in order to analyze the non-autonomous dynamical system generated by (20), it first needs to be stated as an ordinary differential system of equations in an infinite-dimensional space. Indeed, let us denote by C the space Banach space C 0 ([−h, 0]; R 3 ), and for any continuous function g : R → R 3 and t ∈ R we denote by g t the function in C defined by g t (s) = g(t + s), for all s ∈ [−h, 0].
Then, problem (20) can be rewritten as follows. For a given initial time t 0 ∈ R and an initial function φ ∈ C, we consider where y(t) = (E(t), K(t), V (t)) and f (·, ·) is the function defined from R × C by Then, when we particularize f (t, φ) for φ = y t , we obtain   which is our delay system. The asymptotic behavior of this delay system can be carried out in a local or in a global way. For the local analysis, there exists several well-known methods based, for instance, in the construction of certain Lyapunov functionals (see, Kolmanovskii and Shaikhet [14], Caraballo et al. [3], Kuang [15] amongst many others) as well as the use of the Razumikhin-Lyapunov theory which allows for more generality on the variable delay functions (see Hale and Lunel [11], Caraballo et al. [4], etc). We plan to investigate this in a future paper. However, concerning the global asymptotic behavior, as our problem will generate a non-autonomous dynamical system, then we need to choose an appropriate framework for our problem. To this respect we have several options, either the skew-product flow formalism (see e.g. Sell [17]), or the theory of uniform attractor (see Carvalho et al. [5], Kloeden and Rasmussen [13] and the references therein), or the theory of pullback attractor which usually allows for more generality on the delay terms. This problem will be analyzed in a subsequent paper.
We think that this should be the right context to study the consequence of modifying the delay and should be of a real interest for management of PAs to completely clarify this problem. However, we will leave this topic for a future work.
Finally, we would like to mention that there are some other interesting questions which we have not analyzed in this paper as, for instance, the global stability of the positive fixed point and the case sign {h (z * )} = 0 for which the transversality condition does not hold. Another question is to investigate the possible presence of a double Hopf bifurcations that may occur when two critical delays τ − j and τ + j are identical. We plan to investigate these issues in future.