Properties of the Lindemann Mechanism in Phase Space

We study the planar and scalar reductions of the nonlinear Lindemann mechanism of unimolecular decay. First, we establish that the origin, a degenerate critical point, is globally asymptotically stable. Second, we prove there is a unique scalar solution (the slow manifold) between the horizontal and vertical isoclines. Third, we determine the concavity of all scalar solutions in the nonnegative quadrant. Fourth, we establish that each scalar solution is a centre manifold at the origin given by a Taylor series. Moreover, we develop the leading-order behaviour of all planar solutions as time tends to infinity. Finally, we determine the asymptotic behaviour of the slow manifold at infinity by showing that it is a unique centre manifold for a fixed point at infinity.

1. Introduction. A unimolecular reaction occurs when a single molecule undergoes a chemical change. For unimolecular decay (or isomerization) to occur, a certain amount of energy must be supplied externally, namely the activation energy. For some time, there was debate concerning just how the molecules became activated. Frederick Lindemann suggested [15] in 1922 that unimolecular decay involves two steps, namely the activation/deactivation by collision step and the reaction step. Cyril Norman Hinshelwood made further contributions [12] to the Lindemann model in 1926 and, consequently, the Lindemann mechanism is occasionally referred to as the Lindemann-Hinshelwood mechanism. For general references on unimolecular reactions and the Lindemann mechanism, see, for example, [2,8,11,18].
Suppose that the reactant A is to decay into the product P . Then, according to the (nonlinear, selfactivation) Lindemann mechanism, A is activated by a collision with itself producing the activated complex B. This activation can also be reversed. The complex then decays into the product. Symbolically, (1.1)

