Bitstream 1780177 PDF

Download as pdf or txt
Download as pdf or txt
You are on page 1of 15

11 Lorenz equations

In this lecture we derive the Lorenz equations, and study their behavior.

The equations were first derived by writing a severe, low-order truncation of


the equations of R-B convection.

One motivation was to demonstrate SIC for weather systems, and thus point
out the impossibility of accurate long-range predictions.

Our derivation emphasizes a simple physical setting to which the Lorenz


equations apply, rather than the mathematics of the low-order truncation.

See Strogatz, Ch. 9, for a slightly different view. This lecture derives from Tritton, Physical Fluid
Dynamics, 2nd ed. The derivation is originally due to Malkus and Howard.

11.1 Physical problem and parameterization

We consider convection in a vertical loop or torus, i.e., an empty circular


tube:
cold

hot

We expect the following possible flows:

• Stable pure conduction (no fluid motion)


• Steady circulation
• Instabilities (unsteady circulation)

109
The precise setup of the loop:

T0−T1 (external)
T0−T2

z
(T0−T3) (T0+T 3) g
q
a φ

T0+T 2
T0+T 1 (external)

θ = position round the loop.

External temperature TE varies linearly with height:

TE = T0 − T1 z/a = T0 + T1 cos θ (24)

Let a be the radius of the loop. Assume that the tube’s inner radius is much
smaller than a.

Quantities inside the tube are averaged cross-sectionally:

velocity = q = q(θ, t)
temperature = T = T (θ, t) (inside the loop)

As in the Rayleigh-Bénard problem, we employ the Boussinesq approximation


(here, roughly like incompressiblity) and therefore assume
ωδ
= 0.
ωt
Thus mass conservation, which would give � · ψu in the full problem, here
gives
ωq
= 0. (25)
ωθ

110

Thus motions inside the loop are equivalent to a kind of solid-body rotation,
such that
q = q(t).

The temperature T (θ) could in reality vary with much complexity. Here we
assume it depends on only two parameters, T2 and T3, such that

T − T0 = T2 cos θ + T3 sin θ. (26)

Thus the temperature difference is

• 2T2 between the top and bottom, and


• 2T3 between sides at mid-height.

T2 and T3 vary with time:

T2 = T2(t), T3 = T3(t)

11.2 Equations of motion

11.2.1 Momentum equation

Recall the Navier-Stokes equation for convection:


ωψu ψ u = − 1 �p
ψ − ψg�ΓT + ��2ψu
+ ψu · �ψ
ωt δ
We write the equivalent equation for the loop as
ωq 1 ωp
=− + g�(T − T0) sin θ − �q. (27)
ωt δa ωθ
The terms have the following interpretation:

• ψu ∗ q
• ψu · �ψu ∗ 0 since ωq/ωθ = 0.

111
1 γp
• �p ∗ a γπ by transformation to polar coordinates.
• A factor of sin θ modifies the buoyancy force F = g�(T − T0 ) to obtain
the tangential component:

F
Fsinφ
φ φ

The sign is chosen so that hot fluid rises.


• � is a generalized friction coefficient, corresponding to viscous resistance
proportional to velocity.

Now substitute the expression for T − T0 (equation (26)) into the momentum
equation (27):
ωq 1 ωp
=− + g�(T2 cos θ + T3 sin θ) sin θ − �q
ωt δa ωθ
Integrate once round the loop, with respect to θ, to eliminate the pressure
term: � 2�
ωq
2α = g� (T2 cos θ sin θ + T3 sin2 θ)dθ − 2α�q.
ωt 0
The pressure term vanished because
� 2�
ωp
dθ = 0,
0 ωθ
i.e., there is no net pressure gradient around the loop.

The integrals are easily evaluated:


� 2� �2�
1 2 ��
cos θ sin θ dθ = sin θ�
= 0
0 2 0

and � 2�
sin2 θ dθ = α.
0

112

Then, after dividing by 2α, the momentum equation is


