## Abstract

The goal of this paper is to provide a complete picture of the long-term behavior of Lorenz’s maximum simplification equations along with the corresponding meteorological interpretation for all initial conditions and all values of the parameter.

## 1. Introduction

In this paper our goal is to analyze the properties of equilibrium points and long-term behavior of the solution of the Lorenz’s (1960) maximum simplification equations (MSE). Lorenz (1960) introduced the maximum simplification equations and demonstrated that they could be used to understand certain atmospheric phenomena, such as the index cycle. The dynamical equations were simplified to a point where reasonable predictions of atmospheric conditions could not be expected, however, certain physical processes and features of interest were retained. This MSE consists of a set of three nonlinear and coupled system of three ordinary differential equations (ODE) that govern the evolution of the amplitudes of three modes arising from the double Fourier series expansion of the solution of the barotropic vorticity equation (BVE). The importance of MSE stems from the fact that while it is the simplest possible (minimal set) approximation to BVE, using it, Lorenz (1960) was able to capture many of the qualitative features of the flow including the exchange of kinetic energy between the zonal and meridional components. Motivated by this success, the basic spirit of this ground-breaking work of Lorenz (1960) namely approximating the solution of complex realistic problems using a low-order model (LOM) has become a standard arsenal in problem solving in geophysical domain. Refer to Saltzman (1962), Lorenz (1963), De Swart (1989), Wiin-Nielson (1992), Gluhovsky et al. (2002), and the references therein for details.

Thompson (1957) used a set of equations similar to MSE to explore the limits of deterministic predictability that stem from the uncertainty in the initial conditions. Epstein (1969) used MSE to study the stochastic dynamic prediction. Recently, Lakshmivarahan et al. (2003) used MSE to compare the performance of first-order and second-order methods for data assimilation. To our knowledge, long-term behavior of Lorenz’s MSE has not been explored. In this paper we examine the growth of small errors superimposed on equilibrium points of MSE leading to a set of interesting heteroclinic orbits. The necessary conditions for obtaining a solution that approaches an equilibrium point are established. In addition, we provide a complete catalog of the behavior of MSE starting from any arbitrary initial point. The impacts of changing various model parameters on the solutions with respect to the index cycle are also examined.

In section 2, MSE and their key properties are derived. The set of all equilibria for MSE are identified and a set of necessary conditions for the orbits to approach an equilibrium point are derived in section 3. A comprehensive analysis of stability of equilibria is covered in section 4. A complete picture of the evolution of small errors superimposed on equilibria is given in section 5. A complete classification of the behavior of the solution starting from an arbitrary point is covered in section 6. Meteorological implications of the above analysis are given in section 7. The concluding section 8 provides a summary.

## 2. The maximum simplification equations—The model

Our model equations are the Lorenz’s so-called maximum simplification equations given by

where *A* = *A*(*t*), *F* = *F*(*t*), and *G* = *G*(*t*) are the time varying amplitudes of the three spectral components obtained from the standard BVE and *Ȧ*, *Ḟ*, and *Ġ* refer to their time derivatives. See Lorenz (1960) for details of the derivation. The streamfunction *ψ* = *ψ*(*x*, *y*, *t*) and the vorticity ∇^{2}*ψ* = ∇^{2}*ψ*(*x*, *y*, *t*) corresponding to (2.1) are given by

and

where *x* and *y* are the usual space coordinates, ∇^{2} is the standard Laplacian operator, 2*π/l = L* is the distance between zonal wind maximum, and 2*π/k* = *K* is the wavelength of the disturbance. That is, *K* and *L* denote the spatial scales of the flow and the disturbance. As an example, *L* = 10 000 km and *K* = 5000 km, that is, *l* = 6.285 × 10^{−4} and *k* = 1.257 × 10^{−3}. The expression for the easterly wind *u* = *u*(*x*, *y*, *t*) = −∂*ψ*/∂*y* and northerly wind *υ* = *υ*(*x*, *y*, *t*) = ∂*ψ*/∂*x* are given by

and

### a. Dependence on the parameter r

It can be verified that the model (2.1) depends only on the ratio *r* = *k*/*l* = *L*/*K* and can be written as

where

and

A plot of the values of *a* = *a*(*r*), *b* = *b*(*r*), and *c* = *c*(*r*) as a function of *r* > 0 is given in Fig. 1. It can be verified that *b* = *c* when *r* = 41/3 = 0.7598 and *a* = *c* when *r* = 43 = 1.3161. Further, *b* → 0^{+} (from above) as *r* → 0 and *a* → 0^{−} (from below) as *r* increases.

Let ** η** = (

*A*,

*F*,

*G*)

^{T}∈ ℝ

^{3}be a column (real) vector of three spectral components, where T denotes the transpose. That is,

*η*

_{1}=

*A*,

*η*

_{2}=

*F,*and

*η*

_{3}=

*G*. Then (2.6) can be succinctly written as

where **g**: ℝ^{3} → ℝ^{3} is such that **g**(** η**) = [

*g*

_{1}(

**),**

*η**g*

_{2}(

**),**

*η**g*

_{3}(

**)]**

*η*^{T}with

*g*

_{1}(

**) =**

*η**aFG*,

*g*

_{2}(

**) =**

*η**bAG,*and

*g*

_{3}(

**) =**

*η**cAF*.

It should be interesting to note that the Eq. (2.6) is equivalent to the Euler’s equation of motion for a rigid body and is known as the Euler gyrostat (Obukhov 1973). Derivation and analysis of Euler gyrostat can be found in textbooks on classical mechanics (Landau and Lifshitz 1960; Wittenburg 1977; Arnol’d 1989). Recently Gluhovsky and his collaborators (Gluhovsky and Agee 1997; Gluhovsky and Tong 1999; Gluhovsky et al. 2002) have shown that a Volterra gyrostat (of which the Euler gyrostat is a special case) can be used as a building block for generating a wide variety of low-order models of interest in several application domains.