Differential Equations and Common
Approximations. Using the Law of Mass Action, the concentrations of A and B in (1.1) satisfy the planar reduction where τ is time. The traditional initial conditions are a(0) = a 0 and b(0) = 0. However, we will allow the initial condition for b to be arbitrary. Note that dp/dτ = k 2 b and (traditionally) p(0) = 0. Since the differential equations (1.2) do not depend on the differential equation for p, we need only consider the differential equations for a and b.
There are two common approximations for the planar reduction. The Equilibrium Approximation (EA) and the Quasi-Steady-State Approximation (QSSA), which have proved successful for the Michaelis-Menten mechanism of an enzyme-substrate reaction [17], have also been applied to the Lindemann mechanism. See, for example, §2.2 of [11] and pages 122-126 and 313-317 of [18]. These approximations are frequently employed to simplify more complicated networks in chemical kinetics which may involve, for example, inhibition or cooperativity effects. For the EA, one assumes da/dτ ≈ 0 for sufficiently large time. This yields The QSSA, on the other hand, assumes db/dτ ≈ 0 for sufficiently large time. This yields b(τ ) ≈ k 1 a(τ ) 2 k 2 + k −1 a(τ ) .
It will be useful for us to convert the planar reduction to dimensionless form. Define t := k 2 τ, x := k 1 k 2 a, y := k 1 k 2 b, and ε := k −1 k 1 > 0, which are all dimensionless. Thus, t, x, and y are, respectively, a scaled time, reactant concentration, and complex concentration. Moreover, the parameter ε > 0 measures how slow the deactivation of the reactant is compared to the activation. Traditionally, one may want to consider ε to be small. In our analysis, the size of ε does not matter.
It is easy to verify that, with the above rescaling, the planar reduction (1.2) becomeṡ where˙= d/dt. Observe that the system (1.3) is a regular perturbation problem. Occasionally, we will need to refer to the vector field of this planar system. Hence, define g(x) := −x 2 + εxy x 2 − (1 + εx) y , (1.4) where x := (x, y) T . Moreover, we will be working with the scalar reduction y ′ = x 2 − (1 + εx) y −x 2 + εxy , (1.5) where ′ = d/dx, which describes solutions of the planar reduction (1.3) in the xy-plane by suppressing the dependence on time. We will need to refer to the right-hand side of the scalar reduction.
(1.6) Remark 1. 1. The function f (x, y) can be written f (x, y) = g 2 (x, y)/g 1 (x, y), where g 1 (x, y) and g 2 (x, y) are the components of the function g(x) given in (1.4). Note the use of the row vector (x, y) in the arguments of g 1 and g 2 as opposed to the column vector x. To alleviate notational headaches that arise from competing conventions involving row and column vectors, when there will be no confusion we will use the notation appropriate for the given situation.
1.2. Discussion. The Lindemann mechanism has been explored mathematically by others. For example, the planar system (1.3) has been treated as a perturbation problem in [22,24]. Furthermore, Simon Fraser has used the Lindemann mechanism [9,10] as an example in his work on the dynamical systems approach to chemical kinetics. Finally, properties of the Lindemann mechanism have been explored mathematically in [7].
The focus of this paper is the detailed behaviour of solutions to the planar reduction (1.3) in phase space. That is, we perform a careful phase-plane analysis to reveal important details that a common phase-plane analysis would miss. Equivalently, we are studying solutions of the scalar reduction (1.5). It is worth reiterating that our analysis does not depend on the size of the parameter ε (which is traditionally treated as being small). In §2, we present the basic phase portrait in the nonnegative quadrant. Moreover, we establish that the origin is a saddle node and is globally asymptotically stable with respect to the nonnegative quadrant. In §3, we describe the isocline structure which we exploit in later sections. For example, the isocline structure plays an important role in determining the concavity and asymptotic behaviour of solutions. In §4, we prove that there is a unique slow manifold M between the horizontal and vertical isoclines. To this end, we use a nonstandard version of the Antifunnel Theorem. In §5, we determine the concavity of all solutions, excluding the slow manifold, in the nonnegative quadrant by analyzing an auxiliary function. In §6, we use the Centre Manifold Theorem to show that all scalar solutions are given by a Taylor series at the origin. Moreover, we establish the leading-order behaviour of planar solutions as t → ∞. This is nontrivial due to the fact that the origin is a degenerate critical point. In §7, we show that all planar solutions must enter and remain in the region bounded by the horizontal isocline and the isocline for the slope of the slow manifold at infinity. In §8, we single out properties of the slow manifold. These properties include concavity, monotonicity, and asymptotic behaviour at the origin and at infinity. Finally, in §9, we state some open problems.
2. Phase Portrait. A computer-generated phase portrait for the planar reduction (1.3), restricted to the physically relevant and positively invariant nonnegative quadrant S, is given in Figure 2.1. In this paper, we will develop precise mathematical properties of the phase portrait. Equivalently, we develop results on solutions of the scalar reduction (1.5).
The horizontal and vertical isoclines for the planar system (1.3), which are found by respectively The QSSA corresponds to the horizontal isocline (the quasi-steady-state manifold) and the EA corresponds to the vertical isocline (the rapid equilibrium manifold). Observe that H(0) = 0 = V (0), both H and V are strictly increasing, and V (x) > H(x) for all x > 0. It appears from the phase portrait that the region between the isoclines, acts like a trapping region for (time-dependent) solutions of the planar reduction. Moreover, the origin appears to be globally asymptotically stable. Proof: (a) It follows from the definition (1.4) of the vector field g that g • ν < 0 along V and H, where ν is the outward unit normal vector. Thus, solutions cannot exit Γ 0 through the horizontal or vertical isoclines. Furthermore, solutions cannot escape from Γ 0 through the origin since solutions do not intersect. Hence, Γ 0 is positively invariant. (b) We will break the proof into cases.
Case 1: Case 2: x 0 > 0 and y 0 > V (x 0 ). Suppose, on the contrary, that x(t) does not enter Γ 0 . It follows that y(t) > V (x(t)) for all t ≥ 0. Using the differential equation (1.3), we knoẇ x(t) > 0 andẏ(t) < 0 for all t ≥ 0. Now, we see from the definition (1.6) of the function f thatẏ Note that εxy − x 2 > 0 since y > V (x) = x/ε and x, y > 0. Thus, Integrating with respect to s from 0 to t and rearranging, we obtain Let (x 1 , V (x 1 )) be the point of intersection of the vertical isocline y = V (x) and the straight line y = y 0 − (x − x 0 ). Obviously, x 1 > x 0 . Since x(t) is monotone increasing and bounded above by x 1 , we see that there is an x ∈ [x 0 , x 1 ] such that x(t) → x as t → ∞. Similarly, since y(t) is monotone decreasing and bounded below by V (x 0 ), we see that there is a y ∈ [V (x 0 ), y 0 ] such that y(t) → y as t → ∞. Thus, the ω-limit set is ω(x 0 , y 0 ) = {(x, y)}.
Since ω(x 0 , y 0 ) is invariant and (0, 0) is the only equilibrium point of the system, x = 0 and y = 0. This is a contradiction. Case 3: x 0 > 0 and 0 ≤ y 0 < H(x 0 ). This case is proved in a manner similar to Case 2. (c) If x 0 = 0, the solution of (1.3) is x(t) = (0, y 0 e −t ) T . This clearly satisfies x(t) → 0 as t → ∞. Thus, we can assume x 0 > 0 and, by virtue of part (b), we can assume further that (x 0 , y 0 ) ∈ Γ 0 . It follows from the differential equation (1.3) and the fact that Γ 0 is positively invariant thaṫ x(t) ≤ 0 andẏ(t) ≤ 0 for all t ≥ 0. Since both x(t) and y(t) are decreasing and bounded below by zero, by the Monotone Convergence Theorem we know that there are x and y such that x(t) → x and y(t) → y as t → ∞. Thus, the ω-limit set is ω(x 0 , y 0 ) = {(x, y)}. Since ω(x 0 , y 0 ) is invariant and (0, 0) is the only equilibrium point of the system, x = 0 and y = 0.
The Jacobian matrix at the origin for the planar system (1.3) is diag (0, −1). Thus, the origin is a nonhyperbolic fixed point. The Hartman-Grobman Theorem, unfortunately, cannot be applied here. Using Theorem 65 in §9.21 of [1], the origin is a saddle node which consists of two hyperbolic sectors and one parabolic sector. As we will effectively show later, S is contained in the parabolic sector. 3. The Isocline Structure. The horizontal and vertical isoclines, along with all isoclines between them, will be very useful. If we solve f (x, y) = c for y, we obtain y = F (x, c), where That is, y = F (x, c) is the isocline for slope c.  (i) The interior of the region Γ 0 corresponds to 0 < c < ∞ and 0 < K(c) < 1.
(ii) Two exceptional isoclines are y = V (x) (the vertical isocline) and y = 0 which correspond, respectively, to Furthermore, w is concave up for all x > −ε −1 K(c) and satisfies the differential equation Sketch of the isocline structure of (1.5). The isoclines above the vertical isocline have zero slope along the line y = 2x/ε.