dq g�T3
= −�q + (28)
dt 2
where we have written dq/dt instead of ωq/ωt since ωq/ωθ = 0.

We see that the motion is driven by the horizontal temperature difference,


2T3.

11.2.2 Temperature equation

We now seek an equation for changes in the temperature T . The full tem­
perature equation for convection is
ωT ψ = τ�2 T
+ ψu · �T
ωt
where τ is the heat diffusivity.

We approximate the temperature equation by considering only cross-sectional


averages within the loop:
ωT q ωT
+ = K(TE − T ) (29)
ωt a ωθ
Here we have made the following assumptions:

• RHS assumes that heat is transferred through the walls at rate


K(Texternal − Tinternal).
• Conduction round the loop is negligible (i.e., no �2T ).
q γT
• a γπis the product of averages, not (as it should be) the average of a
product; i.e., q is taken to be uncorrelated to ωT /ωθ.

Recall that we parameterized the internal temperature with two time-dependent


variables, T2(t) and T3 (t). We also have the external temperature TE varying
linearly with height. Specifically:
TE = T0 + T1 cos θ
T − T0 = T2 cos θ + T3 sin θ

113
Subtracting the second from the first,
TE − T = (T1 − T2) cos θ − T3 sin θ.
Substitute this into the temperature equation (29):
dT2 dT3 q q
cos θ + sin θ − T2 sin θ + T3 cos θ = K(T1 − T2) cos θ − KT3 sin θ.
dt dt a a
Here the partial derivatives of T have become total derivatives since T 2 and
T3 vary only with time.

Since the temperature equation must hold for all θ, we may separate sin θ
terms and cos θ terms to obtain
dT3 qT2
sin θ : − = −KT3
dt a
dT2 qT3
cos θ : + = K(T1 − T2)
dt a
These two equations, together with the momentum equation (28), are the
three o.d.e.’s that govern the dynamics.

We proceed to simplify by defining


T4(t) = T1 − T2(t),
which is the difference between internal and external temperatures at the top
and bottom—loosely speaking, the extent to which the system departs from
a “conductive equilibrium.” Substitution yields
dT3 qT1 qT4
= −KT3 + −
dt a a
dT4 qT3
= −KT4 +
dt a

11.3 Dimensionless equations

Define the nondimensional variables


q g�T3 g�T4
X= , Y = , Z=
aK 2a�K 2a�K
114
Here
X = dimensionless velocity
Y = dimensionless temperature difference between up and down currents
Z = dimensionless departure from conductive equilibrium

Finally, define the dimensionless time


t∗ = tK.
Drop the prime on t to obtain
dX
= −P X + P Y
dt
dY
= −Y + rX − XZ
dt
dZ
= −Z + XY
dt
where the dimensionless parameters r and P are
g�T1
r = = “Rayleigh number”
2a�K

P = = “Prandtl number”
K

These three equations are essentially the same as Lorenz’s celebrated system,
but with one difference. Lorenz’s system contained a factor b in the last
equation:
dZ
= −bZ + XY
dt
The parameter b is related to the horizontal wavenumber of the convective
motions.

11.4 Stability

We proceed to find the fixed points and evaluate their stability. For now, we
remain with the loop equations (b = 1).

115
The fixed points, or steady solutions, occur where
Ẋ = Ẏ = Ż = 0.
An obvious fixed point is
X � = Y � = Z � = 0,
which corresponds, respectively, to a fluid at rest, pure conduction, and a
temperature distribution consistent with conductive equilibrium.

Another steady solution is



X� = Y � = ± r − 1
Z� = r − 1
This solution corresponds to flow around the loop at constant speed; the ±
signs arise because the circulation can be in either sense. That sgn(X) =
sgn(Y ) implies that hot fluid rises and cold fluid falls.

Note that the second (convective) solution exists only for r > 1. Thus we
see that, effectively, r = Ra/Rac , i.e., the convective instability occurs when
Ra > Rac .