Implicit in the above description are the two sets of coordinate systems. First is the three-dimensional physical space consisting of the two space variables, *x* and *y*, and the time variable *t*, in which the physical quantities such as the streamfunction *ψ*, vorticity ∇^{2}*ψ*, and the wind components *u* and *υ* are described. The second is the three-dimensional spectral space in which the spectral amplitudes *A*, *F*, and *G* evolve as functions of time according to (2.6) or (2.10). The goal is to characterize the long-term behavior of the physical quantities from those of the spectral components.

### b. Special case of r = 1

It can be verified that *a* = −½, b = ½ and *c* = 0 when *r* = 1. In this case, *Ġ* = 0 and *G*(*t*) = *G*(0) for all *t*. Hence, if *G*(0) ≠ 0, (2.6) reduces to a linear harmonic oscillator

whose solution is periodic (Hirsch and Smale 1974). When *G*(0) = 0 and *r* = 1, then ** η**(

*t*) ≡ 0. In the following it is tacitly assumed that

*r*≠ 1.

### c. Conservation of kinetic energy and vorticity

Let 𝗗 = diag(*d*_{1}, *d*_{2}, *d*_{3}) be a diagonal matrix where *d _{i}* > 0,

*i*= 1, 2, and 3. Consider the functional

*f*: ℝ

^{3}→ ℝ defined by

It can be verified that the time derivative of this functional evaluated along the trajectory of (2.6) is given by

Thus, given any parameter vector **p** = (*a*, *b*, *c*)^{T} of the model, for any **d** = (*d*_{1}, *d*_{2}, *d*_{3})^{T} that is orthogonal to **p**, *ḟ*(** η**) ≡ 0 and for this choice of the vector

**d**or matrix 𝗗,

*f*(

**) is conserved by the model. By direct computation it can be verified that the mean kinetic energy**

*η*and the mean square vorticity

are conserved by the model (2.1). That is, *E*(*t*) = *E*(0) and *V*(*t*) = *V*(0) for all *t* ≥ 0. Refer to Lorenz (1960) and Krishnamurti et al. (1998) for more detail.

### d. Solution space and degree of freedom

Let ** η**(0) = [

*A*(0),

*F*(0),

*G*(0)]

^{T}be the given initial condition for (2.10). Then (2.15) implies that the solution

**(**

*η**t*) of (2.10) is confined to the surface of the

*V*ellipsoid, which is centered at the origin and is given by

where

Clearly, *f _{V}*(

**) depends only on**

*η***(0) and not on**

*η**r*. Similarly, from (2.14) it follows that the solution

**(**

*η**t*) of (2.10) is confined to the surface of the

*E*ellipsoid, which is also centered at the origin, and is given by

where

Combining these, it follows that the solution ** η**(

*t*) of (2.10) must lie in the set

which is the intersection of the two concentric ellipsoids. That is, this set 𝒮 is the solution space for the model Eq. (2.10). Refer to Fig. 2 for an example of the two intersecting ellipsoids.

where **∇***f _{V}*(

**) = (1/**

*η**V*(0))(

*A*,

*F*, 2

*G*)

^{T}is the gradient of

*f*(

_{V}**) and**

*η*where **∇***f _{E}*(

**) = [1/2**

*η**E*(0)][(

*A*/

*l*

^{2}), (

*F*/

*k*

^{2}), (2

*G*/

*l*

^{2}+

*k*

^{2})]

^{T}is the gradient of

*f*(

_{E}**). It can be verified that the gradient vectors**

*η***∇**

*f*and

_{V}**∇**

*f*are linearly independent unless

_{E}*r*= 1 and

*G*= 0. Thus, the flow of the dynamics in (2.10) is always orthogonal to the plane defined by these two gradient vectors. Since

**∇**

*f*and

_{V}**∇**

*f*define the outward normal to the surface of the ellipsoids, the flow always lies along the tangent to the points in 𝒮 in (2.20).

_{E}Lorenz (1960) has shown that (2.6) can be solved analytically where the solution *A*(*t*), *F*(*t*), and *G*(*t*) are expressed as elliptic functions of time. In principle, we can use the tables of elliptic functions (Abramowitz and Stegun 1965) to characterize the long-term behavior of the solution of (2.6). However, following Lorenz (1960), we use the standard numerical methods to solve (2.6). While Lorenz (1960) used the central difference scheme that involves three time steps, we use the fourth-order adaptive Runge–Kutta (RK) method.

### e. Boundedness of the trajectories

An immediate consequence of the property described in section 2b is that the trajectories ** η**(

*t*) of (2.10) remains bounded. Let

**(0) = [**

*η**A*(0),

*F*(0),

*G*(0)]

^{T}be the initial state. Then from (2.14) and (2.15) we get

For example, when *A*(0) = 10^{−3} = *G*(0) and *F*(0) = 5, *K* = 5000 and *L* = 10 000, it can be verified that |*A*(*t*)| ≤ 2.5, |*F*(*t*)| ≤ 5.0, and |*G*(*t*)| ≤ 3.5355, for all *t* ≥ 0.

This property in turn implies that the model (2.1) is (externally) stable in the sense of Lyapunov (Hirsch and Smale 1974; Sparrow 1982). However, in view of the nonlinear interaction between modes, there is transfer of energy between modes leading to (internal) instability.

An example illustrating this transfer of energy between modes is given in Fig. 3 where the system (2.6) is solved using the fourth-order adaptive Runge–Kutta method. This figure provides a complete picture of the model trajectories, namely the variation of *A*(*t*), *F*(*t*), and *G*(*t*) with time; two dimensional plots of *A*(*t*) − *F*(*t*), *A*(*t*) − *G*(*t*) and *F*(*t*) − *G*(*t*), and the three-dimensional view of the solution.

## 3. Equilibrium points of the model

Recall that if = **g**(** η**) represents a dynamical system, then the set of points at which

