A Note on the Tikhonov Theorem on an Infinite Interval

In this note we provide a new proof of the Tikhonov theorem for the infinite time interval and discuss some of its applications.


Introduction
Modern modelling dynamical processes with ordinary differential equations usually leads to very large and complex systems with the coefficients that often widely differ in magnitude. These features make any robust analysis of them close to impossible. In particular, the presence of very small and very large coefficients creates a stiffness in the system that renders standard numerical methods unreliable. At the same time, the presence of such coefficients indicates that the process is driven by mechanisms acting on very different time scales. Then one can hope that there is a dominant time scale; that is, the time scale at which the system, obtained by an appropriate aggregation of much faster and/or much slower processes, will have the same main dynamical features as the original one.
The presence of different time scales in a system is revealed if the nondimensionalization with respect to the chosen reference time unit produces coefficients that are significantly larger (or smaller) than the others. In this paper we will be dealing with systems that can be written in the so-called canonical, or Tikhonov, form where ,t denotes the time derivative, f and g are sufficiently smooth functions acting from an open subset of R n × R m into, respectively, R n and R m , and is a small parameter. As we shall see below, many more complex systems can be brought to such a form, see e.g. [4] for a systematic approach to a class of such systems. The interpretation of (1) is that the processes described by g happen much faster than those described by f and thus, if we are interested in larger times, it is plausible to assume that the former reach an equilibrium before any significant change occurs in the latter. Hence, for small , the solution to (1) should be close to the pair consisting of the solution t →v(t) = φ(u, t) to the algebraic equation g(u, v, t, 0) = 0, (2) called the quasi steady state, and the solutionū(t) of the reduced equation The validity of such an approximation is the subject of the Tikhonov theorem, see e.g. [14,15], that was also proven by other methods, such as the Center Manifold Theory, see e.g. [5,7]. The two main assumptions of the Tikhonov theorem are that: a) (2) admits a single isolated solution and, noting that for any fixed (u, t), φ(u, t) is an equilibrium of the initial layer equation v ,τ = g(u, v, t, 0), where τ = t/ but (u, t) are treated as parameters, and b) φ(u, t) is uniformly asymptotically stable in (u, t). The latter assumption is a mathematical expression of the fact that the quasi steady state is indeed practically reached in the fast time t/ ; that is, almost immediately in terms of the slow time t.
A problem with the Tikhonov theorem is that the adopted assumptions only ensure that the convergence is valid on finite intervals of time t and thus the approximation (2), (3) is useless if one wants to investigate the long term dynamics of (1). In other words, within the framework of the Tikhonov theorem, one cannot substitute the quasi steady state into (1) and draw any valid conclusions about the long term dynamics of (1) from the resulting reduced equation. We shall present an example of such a situation in Section 2.
This problem was first addressed in [8] where, assuming additionally that the relevant equilibrium of (3) is uniformly asymptotically stable, the author used the reverse Lyapunov theorem of [11] to construct appropriate Lyapunov functions to push the estimates of the original Tikhonov's proof to infinity. Recently the problem was again picked in [10], where the authors used the ideas of [5] to localize the equations along the quasi steady state and then proved the uniform in time estimates by directly using differential inequalities. The tricky part of this method is that the localization must preserve the properties of g that allow for solving (2) and keep the stability of the localized version of (3). Moreover, the localization must be extended close to t = 0 to a funnel-like domain to encompass initial conditionsẙ that may be far away from the quasi steady state. These, together with the specific form of localization, forced the authors of [10] to consider a restricted form of g which is one-dimensional with a dominant constant coefficient linear part.
The main aim of this paper is to address the restrictions mentioned above. We follow the ideas of [10] but use a different form of the localization that preserves the linear part of f and g. Moreover, to avoid extending the localization to the funnel-like domain (which requires an additional assumption), we use the estimates of the original proof of the Tikhonov theorem close to t = 0 and then employ the differential inequalities only for small initial conditions.
In the present paper we only prove the Tikhonov type result without addressing the order of the convergence, as was done in [10]. The higher order estimates, that can be done using the same approach, are the subject of the follow-up paper [3]. We also note that our assumptions, while being in line with that of [10], are stronger than in the original Tikhonov theorem but they are satisfied in most applications and make the proofs more straightforward.

Example
In this section, following [4,Section 4.3] and [1], we present an example showing that the Tikhonov approximation, being valid on finite time intervals, may not provide any reliable information on the long term dynamics of the original system.