As usual, we determine the stability of the steady-state solutions by deter­


mining the sign of the eigenvalues of the Jacobian.

Let ⎪

⎪ �
X
X�
ψ =

Y ⎞ ,

θ θ� = ⎛ Y � ⎞

Z Z�
Then the Jacobian matrix is
� ⎭

� −P +P 0
ωθ̇i �
� =

r − Z � −1 −X � ⎤

ωθj �

π Y� X � −1
The eigenvalues ε are found by equating the following determinant to zero:

� �
� −(ε + P ) P 0 �
� �
� r − Z � −(ε + 1) −X � � = 0

� �
� �

Y X −(ε + 1) �

116

For the steady state without circulation (X � = Y � = Z � = 0), we have


� �
� −(ε + P ) P 0 �
� �

� r −(ε + 1) 0 � = 0.



0
0 −(ε + 1)

This yields
−(ε + P )(ε + 1)2 + rP (ε + 1) = 0
or � �

(ε + 1)
ε 2 + ε(P + 1) − P (r − 1)
= 0.

There are three roots:


ε1 = −1

−(P + 1) (P + 1)2 + 4P (r − 1)
ε2,3 = ±
2 2

As usual,
Re{ε1, ε2, and ε3} < 0 =∞ stable
Re{ε1, ε2, or ε3} > 0 =∞ unstable
Therefore X � = Y � = Z � = 0 is
stable for 0<r<1
unstable for r>1


We now�calculate the stability of the second fixed point, X � = ± r − 1,
Y � = ± r − 1, Z � = r − 1.

The eigenvalues ε are now the solution of


� �
� −(ε + P ) P 0 �
� � �

� 1 −(ε + 1) −S � = 0,
� S
= ± r − 1.


S
S −(ε + 1) �

(Explicitly,
−(ε + p)(ε + 1)2 − P s2 − S 2 (ε + P ) + P (ε + 1) = 0
(ε + 1)[ε 2 + ε(P + 1)] + εS 2 + 2P S 2 = 0.)

117
The characteristic equation is cubic:
ε 3 + ε 2 (P + 2) + ε(P + r) + 2P (r − 1) = 0
This equation is of the form
ε 3 + Aε 2 + Bε + C = 0 (30)
where A, B, and C are all real and positive.

Such an equation has either

• 3 real roots; or
• 1 real root and 2 complex conjugate roots, e.g.,

or

3 real roots 1 real and 2 complex conjugates

Rearranging equation (30),


ε (ε 2 + B) = −Aε

2
⎜� − C� < 0.
� ⎜� �
positive real negative real

Consequently any real ε < 0, and we need only consider the complex roots
(since only they may yield Re{ε} > 0).

Let ε1 be the (negative) real root, and let


ε2,3 = � ± iλ.
Then
(ε − ε1 )(ε − � − iλ)(ε − � + iλ) = 0
and
A = −(ε1 + 2�)

B = 2�ε1 + �2 + λ 2

C = −ε1(�2 + λ 2)

118
A little trick: � �
C − AB = 2� (ε1 + �)2 + λ 2 .
� ⎜� �
positive real

Since � is the real part of both complex roots, we have

sgn(Re{ε2,3}) = sgn(�) = sgn(C − AB).

Thus instability occurs for C − AB > 0, or

2P (r − 1) − (P + 2)(P + r) > 0, .

Rearranging,
r(2P − P − 2) > 2P + P (P + 2)
and we find that instability occurs for
P (P + 4)
r > rc = .
P −2
This condition, which exists only for P > 2, gives the critical value of r for
which steady circulation becomes unstable.

Loosely speaking, this is analogous to a transition to turbulence.

Summary: The rest state, X � = Y � = Z � = 0, is

stable for 0<r<1


unstable for r > 1.

The convective state (steady circulation), X � = Y � = ± r − 1,
Z � = r − 1, is

stable for 1 < r < rc