**g**(

**) = 0 are called the equilibrium points of the system (Hirsch and Smale 1974; Devaney 1989; Sparrow 1982). Accordingly, we define three sets of equilibria for the system in (2.6). Let**

*η*It can be verified that *ɛ _{A}*,

*ɛ*and

_{F}*ɛ*constitute the set of all equilibrium points of (2.6). It can be seen that the origin is the only point common to all these three sets.

_{G}Given that the solution of the model (2.10) always lies in the set 𝒮 in (2.20), a question of interest at this point is: Which of the above equilibrium points lie arbitrarily close to the solution space 𝒮? Stated in other words, when does the model solution starting from an arbitrary initial point ever approach an equilibrium point?

In answering this question, first consider ** η*** = (

*A**, 0, 0)

^{T}∈

*ɛ*. It can be verified that

_{A}*** ∈ 𝒮**

*η*_{𝒜}when

that is, when

in which case *A** = ±2*V* = ±2*l**E*. Thus, (3.2) is a necessary condition for an orbit to include ** η*** ∈

*ɛ*.

_{A}Considering ** η*** = (0,

*F**, 0)

^{T}∈

*ɛ*and repeating the above argument, it can be verified that

_{F}is a necessary condition for an orbit to include ** η*** ∈

*ɛ*where

_{F}*F** = ±2

*V*= ±2

*k*

*E*. Similarly, considering

*** = (0, 0,**

*η**G**)

^{T}∈

*ɛ*, it readily follows that

_{G}is necessary for the orbit to include ** η*** ∈

*ɛ*, where

_{G}*G** = ±2

*V*= ±2(

*l*

^{2}+

*k*

^{2})

*E*.

Stated in other words, the key to understanding the behavior of orbits of the model in (2.10) is the value of the ratio *V*/*E*. Lorenz (1960) recognized the importance of *V*/*E* (= 2*k*^{2} as the threshold for sub and supercritical disturbances, for example), but did not explain why it was significant. We see here that the necessary conditions for obtaining any solution including an equilibrium point as part of the orbit are set by specific values of the ratio *V*/*E*.

## 4. Stability of equilibria

Recall that the local stability properties of an equilibrium point, say ** η***, is governed by the distribution of the eigenvalues of the Jacobian, 𝗗

**(**

_{g}**) of the mapping**

*η***g**(

**) in (2.10) evaluated at**

*η**** (Hirsch and Smale 1974; Devaney 1989; Lewis et al. 2006, chapter 32). The Jacobian of**

*η***g**(

**) is given by**

*η*### a. Stability of equilibrium points in *ɛ*_{A}

Let ** η*** ∈

*ɛ*Then, the Jacobian (4.1) becomes

_{A}.Let *r* < 1. Then, *a* < 0, *b* > 0, and *c* > 0. It can be verified that the eigenvalues of 𝗗** _{g}**(

***) are given by**

*η*and the corresponding eigenvectors are given by

For latter reference, define the diagonal matrix **Λ**^{(1)} of eigenvalues and 𝗩^{(1)}, the nonsingular matrix of the corresponding eigenvectors as

and

where (it follows from the definition)

Notice that while the eigenvalues *λ*^{(1)}_{2} and *λ*^{(1)}_{3} in (4.3) depend on *A* and the ratio *r*, the eigenvectors **v**^{(1)}_{2} and **v**^{(1)}_{3} are independent of *A*. That is, every point in *ɛ _{A}* shares the same set of eigenvectors that depends only on the value of the ratio

*r*. For purposes of comparison, the eigenvalues and vectors of 𝗗

**(**

_{g}***) when**

*η**** belongs to each of the three equilibria are summarized in Table 1.**

*η*Let *r* > 1. Then *a* < 0, *b* > 0, and *c* < 0. It can be verified that the complex eigenvalues of 𝗗** _{g}**(

***) are given by**

*η*and the corresponding complex eigenvectors are **v**^{(1)}_{1}, **w**^{(1)}_{2} = **u**_{1} + *i***u**_{2} and **w**^{(1)}_{3} = **u**_{1} − *i***u**_{2}, respectively, where

We leave it to the reader to verify that this eigen decomposition in the complex domain can be equivalently expressed in real domain using the standard canonical form (Hirsch and Smale 1974) as

where

and

That is, the eigenvectors that are the columns of 𝗩^{(1)} are given in (4.8). In this paper we use this latter form of eigen decomposition when the eigenvalues are complex. Refer to Table 1 for a summary.

and

Variation of *bc* and *c*/*b* as a function of *r* are given in Figs. 4a,b, respectively. It follows from (4.3)–(4.4) and (4.8)–(4.9) that the eigenvalues and eigenvectors depend on two factors: the value of the first coordinate *A* of the equilibrium point ** η*** = (

*A*, 0, 0)

^{T}and the value of the ratio

*r*. For completeness, we consider several subcases.

#### 1) Case 1: Let *A* > 0 and *r* < 1

In this case, it can be verified that *λ*^{(1)}_{2} > 0 and *λ*^{(1)}_{3} < 0, and the equilibrium point, ** η*** = (

*A*, 0, 0) is called a saddle point, with the eigenvector

**v**

^{(1)}

_{2}corresponding to the unstable (repelling) direction and the eigenvector

**v**

^{(1)}

_{3}to the stable (attracting) direction (Hirsch and Smale 1974). Refer to Table 2 for a summary. In verifying these claims, let

**=**

*η**** +**

*η***where**

*ξ***is a small perturbation. Then, using the first-order Taylor series expansion, we obtain**

*ξ*or

as the local dynamics of evolution of the perturbation. Now change the variables from the standard (*A*, *F*, *G*) orthogonal coordinate system to oblique coordinate system defined by the eigenvectors **v**^{(1)}_{1}, **v**^{(1)}_{2} and **v**^{(1)}_{3} in (4.6). Define a new vector **z** = (*z*_{1}, *z*_{2}, *z*_{3})^{T} where