Example 1
We assume that we have the populations of prey and predators, where the prey can move between two locations, say, grazing grounds and some refuge, while the predators only hunt in the grazing area. The interactions between the predators and the prey are modelled by the mass-action law. We denote by (n 1 , n 2 ) and p, respectively, the prey and the predator populations, and assume that the migrations are fast if compared to the vital processes. Then where, for i = 1, 2, n i denotes the prey density in patch i, m i denotes the migration rate from patch i, r i is the prey population growth rate in patch i, d is the predator death rate, a is the predation rate in patch 1 and b is the biomass conversion rate. We note that (4) is not in the typical Tikhonov form as letting = 0 in the first two equations yields two identical equation and the assumptions of the Tikhonov theorem are not satisfied. However, adding the first two equation and introducing the total prey population n := n 1 + n 2 , we obtain n ,t = n(r 1 − ap) + n 2 (r 2 − r 1 + ap), n 2,t = n 2 r 2 + (m 1 n − n 2 (m 1 + m 2 )), Denoting M i = m i m 1 +m 2 , i = 1, 2, we get the quasi steady staten 2 = M 1n and the reduced systemn which we recognize as the Lotka-Volterra model with the aggregated birth rate for the preȳ r = M 2 r 1 +M 1 r 2 and similarly adjusted predation and biomass conversion rates. We see that the assumptions of the Tikhonov theorem are satisfied as the quasi steady staten 2 = M 1n is a uniformly asymptotically stable equilibrium of the initial layer equation Thus the solution (n,p) (augmented by the initial layer) approximates the solution (n 1 + n 2 , p) of (4) on finite time intervals. On the other hand, the equilibria of (4) are (0, 0, 0) and, for small , The Jacobi matrix of (4), evaluated at (n * 1 , n * 2 , p * ), gives Denoting α = m 1 m 2 /(m 2 − r 2 ), β = ad/b, γ = bp * and using αm 2 − m 1 m 2 = αr 2 , we get the characteristic equation For small > 0 all coefficients are positive and hence e.g. the Hurwitz criterion ensures that real parts of all eigenvalues of J * are negative. Thus the positive equilibrium of the system (4) is asymptotically stable. This is in contrast to (5), for which the positive equilibrium is a centre.

Notation and Assumptions
As mentioned in the introduction, we consider an n × m dimensional system To be consistent with writing systems of equations in the column form, we adopt the convention that any vector x is a column vector and thus for a function t → x(t), x ,t is a column vector. Then, for a scalar function x → h(x), h ,x is the row vector of the first derivatives of h and, for an R m valued vector function . . , m, as its rows, or The estimates for the nonlinear problems depend to large extent on the estimates for their time dependent linearizations. For this we recall some relevant results. Consider an r × r system on R + , where A is a continuous matrix function, and let R + t → Y A (t) be the fundamental matrix for (7). We say that A satisfies the exponential dichotomy property if there are positive constants K, α such that We note that this is a simplified case of the exponential dichotomy discussed in e.g.
(A1) We assume that f : R n+m × R + × I e → R n and g : , e > 0, are C 3 functions with respect to u and v and C 1 with respect to t and , that are bounded together with all existing derivatives on [0, ∞) uniformly for u, v in bounded subsets of R n+m and ∈ I e .
As in the classical Tikhonov theorem, we assume that (A2) 0 = g(u, v, t, 0) admits an isolated solution v(t) = φ(u, t) for any (u, t) ∈ R n × [0, ∞) and we denote byū the solution to We assume that is a bounded differentiable function on [0, ∞).
For any matrix A we denote by σ (A) the spectrum of A and by s(A) := max{ λ; λ ∈ σ (A)} the spectral bound of A.
Finally, we adopt assumption that will ensure that the estimates are valid uniformly on R + .

Remark 3 Direct verification of
, being a constant matrix, has the exponential dichotomy property. Then, arguing as in Remark 2, J f (ū(t)) also has the exponential dichotomy property.