unstable for r > rc .

What happens for r > rc ?

Before addressing that interesting question, we first look at contraction of


volumes in phase space.

119
11.5 Dissipation

We now study the “full” equations, with the parameter b, such that
Ż = −bZ + XY, b > 0.
The rate of volume contraction is given by the Lie derivative
1 dV � ωθ̇i
= , i = 1, 2, 3, θ1 = X, θ2 = Y, θ3 = Z.
V dt i
ωθ i

For the Lorenz equations,


ω Ẋ ω Ẏ ω Ż
+ + = −P − 1 − b.
ωX ωY ωZ
Thus
dV
= −(P + 1 + b)V
dt
which may be solved to yield
V (t) = V (0)e−(P +1+b)t.
The system is clearly dissipative, since P > 0 and b > 0.

The most common choice of parameters is that chosen by Lorenz


P = 10
b = 8/3 (corresponding to the first wavenumber to go unstable).
For these parameters,
41
V (t) = V (0)e− 3 t .
41
Thus after 1 time unit, volumes are reduced by a factor of e− 3 � 10−6. The
system is therefore highly dissipative.

11.6 Numerical solutions

For the full Lorenz system, instability of the convective state occurs for
P (P + 3 + b)
r > rc =
P −1−b
120
For P=10, b=8/3, one has
rc = 24.74.
In the following examples, r = 28.

Time series of the phase-space variables are shown in

Tritton, Fig 24.2, p. 397

• X(t) represents variation of velocity round the loop.


– Oscillations around each fixed point X+� and X−� represent variation
in speed but the same direction.
– Change in sign represents change in direction.
• Y (t) represents the temperature difference between up and downggoing
currents. Intuitively, we expect some correlation between X(t) and Y (t).
• Z(t) represents the departure from conductive equilibrium. Intuitively,
we may expect that pronounced maxima of Z (i.e., overheating) would
foreshadow a change in sign of X and Y , i.e., a destabilization of the
sense of rotation.

Projection in the Z-Y plane, showing oscillations about the unstable convec­
tive fixed points, and flips after maxima of Z:

BPV, Fig. VI.12

A 3-D perspective, the famous “butterfly:”

BPV, Fig. VI.14

Note the system is symmetric, being invariant under the transformation X ∗


−X, Y ∗ −Y , Z ∗ Z.

A slice (i.e., a Poincaré section) through the plane Z = r − 1, which contains

the convective fixed points:

BPV, Fig. VI.15

121
• The trajectories lie on roughly straight lines, indicating the attractor
dimension d ◦ 2.
• These are really closely packed sheets, with (as we shall see) a fractal
dimension of 2.06.
• d ◦ 2 results from the strong dissipation.

Since d ◦ 2, we can construct, as did Lorenz, the first return map

zk+1 = f (zk ),

where zk is the kth maximum of Z(t). The result is

BPV, Fig. VI.16

(These points intersect the plane XY −bZ = 0, which corresponds to Ż = 0.)

The first-return map shows that the dynamics can be approximated by a 1-D
map. It also reveals the stability properties of the fixed point Z = r − 1:

BPV, Fig. VI.17

Finally, sensitivity to initial conditions is documented by

BPV, Fig. VI.18

11.7 Conclusion

The Lorenz model shows us that the apparent unpredictability of turbulent


fluid dynamics is deterministic. Why?

Lorenz’s system is much simpler than the Navier-Stokes equations, but it is


essentially contained within them.

Because the simpler system exhibits deterministic chaos, surely the Navier-
Stokes equations contain sufficient complexity to do so also.

122
Thus any doubt concerning the deterministic foundation of turbulence, such
as assuming that turbulence represents a failure of deterministic equations,
is now removed.

A striking conclusion is that only a few (here, three) degrees of freedom are
required to exhibit this complexity. Previous explanations of transitions to
turbulence (e.g., Landau) had invoked a successive introduction of a large
number of degrees of freedom.

123

You might also like