That is,

and *z*_{2}(*t*) = *z*_{2}(0)*e*^{λ(1)t}_{2} and *z*_{3}(*t*) = *z*_{3}(0)*e*^{λ(1)t}_{3}. Since *λ*^{(1)}_{2} > 0, indeed *z*_{2} grows exponentially and *λ*^{(1)}_{3} < 0 implies that *z*_{3} decreases exponentially as well. Refer to Fig. 5 for a geometric exposition of these ideas. Thus, any perturbation with a nonzero component along **v**^{(1)}_{2} will grow. Hence, in this case ** η*** = (

*A*, 0, 0) is an unstable equilibrium point. Referring to Fig. 5, the angle

*θ*that

**v**

^{(1)}

_{2}makes with the positive

*F*axis (same as the angle

**v**

^{(1)}

_{3}makes with the negative

*F*axis) is given by

Thus, as *r* increases from 0 to 1, *θ* decreases from 90° to 0°.

#### 2) Case 2: Let *A* < 0 and *r* < 1

In this case since *λ*^{(1)}_{2} < 0 and *λ*^{(1)}_{3} > 0, the role of the eigenvalue–eigenvector pairs (*λ*^{(1)}_{2}, **v**^{(1)}_{2}) and (*λ*^{(1)}_{3}, **v**^{(1)}_{3}) are interchanged and the above analysis and conclusions hold. Thus, the equilibrium point ** η*** = (

*A*, 0, 0) is an unstable saddle point. Refer to Fig. 5.

##### (i) Effect of changing A

Thus, keeping r(<1) fixed, the force of attraction along **v**^{(1)}_{3} and the force of repulsion along **v**^{(1)}_{2} as measured by the magnitudes of the eigenvalues |*λ*^{(1)}_{3}| = *A**bc* and |*λ*^{(1)}_{2}| = *A**bc*, respectively, diminish to zero as *A* decreases from 1 to 0. As *A* increases in value in the negative direction, while the magnitudes of these forces increase, their role is reversed, with **v**^{(1)}_{2} being the new direction of attraction and **v**^{(1)}_{3} being the new direction of repulsion. This reversal of the role between **v**^{(1)}_{2} and **v**^{(1)}_{3} occurs when *A* = 0.

##### (ii) Effect of changing r

Let *A* > 0 be fixed. From (4.15) it follows that as *r* → 0^{+}, *c*/*b* → ∞, and *θ* → 90°. Combining this with (4.4) we get a degenerate limit, namely

As *r* increases form 0^{+}, the angle *θ* decreases from 90°. When *r* → 1^{−}, (*c*/*b*) → 0^{+} and *θ* → 0^{+}. Again we get a degenerate limit

When *r* = 1, *c* = 0 and the character of the equilibrium ** η*** = (A, 0, 0)

^{T}changes from being a saddle point to a center. Refer to the next subsection for more details.

#### 3) Case 3: Let *A* > 0 and *r* > 1 or *A* < 0 and *r* > 1

From (4.8) it follows that the nonzero eigenvalues are complex conjugate with zero real part. In this case the equilibrium point is called a center for the corresponding linearized problem (Hirsch and Smale 1974). A geometric view of the two new eigenvectors in (4.9) is shown in Fig. 6. It turns out that while the new eigenvector **v**^{(1)}_{3} = (0, 1, 1)^{T} is independent of *r* and makes a fixed 45° angle with the positive *F* axis, the new eigenvector **v**^{(1)}_{2} = (0, −|(*b*/*c*)|, |(*c*/*b*)|)^{T} makes an angle *θ* with the negative *F* axis where

As *r* increases from one to infinity, |*c*/*b*| increases form 0 to ½ and *θ* from 0 to 26.57°.

##### (i) Change of character of the equilibrium *η** *= (A, 0, 0)*^{T}

*η*

As *r* changes from 1^{−} to 1^{+}, it follows from (4.17) that the old eigenvectors **v**^{(1)}_{2} and **v**^{(1)}_{3} in (4.4) become degenerate and cease to exist. Simultaneously, a new set of eigenvectors **v**^{(1)}_{2} and **v**^{(1)}_{3} in (4.9) are created. Of these, **v**^{(1)}_{3} = (0, 1, 1)^{T} is fixed. When *r* → 1^{+}, |*c*/*b*| → 0 and hence new **v**^{(1)}_{2} → (0, −1, 0)^{T}. Comparing this with (4.17), it looks as though the old **v**^{(1)}_{3} is reincarnated as the new **v**^{(1)}_{2} as *r* changes from 1^{− }*to* 1^{+}. Thus, as *r* changes from 1^{−} to 1^{+}, ** η*** = (

*A*, 0, 0)

^{T}changes its character from being a saddle to a center.

#### 4) Case 4: Either *A* = 0 or *r* = 1

In this case, *λ*^{(1)}_{2} = 0 = *λ*^{(1)}_{3}.

### b. Stability of equilibrium points in *ɛ*_{F}

Let ** η*** = (0,

*F*, 0)

^{T}∈

*ɛ*. Then, the Jacobian in (4.1) becomes

_{F}whose eigenvalues and eigenvectors are given by

with *λ*^{(2)}_{1} = 0, *λ*^{(2)}_{2} = *F**ac*, *λ*^{(2)}_{3} = −*F**ac* and

The product *ac* and the ratio *c*/*a* are given by

and

Variation of *ac* and *c*/*a* as a function of *r* are given in Figs. 7a,b, respectively. From (4.20) and (4.22) it follows that the eigenvalues *λ*^{(2)}_{2} and *λ*^{(2)}_{3} depend on two factors: the value of the second coordinate *F* of the equilibrium point ** η*** = (0,

*F*, 0)

^{T}and the value of the ratio

*r*. A summary of the properties of these eigenvalues for various combinations of values of