Proof:
The proof is straight-forward and omitted.
Remark 3.3. The vertical isocline satisfies the limit (3.3) and the differential equation (3.4). The isocline w(x) = 0 (the isocline for slope −1) also satisfies the differential equation but does not satisfy the limit.
The isocline structure is sketched in Figure 3.2. We will often appeal to the isocline structure. For example, if a scalar solution of (1.5) is above the line y = 0 and below the horizontal isocline y = H(x), we know that −1 < y ′ (x) < 0.

Existence and Uniqueness of the Slow
Manifold. It appears from the given phase portrait, Figure 2.1, that there exists a unique solution to (1.5) that lies entirely in the region between the horizontal and vertical isoclines. To prove this, we will need to use a nonstandard version of the Antifunnel Theorem. See, for example, Chapters 1 and 4 of [13].
and consider the first-order differential equation The curves α and β satisfying (4.1) are, respectively, a lower fence and an upper fence. If there is always a strict inequality in (4.1), the fences are strong. Otherwise, the fences are weak.
Theorem 4.2 (Antifunnel Theorem, p.196 of [13]). Let Γ be an antifunnel with strong lower and upper fences α and β, respectively, for the differential equation Then, there exists a unique solution y(x) to the differential equation

Existence-Uniqueness Theorem.
We want to show that there is a unique scalar solution that lies entirely in the region Γ 0 . However, the vertical isocline is not a strong lower fence since "f (x, V (x)) = ∞." This turns out to be a fortunate obstacle.
Suppose that we want the isocline to be a strong lower fence for the differential equation (1.5) for all x > 0. Note that the condition on r restricts the isocline to being between the horizontal and vertical isoclines. Note also that Since w is concave up and satisfies the limit Since we want the isocline that will give us the thinnest antifunnel, we choose Note that α(x) is the isocline for slope ε −1 and Hence, define the region Proof: (a) We have already established that α is a strong lower fence. To show that H is a strong upper fence, observe Rearranging, we have Thus, using the definition (1.6) of f we have Hence, we can apply the Antifunnel Theorem with r(x) := ε 2 . To see why, consider that We can therefore conclude that there is a unique solution y = M(x) to (1.5) that lies in Γ 1 for all x > 0.
(b) Let y be a solution in Γ 0 lying below M. Since M is the only solution contained in Γ 1 , y must leave Γ 1 through the horizontal isocline. Now, let y be a solution in Γ 0 lying above M. Suppose on the contrary that y never leaves Γ 0 and thus M( Since w is a strong lower fence, the proof of the first part of the theorem can be adapted to show that M is the only solution contained in the region between the isoclines H and w. Thus, there is an a > 0 such that y(a) = w(a).
for all x > a. This is a contradiction since y(x) > V (x) for sufficiently large x.
(i) We are referring to the unique solution between the horizontal and vertical isoclines as the slow manifold. However, all scalar solutions in Γ 0 are technically slow manifolds (and, as it turns out, centre manifolds). This is because, as functions of time, the solutions approach the origin in the slow direction.
To see why this is the case, suppose w( where a > 0. Then, we can extend y(x) and y ′ (x) to say y(0) = 0 and y ′ (0) = 0.
Proof: Observe that the Squeeze Theorem establishes y(0) = 0. Now, Thus, by the Squeeze Theorem again as well as the definition of (right) derivative, we can say y ′ (0) = 0.

Nested Antifunnels.
The region Γ 1 is the thinnest antifunnel for x > 0 with isoclines as boundaries. However, we can find thinner antifunnels than Γ 1 which are valid for different intervals. For an isocline w(x) := F (x, c), where 0 < c < ε −1 , to be a strong lower fence on an interval, we Note that 1 − εc > 0. We will quickly establish a few properties of ξ(c). See  (a) The function ξ(c) satisfies The function ξ(c) is analytic for all c ∈ 0, ε −1 . Furthermore, ξ(c) has analytic inverse ξ −1 (x) defined for all x > 0.
Proof: The proof is routine, tedious, and omitted. (a) The isocline w satisfies .
(b) The slow manifold satisfies Proof: This has two roots, one negative and one positive. The positive root is given by x = ξ(c) with ξ(c) as in (4.3). It is a routine matter to confirm that w ′ (x) < f (x, w(x)) when 0 < x < ξ(c) and that w ′ (x) > f (x, w(x)) when x > ξ(c). (b) It follows from the Antifunnel Theorem.

Concavity.
In this section, we will establish the concavity of all scalar solutions, except for the slow manifold, in the nonnegative quadrant. The concavity of the slow manifold will be established later. These results will be obtained by using an auxiliary function. Moreover, we will construct a curve of inflection points which approximates the slow manifold.

Establishing Concavity.
Let y be a solution to (1.5) and consider the function f given in (1.6). If we differentiate y ′ (x) = f (x, y(x)) and apply the Chain Rule, we obtain y ′′ (x) = p(x, y(x))h(x, y(x)), where p(x, y) := 1 The function p(x, y) is positive everywhere except along the vertical isocline and for x = 0, where it is undefined. We will be considering the functions h(x, y) and p(x, y) along a given solution y(x) so we will abuse notation by writing h(x) := h(x, y(x)) and p(x) := p(x, y(x)).
For a given x > 0 with y(x) = V (x), it follows from (5.1) and the fact that p(x) > 0 that the sign of h(x) is the same as the sign of y ′′ (x). Furthermore, if we differentiate h(x) with respect to x and apply (5.1), we see that the function h has derivative

Region
Concavity of Solutions 0 ≤ y ≤ H concave down H < y < M concave up, then inflection point, then concave down M < y < V concave up y > V concave up, then inflection point, then concave down Table 5.1 A summary of the concavity of solutions of (1.5) in the nonnegative quadrant.
Remark 5.1. The function h cannot tell us anything about the concavity of solutions at x = 0, not even by taking a limit.
Claim 5.2. Let y be a solution to (1.5) and let x 0 > 0 with y(x 0 ) = V (x 0 ). Consider the isocline through the point (x 0 , y(x 0 )), which is given by w(x) := F (x, y ′ (x 0 )). Then, Proof: The first part follows from (3.4) and (5.2). The second part follows from the first.
The concavity of all solutions in all regions of the nonnegative quadrant can be deduced using the auxiliary function h and the following easy-to-verify lemma. Table 5.1 summarizes what we will develop in this section.  H(a). Then, there is a unique x 1 ∈ (0, a) such that y ′′ (x 1 ) = 0. Moreover, y is concave up on (0, x 1 ) and concave down on (x 1 , a]. Proof: Let h be as in (5.2) defined with respect to the solution y. Now, we know y ′ (a) = 0 and, by Proposition 4.6, we can extend y ′ (x) continuously and write y ′ (0) = 0. By Rolle's Theorem, there is an x 1 ∈ (0, a) such that y ′′ (x 1 ) = 0 and hence h(x 1 ) = 0. To show the uniqueness of x 1 , suppose that x 2 ∈ (0, a) is such that h(x 2 ) = 0. Now, since H(x 2 ) < y(x 2 ) < α(x 2 ), by virtue of the isocline structure 0 < y ′ (x 2 ) < ε −1 . Moreover, we see from (5.3) that h ′ (x 2 ) < 0. By Lemma 5.3, we can conclude x 2 = x 1 . Finally, by continuity we can conclude that h(x) > 0 on (0, Since since y(a−) = V (a), we know that J = ∅. We will show separately that h(x) > 0 for each x ∈ I and for each x ∈ J.
If I = ∅ then we are done. Suppose that I = ∅ and let x 3 ∈ J be fixed. We have already shown that h(x 3 ) > 0. We need to show that h(x) > 0 for all x ∈ I. By continuity, it suffices to show that h has no zeros in I. Suppose, on the contrary, that this is not the case. By Lemma 5.3, h has a right-most zero x 1 ∈ I ⊂ (0, x 3 ]. Observe that y ′ (x 1 ) < ε −1 , which follows from the fact that y(x 1 ) < α(x 1 ). Since h(x 1 ) = 0 and y(x 1 ) > 0, (5.3) informs us h ′ (x 1 ) < 0. Consequently, there is an x 2 ∈ (x 1 , x 3 ) such that h(x 2 ) < 0. Since h(x 3 ) > 0, the Intermediate Value Theorem implies that h has a zero in (x 2 , x 3 ), which is a contradiction since x 1 is the right-most zero.  The two thick curves are curves along which solutions of (1.5) have inflection points, for parameter value ε = 0.5. The thin curves are the horizontal, α, and vertical isoclines. We will show later that the slow manifold lies between the lower thick curve and the middle thin curve.
To prove uniqueness, suppose that x 2 ∈ (0, a) is such that h(x 2 ) = 0, where h is as in (5.2) defined with respect to the solution y. Now, we know from the isocline structure that y ′ (x 2 ) < ε −1 . Using (5.3), h ′ (x 2 ) < 0. Since x 2 was an arbitrary zero of h, we can conclude using Lemma 5.3 that x 2 = x 1 . Furthermore, since h ′ (x 1 ) < 0, we can say that y is concave up on (0, x 1 ) and concave down on (x 1 , a). Table 5.1 that solutions to the scalar differential equation (1.5) can only have inflection points between H and M or above V . We can construct a curve of inflection points, between H and M, which is close to the slow manifold.

It is easily verified that
where h is as in (5.2). Thus, there are three curves along which solutions have zero second derivative, given implicitly by One curve lies below the x-axis and is discarded. The other two curves, as expected, are in the positive quadrant. See Figure 5.1.
Recall that, for a fixed c ∈ 0, ε −1 , the isocline w(x) := F (x, c) switches from being a strong lower fence to being a strong upper fence at x = ξ(c) and y = F (ξ(c), c), where F is defined in (3.1) and ξ is defined in (4.3). As it turns out, will be a curve of inflection points between H and M. Note that H(x) < F (x, c) < α(x) for all x > 0 and c ∈ 0, ε −1 , which follows from the isocline structure. Moreover, note 0  F (x, c). Then, the isocline w satisfies .
Proof: Note that 0 < c < ε −1 and y 0 = w(x 0 ). We will only show the third case since the other two cases are similar. Assume that Y(x 0 ) < y 0 < α(x 0 ). Appealing to the isocline structure, we .
Proposition 5.11. The curve y = Y(x) satisfies Proof: We know already that H(x) < Y(x) < α(x) for all x > 0. We know from our results on concavity (see Table 5.1) that h(x, y) > 0 if x > 0 and M(x) < y < α(x), where h is the function defined in (5.2). By continuity, we can conclude h(x, M(x)) ≥ 0 for all x > 0. It follows from Proposition 5.10 that Y(x) ≤ M(x) for all x > 0.
To establish a strict inequality, let h be defined along the solution y = M(x). Assume, on the contrary, that there is an x 0 > 0 such that h(x 0 ) = 0. Using (5.3), h ′ (x 0 ) < 0. This contradicts the fact that h(x) ≥ 0 for all x > 0.

Slow Tangent Manifold.
The curve y = Y(x) can be referred to as a slow tangent manifold (or as an intrinsic low-dimensional manifold) since it consists of the points for which the tangent vector for the planar system (1.3) points in the slow direction. See, for example, [14,16,19,20,21,23]. To see why the curve of inflection points and the slow tangent manifold are equivalent, first consider the general planar systemẋ = g(x), where˙= d/dt and g ∈ C 1 R 2 , R 2 , along with the corresponding scalar system y ′ = g 2 (x, y)/g 1 (x, y), where ′ = d/dx. Consider the linearization matrix A(x) := {g ij (x)} 2 i,j=1 , where g ij (x) := ∂g i (x)/∂x j , which has characteristic equation λ 2 − τ λ + ∆ = 0, where τ := g 11 + g 22 and ∆ := g 11 g 22 − g 12 g 21 .
For notational brevity, we are suppressing the dependence on x.
Proof: The proof is routine.
Proposition 5.13. Suppose that g 1 = 0, g 12 = 0, and 4∆ < τ 2 at some fixed point (a, b) and let y(x) be the scalar solution through (a, b). Then, y ′′ (a) = 0 if and only if g v + or g v − at (a, b).

Proof:
First, note that σ ± is the slope of the eigenvector v ± . If we differentiate the scalar differential equation and manipulate the resulting expression, we obtain The conclusion follows.
Proposition 5.14. Consider the planar system (1.3) and the scalar system (1.5). The curve of inflection points y = Y(x) is a slow tangent manifold.
6. Behaviour of Solutions Near the Origin. In this section, we establish the full asymptotic behaviour of scalar solutions y(x) as x → 0 + . Moreover, we will obtain the leading-order behaviour of planar solutions x(t) as t → ∞.
6.1. Scalar Solutions. We will begin by attempting to find a Taylor series solution. Consider the differential equation (1.5), which can be rewritten Assume that y(x) is a solution in Γ 0 of the form for undetermined coefficients {b n } ∞ n=0 . If we substitute the series (6.2) into (6.1) and then solve for the coefficients, we obtain We will use centre manifold theory to show that the series (6.2) is fully correct for each solution inside the trapping region Γ 0 . However, we must first show that each solution is a centre manifold.
That is, we must show that each solution y(x) satisfies y(0) = 0 and y ′ (0) = 0. Proposition 4.6 already established that this is true for y(x) inside Γ 1 .
Proof: These limits have already been established if y is the slow manifold or if y lies below the slow manifold M. Hence, we will assume that M(x) < y(x) < V (x) for all x ∈ (0, a).
Let c := y ′ (a/2). We know from Table 5.1 that y is concave up on (0, a/2). Since M(x) > 0 for x ∈ (0, a/2), we thus have 0 < y(x) < F (x, c) for all x ∈ 0, 1 2 a , where F is the function given in (3.1). Note that It follows from the Squeeze Theorem that we can take y(0) = 0. Now, observe that Again by the Squeeze Theorem, we see that we can take y ′ (0) = 0. b n x n as x → 0 + .
Note that the Taylor coefficients of the series for u(x) must be {b n } ∞ n=2 since they are generated uniquely by the differential equation. Since y(x) is a centre manifold, it follows from centre manifold theory that for any k ∈ {2, 3, . . .}. See, for example, Theorem 1 on page 16, Theorem 3 on page 25, and properties (1) and (2) on page 28 of [5]. The conclusion of the theorem follows.
Remark 6.3. For analytic systems of ordinary differential equations for which the Centre Manifold Theorem applies, if the Taylor series for a centre manifold has a nonzero radius of convergence, then the centre manifold is unique. Since all solutions y to (1.5) lying inside Γ 0 are centre manifolds, we can conclude that the Taylor series ∞ n=2 b n x n has radius of convergence zero.

Planar Solutions.
We can use the isoclines to extract the leading-order behaviour of planar solutions as time tends to infinity.
Proof: Let c > 0 be fixed and arbitrary. We know from Theorem 2.1, Table 5.1, and the isocline structure that there exists a T ≥ 0 such that where F is given in (3.1). Using (1.3), (2.1), (3.1), and (6.4), we can see that for all t ≥ T, (6.5) where b := ε/K(c) and K is the function defined in (3.2). Note that the solution of the initial value problemu where a, u 0 > 0 and t 0 ≥ 0 are constants, is where W is the Lambert W function [6]. A simple comparison argument applied to (6.5) establishes ϕ(t; ε, T, x(T )) ≤ x(t) ≤ ϕ(t; b, T, x(T )) for all t ≥ T. (6.6) A standard property of the Lambert W function is W (t) = ln(t) − ln(ln(t)) + o(ln(ln(t))) as t → ∞.

Consequently, it can be shown
With a little manipulation, it can be verified that for any a, u 0 > 0 and t 0 ≥ 0. It follows from (6.6) that Since c was arbitrary with K(c) → 1 and b → ε as c → 0 + , This yields the desired conclusion for x(t). The conclusion for y(t) follows from Theorem 6.2.
It is possible to derive the expression for x(t) in Proposition 6.4 without appealing to the isocline structure and concavity. To achieve this, we will note that x(t) satisfies the integral equation and then twice utilize the following easy-to-verify lemma.
Lemma 6.5. Let a ∈ R be a constant and let f, g : [a, ∞) → R be nonnegative, integrable functions such that We know from Theorems 2.1 and 6.2 that It follows from Lemma 6.5 that t 0 y(s) x(s) By virtue of the integral equation (6.7), To take this one step further, observe now that which follows from Theorem 6.2, and so by Lemma 6.5 we have t 0 y(s) x(s) ds = ln(t) + o(ln(t)) as t → ∞.
By virtue of the integral equation (6.7) once again, 7. All Solutions Must Enter the Antifunnel. Earlier, in Theorem 2.1, we showed that all solutions x(t) to the planar system (1.3), except for the trivial solutions, eventually enter the trapping region Γ 0 . Here, we show that Γ 1 is itself a trapping region.
(a) There is a t * ≥ 0 such that x(t) ∈ Γ 1 for all t ≥ t * . (b) Define the region Then, there is a t * ≥ 0 such that x(t) ∈ Γ 2 for all t ≥ t * . Proof: (a) We know from Theorem 2.1 that x(t) eventually enters and stays in Γ 0 . Let y(x) be the corresponding scalar solution to (1.5). Then, we can say y ′ (0) = 0. Appealing to the isocline structure, this means that x(t) has entered Γ 1 . Furthermore, since g • ν < 0 along the horizontal and α isoclines which form the boundaries of the region in question, we see that Γ 1 is positively invariant. (b) It follows from Table 5.1, Proposition 5.10, and the previous part of the theorem.
8. Properties of the Slow Manifold. In this section, we will highlight some properties of the slow manifold.   Proof: It follows from Propositions 5.10 and 8.1 that h(x, M(x)) > 0 for all x > 0, where h is the function defined in (5.2). Since sgn(M ′′ (x)) = sgn(h(x, M(x))), it must be that the slow manifold is concave up for all x > 0. Proof: The first part is a consequence of Proposition 8.1 and the isocline structure. The first limit is a special case of Proposition 6.1. To prove the second limit, let c ∈ 0, ε −1 . It follows from Proposition 4.8 and the isocline structure that where ξ is the function defined in (4.3). Applying (4.4) and the Squeeze Theorem gives the second limit.
Remark 8.4. The justification which Fraser provides in [9] (just before Theorem 1) that M ′ (x) → ε −1 as x → ∞ is incorrect. The error is that the distance between the horizontal and vertical isoclines does not tend to zero as x tends to infinity. Thus, the asymptotic behaviour of M ′ (x) need not be the same as the asymptotic behaviour of H ′ (x) and V ′ (x). Proof: Since the slow manifold is contained entirely in Γ 0 , we can apply Theorem 6.2.
Corollary 8.6. The slow manifold satisfies Moreover, this statement would not be true if we replace H(x) with any other isocline F (x, c).
Proof: It follows from a comparison of the asymptotic expansions for M(x), H(x), and F (x, c). Now we will establish the full asymptotic behaviour of the slow manifold at infinity. First, we will extract as much information as possible from the isoclines. Second, we will attempt to find a series in integer powers of x. Third, we will prove definitively that the resulting series is indeed fully correct.
Let c ∈ 0, ε −1 . We know from Proposition 4.8 that F (x, c) < M(x) < α(x) for all x > ξ(c), Proof: To prove the result, we will apply the Centre Manifold Theorem to a fixed point at infinity. Consider the change of variables X := x −1 and Y := y − r(x), where r(x) := ρ −1 x + ρ 0 + ρ 1 x −1 , with the coefficients ρ −1 , ρ 0 , and ρ 1 being given in (8.3). Differentiate the new variables with respect to time and use the differential equation (1.3) to obtain the systeṁ X = −X 2 g 1 X −1 , r X −1 + Y , where g 1 and g 2 are as in (1.4). This system is not polynomial but there is no harm in considering the systemẊ = −X 3 g 1 X −1 , r X −1 + Y , Y = X −r ′ X −1 g 1 X −1 , r X −1 + Y + g 2 X −1 , r X −1 + Y , (8.4) which is polynomial. This is because the resulting scalar differential equation will be the same. The system at hand, while messy, is in the canonical form for the Centre Manifold Theorem. Note that the eigenvalues of the matrix for the linear part of this system are 0 and −(1 + ε). We know from centre manifold theory that there is a C ∞ centre manifold Y = C(X) which, we claim, must be the slow manifold. ρ n X n as X → 0 + , for some coefficients { ρ n } ∞ n=2 . Upon reverting back to original coordinates and observing that the coefficients in (8.3) are generated uniquely from the differential equation, the conclusion follows. Moreover, this statement would not be true if we replace α(x) with any other isocline F (x, c). 9. Open Questions. It would be nice to extend Proposition 6.4 to include more terms. In particular, it is desirable to have the lowest-order term which depends on the initial condition. For x(t), it is expected that the initial condition x 0 first appears in the 1/t 2 term since this is the case when ε = 0, which has x n−1 0 t n as t → ∞.
More generally, we would like to develop an iterative procedure to extract as many terms as possible from the asymptotic expansion of a solution of a nonlinear differential equation which approaches a degenerate critical point in the direction of a centre manifold.