Error Estimates
As in the classical proof of the Tikhonov theorem, the estimates are split into estimates in the initial layer and in the bulk part. The first part is done as in [

Estimates in the Initial Layer
Let (u , v ) be the solution to (6).
Note that the solution of the last equation isv 0 (τ ). Thus there is ρ such that for any 0 < < ρ and τ ∈ [0, τ ρ ] we have or, returning to the original variable, uniformly for t ∈ [0, τ ρ ]. Using (19a) and the continuity of φ, we may take 0 small enough to ensure that Next, sinceū clearly is a continuous function, for sufficiently small we have ū(t ρ ) −ů ≤ 5ρ 6 and hence, by (19a), Using again (19a) and the continuity of φ, for sufficiently small we have Corollary 1 For any ρ > 0 there are ρ and τ ρ , such that for t ∈ [0, τ ρ ] and < ρ .
Proof Using (15) as the middle term on the right hand side of the inequality is 0 by (15).

ū(t), φ(ū(t), t), t, ).
Here H # g = H g + J g , where H g is the second order reminder of the expansion of g(u, v, t, 0) with respect to (u, v) aroundῡ(t) (see e.g. [2, Proof of Theorem 4.17]) and J g the first order reminder of the expansion of g(u, v, t, ) in around = 0. In particular, by op. cit, H g is of order of δ 2 , while J g is of order of .

), η(t),ū(t), φ(ū(t), t), t, )
where, as before, H # f = H f + J , H f is the second order reminder with respect to ζ and η evaluated at = 0, while J f is the first order reminder with respect to . Thus, M i , i = 1, 2, 3, depend only on the derivatives of f , g in E δ (ῡ(t)) up to second order and are finite irrespectively of u, v, , t; in particular, Then we consider the localized problem for some t 0 ≥ 0. We have the following lemma.
b) The equation forv can be written as where, by the definition ofg, h (t) ≤ M for some M independent of t and . Hence, by Similarly, by (28) and the definition of ζ , we find that Thus, using assumption (A5), for some constants C 1 , where η ∞ = sup 0< < ρ , t≥t 0 η (t) < ∞. Thus,ȗ (t) is bounded uniformly in t and .
Proof Using Lemma 1 with arbitrary ρ <δ 2 we consider (29) on [t ρ , ∞); that is, with the initial conditionsȗ(t ρ ) = u (t ρ ) andv(t ρ ) = v (t ρ ) for arbitrary fixed < ρ . The initial conditions belong to T ρ ⊂ Tδ /2 . As the first step, we consider a modified approximation foȓ v whose error is defined by whereȗ is the exact solution. We have, where v * is some point betweenv and the approximation, * is an intermediate point between 0 and and we usedg φ(ȗ (t), t), t, 0) ≡ 0.
(32) Before we move on, we observe that by, say [9, Theorem 2.6], the exponential dichotomy (16) is satisfied in some T δ , possibly with different constants K 2 , α 2 . Then, noting that we can decrease the tube Tδ without changing the constants (that can be left as they were for the larger set), we take ρ, ρ andδ small enough that for < ρ that is,ζ stays in the region wheref ,u +f ,vφ,u has the exponential dichotomy property.

An Application to Derive an Allee Type Dynamics
A population displays the so-called Allee type dynamics if it has some carrying capacity to which it monotonically increases if it is large enough but goes extinct if it is too small. Mathematically, the equation describing the evolution of the population should have three equilibria: asymptotically stable 0 as the extinction equilibrium, the repelling threshold equilibrium and the attractive carrying capacity. One of the ways to derive equations of this type is to look at populations interacting with each other according to the mass action law and exploiting multiple time scales occurring in such models. We present an example introduced in [13] and further analysed in [4]. In this model we consider a population N of females subdivided into subpopulations N 1 of females who recently have mated and N 2 of females who are searching for a mate. We assume that there is an equal number N of males. If the females reproduce in a very short time after mating, then where β is the reproduction rate of the recently mated females, μ + νN is the mortality rate of the recently mated females, μ + λ + νN denotes the (increased) mortality rate of the searching females, σ denotes the rate at which the females switch from the reproductive stage to the searching stage and ξN denotes the per capita rate at which a searching female finds one out of N potential mates.
To nondimensionalize the system, we rescale the time in units of the natural life expectancy 1/μ, s = μt, and, assuming β − μ > 0, we introduce the carrying capacity K = β − μ ν and setting N 1 = xK, N 2 = yK and z = x + y, we obtain our system in dimensionless form, where s is the rescaled time. Let us denote ε = μ σ ; that is, ε is the ratio of the average time of satiation to the average life span. In many cases it is a small parameter (for instance, for wild rabbits the average lifespan is 4 years and they breed 6-7 times per year, giving ≈ 0.04). For the population not to become extinct, we can argue that the rate at which a female finds a mate should be comparable with the rate she switches to a searching mode after reproduction, see [4] for a discussion of other cases. Thus, writing ξ/μ =ξ/ and denoting R 0 = β μ , we consider The right-hand side of the first equation can be simplified as We obtain the quasi steady state as and hence the reduced equation is given bȳ The auxiliary equation (13) here takes the form dŷ dτ = −ŷ(1 +ξKz) +z.
We are interested inz ∈ [0, ∞). Then the only equilibrium of (39) isŷ = φ(z) ≥ 0 and the right-hand side is decreasing forz > −1/ξK. So, (11) is satisfied. Further, by the above, any nonnegativeẙ belongs to the domain of attraction of the equilibrium solution.