*F*and

*r*are given in Table 3.

Again it can be verified that ** η*** = (0,

*F*, 0)

^{T}is an unstable saddle point when

*r*> 1 and

*F*≠ 0. However, when

*r*< 1 and F ≠ 0, perturbations around this equilibrium will have a local periodic behavior. It is interesting to note that the behavior of the equilibria in

*ɛ*for

_{A}*r*< 1 are quite similar to those in

*ɛ*for

_{F}*r*> 1 and vice versa.

### c. Stability of equilibrium points in *ɛ*_{G}

Let ** η*** = (0, 0,

*G*)

^{T}∈

*ɛ*Then the Jacobian in (4.1) becomes

_{G}.whose complex eigenvalues are given by *λ*^{(3)}_{1} = 0, *λ*^{(3)}_{2} = *iG*|*ab*|, *λ*^{(3)}_{3} = −*iG*|*ab*|. Using the standard canonical form, the eigen decomposition becomes

where

and

The product *ab* and ratio *b*/*a* are given by

From (4.25) it follows that the eigenvalues *λ*^{(3)}_{2} and *λ*^{(3)}_{3} are always complex conjugates (with zero real parts) for all values of *r*. Variation of *ab* and *b*/*a* as a function of *r* are given in Figs. 8a,b, respectively. Refer to Table 4 for a summary of the properties of these eigenvalues.

## 5. Small perturbation of equilibrium

In this section our goal is to quantify the growth of small errors superimposed on equilibrium points.

Let ** η*** be an equilibrium point and

**∈ ℝ**

*ξ*^{3}denote a (small) perturbation or error vector. Let

**=**

*η**** +**

*η***be the perturbed initial condition for the model (2.6). Since * =**

*ξ***g**(

***) = 0 by definition, we get**

*η*as the exact dynamics of the perturbation. Following the developments in section 4, we consider three cases, one for each of the three sets of equilibria.

### a. Small perturbation of equilibrium points in *ɛ*_{A}

Let ** η*** = (

*A*, 0, 0)

^{T}∈

*ɛ*. Depending on the values of

_{A}*A*and

*r*, consistent with section 4, we again consider four subcases.

#### 1) Case 1: Heteroclinic orbits

Let *A* > 0 and *r* < 1. In this case, ** η*** = (

*A*, 0, 0)

^{T}is a saddle point with

**v**

^{(1)}

_{2}= (0, 1,

*c*/

*b*)

^{T}denoting the local unstable direction and

**v**

^{(1)}

_{3}= (0, 1, −

*c*/

*b*)

^{T}denoting the local stable direction. Refer to Tables 1 and 2 for details.

These eigenvectors **v**^{(1)}_{2} and **v**^{(1)}_{3} partition the *F*–*G* plane (where *A* is constant) into four regions *R*_{1} to *R*_{4}. Refer to Fig. 9 for details. It turns out that there is an intrinsic difference in the qualitative behavior of the trajectories starting from each of these four regions but close to the equilibrium point ** η***.

##### (i) Heteroclinic orbits

The orbit or the trajectory starting from ** η** ∈ ℛ

_{1}that is close to the equilibrium point

*** = (1, 0, 0)**

*η*^{T}moves away from it along the unstable direction

**v**

^{(1)}

_{2}. In due course, this orbit bends with

*A*decreasing but

*F*remaining positive. It soon approaches the equilibrium point

**** = (−1, 0, 0)**

*η*^{T}, which is the symmetric dual of

***, along**

*η***v**

^{(1)}

_{2}, which is the local direction of attraction for

****. As this trajectory gets very close to**

*η*****, it gets repelled along**

*η***v**

^{(1)}

_{3}, which is the local unstable direction for

****. Again it bends with**

*η**A*increasing but

*F*still remaining positive. It now approaches

*** along**

*η***v**

^{(1)}

_{3}, which is its direction of attraction. As it gets closer to

***, it again gets repelled along**

*η***v**

^{(1)}

_{2}and this cycle repeats. Figure 10 shows the results of integration over a time period

*T*= 200 units using the adaptive fourth-order Runge–Kutta scheme.

The class of orbits so generated is called a heteroclinic orbit (Sparrow 1982). A distinguishing feature of this orbit is that one part of it connecting ** η*** and

**** lies at the intersection of unstable manifold of**

*η**** and stable manifold of**

*η*****. The other part of the orbit, from**

*η***** to**

*η****, lies at the intersection of the unstable manifold of**

*η***** and the stable manifold of**

*η****.**

*η*##### (ii) Four distinct heteroclinic orbits

It turns out that each equilibrium point ** η*** ∈

*ɛ*is endowed with four distinct classes heteroclinic orbits—one for each of the four regions ℛ

_{A}_{1}to ℛ

_{4}. Figures 11 –13 contain illustrations of heteroclinic orbits starting from initial points in regions ℛ

_{2}to ℛ

_{4}, respectively. The key distinguishing features of these four heteroclinic orbits are summarized in Table 5.

##### (iii) Effect of changing A

From Tables 1 and 2, it is clear that eigenvectors **v**^{(1)}_{2} and **v**^{(1)}_{3} of 𝗗** _{g}**(

***) for**

*η**** ∈**

*η**ɛ*are independent of

_{A}*A*. But their associated eigenvalues

*λ*

^{(1)}

_{2}=

*A*

*bc*and

*λ*

^{(1)}

_{3}= −

*A*

*bc*, which are measures of the force of repulsion and attraction are directly proportional to

*A*. Thus, as

*A*increases in magnitude, the resulting increases in the force of repulsion and attraction, while do not alter the basic properties of the heteroclinic orbits, do affect the frequency of oscillations of the components of

**. Examples of heteroclinic orbits for**

*η**A*= 5 and 10 with

**= (0.0, 0.001, 0.0)**

*ξ*^{T}and

*r*= 0.5 are given in Figs. 14 and 15.

##### (iv) Effect of changing r

From Fig. 4, it follows that *c*/*b* is a decreasing function of *r* with *c*/*b* being positive for 0 < *r <* 1 and negative for *r* > 1. We restrict this discussion to the case when 0 < *r* < 1 (*r* > 1 is considered in case 3 below). Referring to Fig. 9, it follows that (*c*/*b*) determines the slope of the eigenvectors **v**^{(1)}_{2} and **v**^{(1)}_{3} in the *F*–*G* plane. To understand the impact *r*, consider the example of the heteroclinic orbit in Fig. 10. In this example, since *F* remains positive, combining it Fig. 9, it can be verified that the heteroclinic orbit approaches ** η*** = (1, 0, 0)

^{T}along

**v**

^{(1)}

_{3}(that separates region ℛ

_{1}and ℛ

_{4}) and leaves of

*** along**

*η***v**

^{(1)}

_{2}(that separates region ℛ

_{1}and ℛ

_{2}). That is, this orbit makes a turn (to the right) by an angle

*ϕ*= 180 − 2

*θ*where

*θ*is given by (4.15). Clearly, this angle

*ϕ*is small or large depending on if

*r*is small or large.

##### (v) Effect of changing L and K

#### 2) Case 2: *A* < 0 and *r* < 1

In this case, the equilibrium point ** η*** = (−

*A*, 0, 0)

^{T}∈

*ɛ*is again a saddle point, with

_{A}**v**

^{(1)}

_{2}being the stable direction and

**v**

^{(1)}

_{3}the unstable direction. Hence all the comments and conclusions of case 1 carry over to this case.

#### 3) Case 3: *A* > 0 and *r* > 1 or *A* < 0 and *r* > 1

In this case *b* > 0 but *c* < 0 and recall from Table 1 that the eigenvalues *λ*^{(1)}_{2} = *iA*|*bc*| and *λ*^{(1)}_{3} = −*iA*|*bc*|. This in turn implies that ** η*** = (

*A*, 0, 0)

^{T}is a center and the components

*F*(

*t*) and

*G*(

*t*) oscillate where the frequency of oscillation depends on both

*A*and

*r*. In verifying this, assume that

*A*is nearly constant. In this case (2.6) becomes

Referring to Hirsch and Smale (1974), we get the solution of (5.2) as

Thus, *T* is a decreasing function of *A* and *r*.

The question now is: when will *A* be nearly a constant? Recall that when *r* is large, *a* < 0 but close to zero (*a* → 0^{−}). Hence *Ȧ* is small in magnitude but its sign depends on those of *F* and *G*. Refer to Fig. 16. Thus for large *r*, *A* is very slowly varying and the orbits of (2.6) look quite similar to a merry-go-round where the cars go round and round in closed circuits but they reach a maximum when they cross the *G* axis and a minimum when they cross the *F* axis. Refer to Fig. 17 for an example of the trajectory in this case.

#### 4) Case 4: *A* = 0 and *r* = 1

### b. Small perturbation of equilibrium points in *ɛ*_{F}

Comparing entries in Table 1, it turns out that there is a one to one correspondence between the stability properties of equilibria in *ɛ _{A}* and

*ɛ*. In particular, the properties of the eigenvalues and eigenvectors of the Jacobian of equilibria in

_{F}*ɛ*for

_{F}*r*> 1 (

*r*< 1) are very similar to those of

*ɛ*for

_{A}*r*< 1 (

*r*> 1). This in turn implies that the qualitative properties of the trajectories starting close to

*ɛ*when

_{F}*r*> 1 (

*r*< 1) are quite similar to those starting close to

*ɛ*when

_{A}*r*< 1 (

*r*> 1). In view of this similarity, to save space, we observe that all the remarks and conclusion of section 5a carry over to equilibria in

*ɛ*.

_{F}### c. Small perturbations of equilibrium points in *ɛ*_{G}

From (2.7) and (2.8) we have *a* < 0 and *b* > 0 for all *r* > 0. In view of this, from Table 1 it turns out that the eigenvalues of the Jacobian at the equilibrium point in *ɛ _{G}* are purely imaginary for all

*r*> 0. This in turn implies that the qualitative properties of trajectories starting close to

*ɛ*are quite similar to those starting close to

_{G}*ɛ*when

_{A}*r*> 1 or

*ɛ*for

_{F}*r*< 1. Hence all the relevant conclusions from section 5a carry over to this case as well.

## 6. Properties of orbits starting from an arbitrary point

### a. Case 1

This ratio depends on the initial point ** η** = (

*A*,

*F*,

*G*)

^{T}and the value of

*r*. Following Lorenz (1960) we define three types of initial conditions:

*R*

_{1}(

*l*) = 1 is called critical and

*R*

_{1}(

*l*) < (>)1 is called sub (sup) critical cases. From (6.1) it can be verified that

Let

Since *G*^{2}/*F*^{2} is nonnegative, we consider two subcases.

#### 1) 0 < *r <* 1

In this case *β* is positive. Hence, 0 < *R*_{1}(*l*) < 1 for all −*β* < *δ* < 0 and *R*_{1}(*l*) > 1 when *δ* > 0 with *R*_{1}(*l*) = 1 when *δ* = 0. As an example, if *r* = ½, then *β* = 15. Thus, *R*_{1}(*l*) is subcritical for all −15 < *δ* < 0, and supercritical for all *δ* > 0.

The cross sections of the *V* and *E* ellipsoids by the condition (6.3) are given by

and

Combining these with (6.1) we get

Thus, under the condition (6.3), *R*_{1}(*l*) = 1 when *δ* = 0. To check if an equilibrium point in *ɛ _{A}* is reachable, let

*G*= 0 in (6.3). Then

*F*= 0 and from (6.4) and (6.5) we get

which from (6.6) is consistent only when *δ* = 0. This argument leads to the following inescapable conclusion: orbits starting from initial points lying outside of the pair of straight lines given by

cannot approach points in *ɛ _{A}*.

A number of observations are in order.

- Critical points and eigenvectors of 𝗗
(_{g}*),*η** ∈*η**ɛ*. Combining (4.25) with (6.8), the latter can be written as Thus, the pair of straight lines given by_{A}*G*= ±(*c*/*b*)*F*in the*F*–*G*plane (refer to Fig. 9) indeed corresponds to the eigenvectors**v**^{(1)}_{2}= [0, 1, (*c*/*b*)]^{T}and**v**^{(1)}_{3}= [0, 1, −(*c*/*b*)]^{T}. Thus, points along these two eigenvectors in the*F*–*G*plane in Fig. 9 correspond to the critical case. Similarly, the region ℛ_{1}and ℛ_{3}correspond to the subcritical and region ℛ_{2}and ℛ_{4}to the supercritical cases. Critical points and heteroclinic orbits. Since the orbits starting from points along

**v**^{(1)}_{2}and**v**^{(1)}_{3}approach points in*ɛ*, these orbits are necessarily heteroclinic orbits corresponding to the unstable case in section 5a. However, orbits starting from the sub and super critical regions can never approach an equilibrium point in_{A}*ɛ*._{A}

#### 2) *r* > 1

Since

it follows from (6.1) that

Also, since *β* < 0, from (6.3) it follows that *δ* ≥ |*β*| > 0 for 2*G*^{2}/*F*^{2} to be nonnegative. Hence this corresponds to the supercritical case and in this case orbits starting from any ** η** = (

*A*,

*F*,

*G*)

^{T}∈ ℝ

^{3}cannot reach an equilibrium point in

*ɛ*. This is to be anticipated since when

_{A}*r*> 1, any equilibrium point in

*ɛ*is a center.

_{A}### b. Case 2

Analogous to (6.2) we readily see that

Let

Since the lhs of (6.13) is nonnegative, we consider two subcases.

#### 1) *r* > 1

In this case, since *α* is positive, we get 0 < *R*_{2}(*k*) < 1 (subcritical) when −*α* < *δ* < 0 and *R*_{2}(*k*) > 1 (supercritical) when *δ* > 0. In the critical case (*δ* = 0), using (4.23), we obtain

Thus, the pair of straight lines corresponding to (6.14) are given by *G* = ±(*c*/*a*)*A*. These straight lines correspond to the eigenvectors **v**^{(2)}_{2} = [1, 0, (*c*/*a*)]^{T} and **v**^{(3)}_{2} = [1, 0, −(*c*/*a*)]^{T}. Thus, the orbits starting from points along **v**^{(2)}_{2} and **v**^{(2)}_{3} approach the equilibrium points ** η*** = (0,

*F**, 0)

^{T}where

*F** = ±2

*υ*= ±2

*k*

*E*. As in case 1, orbits starting from sub (super) critical regions cannot approach an equilibrium points in

*ɛ*.

_{F}#### 2) *r* < 1

Since

it follows that *R*_{2}(*k*) > 1 for all 0 < *r* < 1. Hence this corresponds to the supercritical case and orbits starting from any initial point ** η** ∈ ℝ

^{3}cannot approach an equilibrium point in

*ɛ*. Again this is to be anticipated since

_{F}*** ∈**

*η**ɛ*is a center when

_{F}*r*< 1.

### c. Case 3

Since

for all *r*, *R*_{3}(*l*, *k*) < 1. Again it can be verified by following the above arguments that equilibrium ** η*** ∈

*ɛ*can never be reached by orbits starting from any point

_{G}**∈ ℝ**

*η*^{3}.

## 7. Index cycle

Lorenz (1960) demonstrated that the set of three ODEs in (2.1) were capable of portraying the tendency of the average strength of the midlatitude westerly winds to undergo considerable variations in intensity, a phenomenon known as the index cycle (Namias 1950). The index cycle represents a cyclical transformation of the general circulation from a pattern of strong, zonal flow in the midlatitudes along with weak subtropical flow (high index) to strong westerly flow in the subtropical latitudes combined with weak midlatitude flow (low index) and back again, a variation lasting over an extended period of several weeks. The first term in the *u* component of the wind (2.4) is independent of *x*, representing the zonal flow. The latitudes of the zonal wind maxima are fixed (*y* = 3*π*/2*l*, 7*π*/2*l*, etc.) and the intensity (*A*/*l*) varies, therefore Lorenz (1960) noted that the variable *A* could represent the zonal index. Large positive values of *A* correspond to strong west-to-east flow at locations of *u*-wind maxima, while large negative values of *A* result in a shift in the location of the maximum zonal flow to the south by a distance *π*/*l*. An index cycle in this model will therefore be indicated by the value of *A* changing from *A*_{max} to −*A*_{max} and back to *A*_{max}. Physically plausible values of *A* are determined by assignment of an arbitrary standard unit of vorticity, for example, Lorenz (1960) set 1 unit of vorticity = 1/3 h^{−1}. Given this, a value of *A* = 0.5 would correspond to a maximum zonal wind of approximately 75 m s^{−1}.

### a. Index cycles associated with heteroclinic orbits

In this work, we define the start of an index cycle to be the time when *A* first drops below 99% of *A*_{max} and the end of the cycle to be when *A* rises above 99% of *A*_{max}. The period of the model index cycle is defined as the length of time between the start and end of the cycle. Given this definition, we can examine the sensitivity of the index cycle period to changing initial conditions and values of *r*. Since we use the variable *A* to define the index cycle, we focus our analysis around the *A* equilibrium. The results of this sensitivity analysis are summarized in Tables 6 –8

An example of an index cycle portrayed by this model is provided in Fig. 10, given by a small *F* perturbation imposed upon *ɛ _{A}*. As noted in Fig. 9 this is a perturbation in the ℛ

_{1}region of the

*F*–

*G*domain, and results in a heteroclinic orbit with properties summarized in Table 5. For several units of time,

*A*remains near its initial value (

*A*= +1) then at approximately

*t*= 15, the

*F*and

*G*variables quickly increase while the value of

*A*rapidly switches from +1 to −1. As

*A*approaches −1, the values of

*F*and

*G*decrease to zero. The solution remains nearly the same for several units of time (around

*t*= 50). The solution suddenly shifts again at approximately

*t*= 85, when

*F*increases and

*G*decreases rapidly, and the value of

*A*switches from −1 to +1. The period of the index cycle in this case is 69.8 units of time.

A similar example is provided in Fig. 11, given by a small *G* perturbation imposed upon *ɛ _{A}*. This is a perturbation in the ℛ

_{2}region of the

*F*–

*G*domain, which results in a heteroclinic orbit with properties described in Table 5. As with the ℛ

_{1}example, the evolution of

*A*is consistent with the index cycle concept: A remains near its initial value (

*A*= +1) for several units of time, switches quickly from +1 to −1, remains near the minimum value for some time, then quickly changes again from −1 to +1. However, the behavior of the

*F*evolution in this example is similar in character to the

*G*evolution in the previous (ℛ

_{1}) example,

*F*changes sign as

*A*switches from +1 to −1 and back again, while

*G*remains positive throughout the evolution. In this case, the period of the index cycle is slightly longer (77.1 units of time) than in the previous example.

#### 1) Effect of changing *A*

In the unstable situation (*r* < 1) near *ɛ _{A}*, the index cycle period decreases as the initial value of

*A*increases (Table 6). It appears that relatively strong zonal flow (high initial

*A*) will result in a shorter index cycle. This makes physical sense, in the cases with stronger initial zonal flow, more energy will be available for transfer to the perturbations, allowing for a faster transformation of the flow. This is also consistent with the stability analysis of

*ɛ*described in section 4a(2). Larger values of

_{A}*A*produce a larger force of attraction/repulsion about the saddle point, which is consistent with a short index cycle period.

#### 2) Effect of changing *r*

The sensitivity of the index cycle period to variations in *r* is summarized in Table 7. Similarly to the previous example, for a given value of *A*, one expects the minimum index cycle period to occur when the magnitude of the eigenvalue is maximized. Since *λ*^{(1)}_{2} = *A**bc*, this would occur when the product *bc* is maximized. As shown in (4.10), the product *bc* is maximized for *r* = 0.6436. This is confirmed experimentally in Table 7 for values of *r* < 0.6 the index cycle period increases with decreasing *r*, and for values of *r* > 0.6 the index cycle increases with increasing *r*.

#### 3) Effect of changing perturbation size

The size of the initial perturbation also affects the length of the index cycle (Table 8). Larger perturbations result in shorter index cycles. Since the evolution of *F* in the ℛ_{1} class of solutions is characterized as an oscillation between the small initial value and some much larger maximum, the initial perturbation fixes the minimum value of *F*. In order for *A* to remain near its maximum/minimum value (as in a long period index cycle), *F* and *G* must remain near zero for a relatively long time. This will not happen if the initial *F* perturbation is not near zero. The same is true for *G* perturbations (ℛ_{2} class of solutions). It is interesting to note that, in the unstable case, identical index periods are obtained for *G* perturbation =(*c*/*b*) × *F* perturbation (for *r* = 0.5, (*c*/*b*) = 2.7386).

### b. Index cycles in the stable case

Lorenz (1960) stated that cycling of *A* was possible in the stable (*r* > 1) situation if the perturbations were large and *V*/*E* was supercritical. As shown in section 6, the critical point for *r* > 1 is determined by *G*^{2}/*A*^{2} = (*r*^{4} − 1)/2. However, in both the subcritical and supercritical cases, the evolution of *A* is not consistent with the definition of an index cycle. An example of a subcritical case is shown in Fig. 3. While *A* varies considerably in this example, it remains positive and does not remain near its maximum or minimum value for any considerable period of time. On the other hand, the *F* evolution is consistent with the index cycle definition. In the supercritical situation, the behavior of *A* and *G* are switched, *A* varies from *A*_{max} to −*A*_{max} and *G* remains positive for initial *G* > 0 (negative for initial *G* < 0). While this is consistent with the index cycle definition, *A* remains at its maximum or minimum value for only an instant, therefore the *A* evolution is inconsistent with the index cycle definition in the supercritical case (for *r* > 1) as well.

## 8. Summary

In this paper we provide a complete portrayal of the long-term behavior of the solution of the Lorenz’s (1960) maximum simplification equations. While the solution is always periodic, it can be broadly classified into two groups—the unstable case leading to heteroclinic orbits and the stable case of periodic orbits. It turns out that this dichotomy critically depends on the value of the ratio *V*/*E*. The unstable case arise when *V*/*E* = 2*l*^{2} or 2*k*^{2} and the stable case for all other values of this ratio. We also examine the effect of various parameters on the index cycle for both the unstable and stable cases.

In the light of this analysis it should be interesting to refine the data assimilation experiment in Lakshmivarahan et al. (2003) and the properties of the moment dynamics in Epstein (1969) to include the unstable heteroclinic orbits and stable periodic orbits.

## Acknowledgments

We wish to thank John Lewis for stimulating our interest in the lower-order models. When this was in progress, Dan Dobrovolschi of the National Meteorological Administration in Bucharest, Romania, suggested that interpretation of heteroclinic orbits is an interesting problem, which provided further impetus for this work. We wish to record our sincere thanks to Alex Gluhovsky, Yunheng Wang, and the two anonymous reviewers whose comments greatly improved the presentation.

## REFERENCES

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

## Footnotes

*Corresponding author address:* S. Lakshmivarahan, School of Computer Science, University of Oklahoma, Norman, OK 73019. Email: E-mail:varahan@ou.edu