7.6 General Hamiltonian flows
8.1.3 Continuous dependence on initial data
A priori estimate
1. We recall in logistic equation y0= ry(1 − y/K), if a solution y(·) with 0 < y(0) < K exists, then it always satisfies 0 < y(t) < K for all t. Such kind of estimate is called an a priori estimate.
2. In spring-mass model m¨x − kx = 0, if the solution exists, then it always satisfies 1
2x˙2+ kx2 = E
for some constant E > 0. This automatically gives boundedness of (x(t), ˙x(t)), as long as it exists. The estimate is called a priori estimate.
A vector field f (t, y) is said to grow at most linearly as |y| → ∞ is that
|f (t, y)| ≤ a|y| + b (8.6)
for some positive constants a, b.
Theorem 8.12. If f (t, y) is smooth and grows at most linearly as |y| → ∞, then all solutions of ODEy0= f (t, y) can be extended to t = ∞.
Proof. Suppose a solution exists in [0, T ), we give a priori estimate for this solution. From the grow condition of f , we have
|y(t)|0≤ |y0(t)| ≤ a|y(t)| + b.
Multiplying e−aton both sides, we get
e−at|y(t)|0
≤ e−atb.
8.1. WELL-POSHNESS 161 Integrating t from 0 to T , we obtain
e−aT|y(T )| − |y(0)| ≤ Z T
0
e−atb dt = b
a 1 − eaT . Hence
|y(T )| ≤ |y(0)|eaT + b aeaT.
Such an estimate is called a priori estimate of solutions. It means that as long as solution exist, it satisfies the above estimate.
Now suppose our solution exists in [0, T ) and cannot be extended. From the above estimate, the limit y(T −) exists. This is because y(·) is bounded, hence y0(t) = f (t, y(t)) is bounded and continuous for t ∈ [0, T ). Hence the limit
y(T −) = lim
t→T −y(0) + Z t
0
y0(s) ds
exists. We can extend y(·) from T with the y(T +) = y(T −). By the local existence theorem, the solution can be extended for a short time. Now, we have a solution on two sides of T with the same data y(T −), we still need to show that it satisfies the equation at t = T . To see this, on the right-hand side
t→T +lim y0(t) = lim
t→T +f (t, y(t)) = f (T, y(T −)).
On the left-hand side, we also have
t→T −lim y0(t) = lim
t→T −f (t, y(t)) = f (T, y(T −)).
Therefore y0(t) is continuous at T and y0(T ) = f (T, y(T )). Hence we get the extended solution also satisfies the equation at T . This is a contradiction.
Remarks.
1. We can replace the growth condition by
|f (t, y)| ≤ a(t)|y| + b(t) (8.7)
where a(t) and b(t) are two positive functions and locally integrable, which means Z
I
a(t) dt, Z
I
b(t) dt < ∞ for any bounded interval I.
2. In the proofs of the uniqueness theorem and the global existence theorem, we use so called the Gronwall inequality, which is important in the estimate of solutions of ODE.
Lemma 8.1 (Gronwall inequality). If
η0(t) ≤ a(t)η(t) + b(t) then
η(t) ≤ e
Rt
0a(s) dsη(0) + Z t
0
e
Rt
sa(τ ) dτb(s) ds (8.8)
Gronwall inequality can be used to show that the continuous dependence of solution to its initial data.
Lyapunov functional
Theorem 8.13. Consider the ODE in Rn:
y0= f (y), y(0) = y0. Suppose there exists a functionΦ such that
(i) ∇Φ(y) · f (y) ≤ 0, and (ii) Φ(y) → ∞ as |y| → ∞.
Then the solution exists on[0, ∞).
Proof. Consider Φ(y(t)). It is a non-increasing function because d
dtΦ(y(t)) = ∇Φ(y(t)) · f (y(t)) ≤ 0 Thus,
Φ(y(t)) ≤ Φ(y(0)) Since Φ(y) → ∞ as y → ∞, the set
{y|Φ(y) ≤ Φ(y0)}
is a bounded set. If the maximal existence of interval is [0, T ) with T < ∞, then y(·) is bounded in [0, T ) and can be extended to T . By the local existence of ODE, we can always extend y(·) to T + . This is a condiction. Hence T = ∞.
As an example, let us consider a damping system
¨
x + γ ˙x = −V0(x)
where V is a trap potential, which means that V (x) → ∞ as |x| → ∞. By multiplying ˙x both sides, we obtain
dE
dt = −γ| ˙x|2 ≤ 0
8.2. SUPPLEMENTARY 163 Here,
E(t) := 1
2| ˙x|2+ V (x)
is the energy. The term γ| ˙x|2is called the energy dissipation rate. We integrate the above equation from 0 to t, drop the dissipation term to get
E(t) ≤ E(0), for all t > 0.
This gives a priori estimate of solution 1
2| ˙x(t)|2+ V (x(t)) ≤ E(0).
This implies both ˙x(t) and x(t) are bounded, because of the property of V .
8.2 Supplementary
8.2.1 Uniform continuity
Pointwise continuity. The concept of continuity is a local concept. Namely, y is continuous at t0 means that for any > 0 there exists δ > 0 such that |y(t) − y(t0)| < as |t − t0| < δ. The continuity property of y at t0is measured by the relation δ(). The locality here means that δ also depends on t0. This can be read by the example y = 1/t for t0 ∼ 0. For any , in order to have
|1/t − 1/t0| < , we can choose δ ≈ t20(Check by yourself). Thus, the continuity property of y(t) for t0 near 0 and 1 is different. The ratio /δ is of the same magnitude of y0(t0), in the case when y(·) is differentiable.
Uniform continuity
Theorem 8.14. When a function y is continuous on a bounded closed interval I, the above local continuity becomesuniform. Namely, for any > 0, there exists a δ > 0 such that |y(t1) − y(t2)| <
whenever |t1− t2| < δ.
Proof. For any > 0, any s ∈ I, there exists δ(, s) > 0 such that |y(t) − y(s)| < whenever
|t − s| < δ(, s). Let us consider the open intervals U (s, δ(, s)) := (s − δ(, s), s + δ(, s)). The union ∪s∈IU (s, δ(, s)) contain I. Since I is closed and bounded, by so called the finite convering lemma, there exist finite many U (si, δ(, si)), i = 1, ..., n such that I ⊂ ∪ni=1U (si, δ(, si)). Then we choose
δ :=minn
i=1δ(, si)
then the distances between any pair siand sjmust be less than δ. For any t1, t2 ∈ I with |t1− t2| <
δ, Suppose t1 ∈ U (sk, δ(, sk)) and t2 ∈ U (sl, δ(, sl)), then we must have |sk− sl| < δ.
|y(t1) − y(t2)| ≤ |y(t1) − y(sk)| + |y(sk) − y(sl)| + |y(sl) − y(t2)| < 3.
This completes the proof.
The key of the proof is the finite convering lemma. It says that a local property can be uniform through out the whole interval I. This is a key step from local to global.
8.2.2 C(I) is a normed linear space
If this distance is zero, it implies y1≡ y2in I. Also, kayk = |a|kyk for any scalar a. Moreover, we have
ky1+ y2k ≤ ky1k + ky2k.
If we replace y2 by −y2k, it says that the distance between the two function is less than ky1k and ky2k. This is exactly the triangular inequality. To show this inequality, we notice that
|y1(t)| ≤ ky1k, |y2(t)| ≤ ky2k, for all t ∈ I Hence,
|y1(t) + y2(t)| ≤ |y1(t)| + |y2(t)| ≤ kbf y1k + ky2k.
By taking maximal value on the L.H. side for t ∈ I, we obtain ky1+ y2k ≤ ky1k + ky2k.
The function space C[I] with the norm k · k is called a normed vector space.
8.2.3 C(I) is a complete
Such a space is called a Banach space.
Definition 8.5. A sequence {yn} is called a Cauchy sequence if for any > 0, there exists an N such that for anym, n ≥ N , we have
kyn− ymk < .
Theorem 8.15. Let {yn} be a Cauchy sequence in C(I). Then there exist y ∈ C(I) such that kyn− yk → 0 as n → ∞.
To prove this theorem, we notice that for each t ∈ I, {yn(t)} is a Cauchy sequence in R. Hence, the limit limn→∞yn(t) exists. We define
y(t) = lim
n→∞yn(t) for each t ∈ I.
We need to show that y is continuous and kyn− yk → 0. To see y is continuous, let t1, t2 ∈ I.
At these two points, limnyn(ti) = y(ti), i = 1, 2. This means that for any > 0, there exists an N > 0 such that
|yn(ti) − y(ti)| < , i = 1, 2, for all n ≥ N.
8.2. SUPPLEMENTARY 165 With this, we can estimate |y(t1) − y(t2)| through the help of ynwith n ≥ N . Namely,
|y(t1) − y(t2)| ≤ |y(t1) − yn(t1)| + |yn(t1) − yn(t2)| + |yn(t2) − y(t2)|
≤ 2 + |yn(t1) − yn(t2)| ≤ 3
In the last step, we have used the uniform continuity of ynon I. Hence, y is continuous in I.
Also, from the Cauchy property of yn in C(I), we have for any > 0, there exists an N > 0 such that for all n, m > N , we have
kyn− ymk < But this implies that for all t ∈ I, we have
|yn(t) − ym(t)| < Now, we fix n and let m → ∞. This yields
|yn(t) − y(t)| ≤
and this holds for n > N . Now we take maximum in t ∈ I. This yields kyn− yk ≤
Thus, we have shown lim yn= y in C(I).
Chapter 9
Stability Analysis
9.1 Damping
In this section, we consider dissipative nonlinear oscillators. The dissipation is due to friction. The friction force is usually a function of velocity. Let us call it b( ˙y). In general, we consider
¨
y = F (y) + b( ˙y) (9.1)
The friction force has the property:
b( ˙y) · ˙y < 0 and b(0) = 0.
This means that the direction of the frictional force is in the opposite direction of the velocity and the friction force is zero if the particle is at rest. Here are two concrete examples of damping.
Simple pendulum with damping The equation for simple pendulum is
ml ¨θ = −mg sin θ.
A simple damping force is portortional to the angular speed β ˙θ, provided the damping comes from the friction at the fixed point. Here, β > 0. Thus the model for simple pendulum with friction reads
ml ¨θ = −β ˙θ − mg sin θ. (9.2)
An active shock absorber
In the mass-spring model, the friction force may depend on the velocity nonlinear, say β(v), say β(v) = v4. Then the corresponding oscillation is nonlinear:
m¨y = −β( ˙y)¨y − ky, β(v) = v4, (9.3) 167
9.1.1 Stability and asymptotic stability
Definition 9.1. An equilibrium ¯y of the ODE ˙y = f (y) is said to be stable if for any > 0, there exists aδ > 0 such that for any solution y(·) with |y(0) − ¯y| < δ, we have |y(t) − ¯y| < ε.
Definition 9.2. An equilibrium ¯y of the ODE ˙y = f (y) is said to be asymptotically stable if there exists aδ > 0 such that any solution y(·) with |y(0) − ¯y| < δ satisfies y(t) → ¯y as t → ∞.
Definition 9.3. An equilibrium ¯y of the ODE ˙y = f (y) is said to be exponentially stable if there exist anα > 0 and a δ > 0 such that any solution y(·) with |y(0) − ¯y| < δ satisfies y(t) − ¯y = O(e−αt) as t → ∞.
Remarks.
• It is clear that asymptotical stability implies stability.
• For linear systems, centers are stable whereas sinks and spiral sinks are asymptotically stable.
• For hamiltonian system, the minimum of a hamiltonian H is a stable center.
Now we shall perturb a hamiltonian system in a special way, called dissipative perturbation.
In this case, the center becomes an asymptotical stable equilibrium. We shall use the following example for demonstration.
Consider the Hamiltonian H(x, p) = p2/2 + V (x). Let us assume
• lim|x|→∞V (x) = ∞
• 0 is the unique minimum of V with V (0) = 0.
Let us perturb this mechanical system by some damping. This means that we the particle exerted a friction force b(p), where p = ˙x is the velocity. This term is a friction if
b(0) = 0, b(p)p < 0
The latter simply means that the force is in the opposite direction of the velocity. Thus, the system reads
x˙ = p
˙
p = −V0(x) + b(p) You can check that along any trajectory (x(t), v(t)),
d
dt(H(x(t), p(t)) = Hxx + H˙ pp = −V˙ 0(x)p + p(−V0(x) − b) = bp < 0.
Thus, the Hamiltonian h decreases along any trajectories until p = 0. Such a perturbation is called a dissipative perturbation. As a result, we can see that (0, 0) becomes asymptotic stable. Indeed, we shall show in the section of Liapunov function that (x(t), p(t)) → (0, 0) as t → ∞ for any trajectories. Here, we just show that, from the linear analysis, the center becomes a spiral sink for
9.1. DAMPING 169 a Hamiltonian system with dissipative perturbation. We shall assume b0(0) 6= 0. The variational matrix at (0, 0) is
0 1
−Hxx −Hxp+ b0(0)
=
0 1
−V00(0) −b0(0)
Its eigenvalues are
λ± = b0(0) ± ip V00(0)
Now, the force b(p) is a friction which means that b(p)p < 0. But b(0) = 0. We get that 0 > b(p)p ∼ b0(0)p · p
Thus, if b0(0) 6= 0, then b0(0) < 0. Hence, (0, 0) becomes a spiral sink.
9.1.2 Stability and Lyapunov method Let us go back to consider the general problem
¨
y = F (y) + b( ˙y), with b( ˙y) · ˙y < 0, b(0) = 0. (9.4) Suppose F (¯y) = 0. Then steady state y(t) ≡ ¯y is a solution. Here, we have used b(0) = 0. We call the state (¯y, 0) an equilibrium on the phase plane.
Recall in the case of linear spring with damping, we have seen in Chapter 2 that every solution tends to the zero state as t tends to infinity. This zero state is an equilibrium and it is globally stable(i.e. any solution tends to this equilibrium as t → ∞). Let us study the same global stability problem for the damped nonlinear equation(9.4). Let us suppose that the potential corresponding to the force F (y) is V (y), i.e. V0(y) = −F (y). We consider the following damped system:
¨
y = −V0(y) + b( ˙y). (9.5)
The following theorem give a sufficient condition for global stability of the equilibrium.
Theorem 9.1. Suppose V (y) → ∞ as |y| → ∞ and V (y) has only one minimum ¯y. Then any solutiony satisfies
y(t) → ¯y and ˙y(t) → 0 as t → ∞.
Proof. For simplicity, we assume b( ˙y) = − ˙y. It should be clear if this case is done, how it can be generalized to more general cases. Without loss of generality, we may also assume ¯y = 0 and V (0) = 0. Otherwise, we may just replace y by y − ¯y and V (y) by V (y) − V (¯y), which does not alter F (y) in the original problem.
We use energy method: multiplying ˙y on both sides of (9.4), we obtain
˙
y ¨y = −V0(y) ˙y − ˙y ˙y dE
dt = − ˙y2, (9.6)
where
E(y, ˙y) := y˙2
2 + V (y). (9.7)
The strategy is to prove (i) E(t) → 0, and (ii) E(y, ˙y) = 0 if and only if (y, ˙y) = (0, 0), and (iii) (y(t), ˙y(t)) → (0, 0). We divide the proof into the following steps.
Step 1. From (9.6), E(t) := E(y(t), ˙y(t)) is a decreasing function along any trajectort (y(t), ˙y(t)).
Further, it has lower bound, namely, E(y, ˙y) ≥ 0. we get E(t) & α as t → ∞ for some number α.
Step 2. Let us call the limiting set of (y(t), ˙y(t)) by Ω+. That is
Ω+= {(y, ˙y)|∃tn, tn→ ∞ s.t. (y(tn), ˙y(tn)) → (y, ˙y)}.
Such a set is called an ω-limit set. We claim that any trajectory (˜y(·), ˙˜y(·)) with initial data (˜y(0), ˙˜y(0)) ∈ Ω+lies on Ω+forever. The proof of this claim relies on the continuity theorem on the initial data. Namely, the solution of an ODE depends on its initial data continuously. Let us accept this fact. Suppose (˜y(0), ˙˜y(0)) ∈ Ω+, we want to prove that for any fixed s > 0, (˜y(s), ˙˜y(s)) ∈ Ω+. Given fixed s > 0, by the continuous dependence of initial data, we have for any > 0, there exists a δ > 0 such that if |(y1, ˙y1) − (˜y(0), ˙˜y(0))| < δ, then the solution y(·) with initial data (y1, ˙y1) is in an neighborhood of ˜y(s). Now, since (˜y(0), ˙˜y(0)) ∈ Ω+, with this δ > 0, there exist tn such that |(y(tn), ˙y(tn)) − (˜y(0), ˙˜y(0))| < δ. Let us consider two solutions, one has initial data (y(tn), ˙y(tn)), the other has initial data (˜y(0), ˙˜y(0)). By the continuity dependence of the initial data, we get (y(tn+ s), ˙y(tn+ s)) − (˜y(s), ˙˜y(s))| < . This yields that ∀ > 0, there exists an n such that |(y(tn+ s), ˙y(tn+ s)) − (˜y(s), ˙˜y(s))| < . Thus, (˜y(s), ˙˜y(s)) ∈ Ω+.
Step 3. We claim that E(˜y(s), ˙˜y(s)) = α for any s ≥ 0. This is because (1) for any fixed s, there exist tn → ∞ such that (y(tn+ s), ˙y(tn+ s)) → (˜y(s), ˙˜y(s)) as n → ∞, and (2), from step 1, E(y(t), ˙y(t)) & α as t → ∞. Thus, we get E(y(tn+ s), ˙y(tn+ s)) → α for any s. This implies
d
dsE(˜y(s), ˙˜y(s)) = 0.
On the other hand, dsdE(˜y(s), ˙˜y(s)) = − ˙˜y2(s). Hence, we get ˙˜y(s) ≡ 0. This again implies
˜
y(s) ≡ ˆy for some constant ˆy. Thus, (ˆy, 0) is an equilibrium state of the damping oscillation system (9.5). However, the only equilibrium state for (9.5) is (0, 0) because V has a unique minimum and thus the only zero of F := −V0 is 0. This implies
E(˜y(s), ˙˜y(s)) = α = 0.
We conclude that
E(y(t), ˙y(t)) → α = 0 as t → ∞.
Step 4. From step 3,
E(y(t), ˙y(t)) = 1
2y(t)˙ 2+ V (y(t)) → 0 as t → ∞.
and V (y) ≥ 0, we get
˙
y(t) → 0 and V (y(t)) → 0, as t → ∞.
Since 0 is the unique minimum of V , we get that V (y) → 0 forces y → 0.
9.2. LIAPUNOV FUNCTION AND GLOBAL STABILITY 171 The above method to show global stability is called the Lyapunov method. The energy function E above is called a Lyapunov function. Thus, the effect of damping (dissipation) is a loss of energy.
In the active shock absorber:
m¨y = −β( ˙y) ˙y − ky, β(v) = v4,
the equilibrium state is (0, 0). From Lyapunov method, we see that this equilibrium is globally stable.
For the simple pendulum, we see that V (θ) = −g/l cos θ has infinite many minima: θ = 2nπ.
The function E(y, ˙y) has local minima (2nπ, 0). The local minimum (2nπ, 0) sits inside the basin Bn= {(y, ˙y) | E(y, ˙y) < g/l}.
The equilibrium (2nπ, 0) is the only minimum of E in the basin Bn. Suppose a solution starts from a state (y(0), ˙y(0)) ∈ Bn, then by using the Lyapunov method, we see that (y(t), ˙y(t)) → (2nπ, 0) as t → ∞.
What will happen if E(0) ≥ g/l initially? From the loss of energy we have E(t) will eventually go below g/l. Thus, the trajectory will fall into some basin Bn for some n and finally goes to (2nπ, 0) as t → ∞.
Homeworks
1. Plot the phase portrait for the damped simple pendulum (9.2).
2. Consider a simple pendulum of length l with mass m at one end and the other end is attached to a vibrator. The motion of the vibrator is given by (x0(t), y0(t). Let the angle of the pendulum to the verticle axis (in counterclockwise direction) is θ(t).
(a) Show that the position of the mass m at time t is (x(t), y(t)) = (x0(t)+l sin θ(t), y0(t)−
cos θ(t)).
(b) Find the velocity and acceleration of m.
(c) Suppose the mass is in the uniform gravitational field (0, −mg). Use the Newton’s law to derive the equation of motion of m.
(d) Suppose (x0(t), y0(t)) is given by (0, α sin(ω0t)). Can you solve this equation?
3. B-D, pp. 502: 22
9.2 Liapunov function and global stability
We recall that when the perturbation of a hamiltonian system is dissipative, we observe that the hamiltonian H decreases along any trajectory and eventually reaches a minimum of H. If there is only one minimum of H, then this minimum must be globally asymptotically stable. That is, every trajectory tends to this minimum as t → ∞. So, the key idea here is that the globally asymptotical
stability of an equilibrium is resulted from the the decreasing of H. This idea can be generalized to general systems. The dissipation is measured by so called the Liapunov function Φ, which decreases along trajectories. More precisely, let consider the general system
x˙ = f (x, y)
˙
y = g(x, y) (9.8)
Suppose (0, 0) is an equilibrium of this system. We have the following definition.
Definition 9.4. A C1-functionΦ(x, y) is called a Liapunov function for (9.8) if (i) Φ(0, 0) = 0, Φ(x, y) > 0 for (x, y) 6= (0, 0).
(ii) Φ(x, y) → ∞ as |(x, y)| → ∞.
(iii) ˙Φ := Φx(x, y)f (x, y) + Φy(x, y)g(x, y) < 0 for (x, y) 6= (0, 0).
Remark
• Condition (i) says that (0, 0) is the only isolated minimum of Φ.
• Condition (ii) says that the region Φ(x, y) ≤ E is always bounded.
• Condition (iii) implies that along any trajectory dΦ(x(t), y(t))
dt < 0, (9.9)
unless it reaches the equilibrium (0, 0). Thus, Φ(x(t), y(t)) is a decreasing function.
Theorem 9.2. Consider the system (9.8). Suppose (0, 0) is its equilibrium. Suppose the system possesses a Liapunov functionΦ, then (0, 0) is globally and asymptotically stable. That is, for any trajectory, we have
t→∞lim(x(t), y(t)) = (0, 0).
Proof. We shall use the extremal value theorem to prove this theorem. The extremal value theorem states that
a continuous function in a bounded and closed domain in Rnattains its extremal value.
Along any trajectory (x(t), y(t)), we have that Φ(x(t), y(t)) is decreasing (condition (iii)) and bounded below (condition (i)). Hence it has a limit as t tends to infinity. Suppose limt→∞Φ(x(t), y(t)) = m > 0. Then the orbit (x(t), y(t)), t ∈ (0, ∞) is confined in the region S := {(x, y)|m ≤ Φ(x, y) ≤ Φ(x(0), y(0))}. From condition (ii), this region is bounded and closed. Hence dΦ(x(t), y(t))/dt can attain a maximum in this region (by the extremal value theorem). Let us call it α. From (9.9), we have α < 0. But this implies
Φ(x(t), y(t)) = Z t
0
dΦ(x(t), y(t))
dt dt ≤ αt → −∞ as t → ∞.
9.2. LIAPUNOV FUNCTION AND GLOBAL STABILITY 173 This is a contradiction. Hence limt→∞Φ(x(t), y(t)) = 0.
Next, we show (x(t), y(t)) → (0, 0) as t → ∞. Let ρ(t) = x(t)2+ y(t)2. Suppose ρ(t) does not tend to 0. This means that there exists a sequence tnwith tn → ∞ such that ρ(tn) ≥ ρ0 > 0.
Then the region
R := {(x, y)|x2+ y2 ≥ ρ0and Φ(x, y) ≤ Φ(x(0), y(0))}
is bounded and closed. Hence, by the extremal value theorem again that Φ attains a minimum in this region. Since Φ > 0 in this region, we have
minR Φ(x, y) ≥ β > 0.
and because (x(tn), y(tn)) ∈ R, we obtain mintn
Φ(x(tn), y(tn)) ≥ β > 0.
This contradicts to limt→∞Φ(x(t), y(t)) = 0. Hence, x2(t) + y2(t) → 0 as t → ∞. Thus, we obtain that the global minimum (0, 0) is asymptotically stable.
If the Lyapunov function Φ satisfies additional conditions:
(iv) The condition (iii) is replaced by ˙Φ(x, y) ≤ −αΦ(x, y) for some positive constant α for all (x, y),
(v) Φ ∈ C2 andΦ(x, y) ≥ c(x2+ y2) in a neighborhood of (0, 0).
Theorem 9.3. Under the assumptions (i),(ii),(iv) (v), the state (0, 0) is asymptotically stable and any solution|(x(t), y(t))| = O(e−αt) as t → ∞.
If the Lyapunov only satisfies the following weaker condition, then we can only have stability result, not asymptotic stability.
Definition 9.5. A C1-functionΦ(x, y) is called a (weak) Liapunov function for (9.8) if (i) Φ(0, 0) = 0, Φ(x, y) ≥ 0 for (x, y) 6= (0, 0).
(ii) Φ(x, y) → ∞ as |(x, y)| → ∞.
(iii) Φx(x, y)f (x, y) + Φy(x, y)g(x, y) ≤ 0 for (x, y) 6= (0, 0).
Theorem 9.4. Consider the system (9.8). Suppose (0, 0) is its equilibrium. Suppose the system possesses a weak Liapunov functionΦ, then (0, 0) is stable.
Example. Damped simple pendulum.
θ =¨ g
l sin θ − b ˙θ
Here, b > 0 is the damping coefficient. In the form of first order equation, it reads (x, y) = (0, 0). So, it only satisfies condition of the weak Lyapunov function. Thus, we can only get a stability result, not asymptotical stability result. However, suppose we consider the linear problem, say the spring-mass system with a linear damper. We know the solutions decay to (0, 0) state exponentially fast from explicit solution formula. Such result cannot be obtained via the weak Lyapunov function. There are two ways to solve this. One is we modify the Lyapunov function.
The other is we provide another linear stability theory based on perturbation theory.
9.2.1 Construction of Lyapunov functions
Lyapunov function for Linear Stable System Consider the linear system
˙x = Ax.
9.2. LIAPUNOV FUNCTION AND GLOBAL STABILITY 175 We also want to show that
xTQx ≤ −αxTPx.
This is equivalent to
Φ ≤ −αΦ,˙ and leads to
Φ(x(t)) ≤ Φ(x(0))e−αy. This gives exponential convergence of all solutions to (0, 0).
Homework For the damped spring-mass system
m¨y + γ ˙y + ky = 0.
Find a Lyapunov function Φ that gives exponential convergence result.
9.2.2 Stability result for perturbed equations Theorem 9.6. Consider the nonlinear equation
y0= f (y)
Suppose ¯y is an equilibrium of the nonlinear equation, i.e. f (¯y) = 0. If ¯y is an exponentially stable equilibrium for the linearized ODE:
y0= f0(¯y)(y − ¯y), that is
Re(λ(f0(¯y))) < 0.
then ¯y is also an exponentially stable equilibrium of the nonlinear squation.
Proof. Let us write u = y − ¯y. Then u satisfies the perturbed equation u0 = Au + g(u),
where A = f0(¯y) and
g(u) = f (¯y + u) − f (¯y) − f0(¯y)u = O(|u|2).
Since Re(λ(A)) ≤ −α < 0, for any positive definite matrix Q, there exists a positive definite matrix P such that
ATP + PA = −Q.
Now, we define the Lyapunov function
Φ(u) := uTAu
If u(t) is a solution of the perturbed equation, then Φ(u(t)) = ˙˙ uTPu + uTP ˙u
= (uTAT + g(u)T)Pu + uTP(Au + g(u))
= uT(ATP + PA)u + g(u)TPu + uTPg(u)
= −uTQu + g(u)TPu + uTPg(u)
≤ −uTQu + O(|u|)uTPu
The last line is due to g(u) = O(|u|2). Now, we choose any α0 with 0 < α0 < α. We then choose δ > 0 such that α0+ δ ≤ α. With these, we have
uTQu ≥ αuTPu ≥ (α0+ δ)uTPu.
Thus, when |u|∞< δ, we have
Φ(u(t)) ≤ −u˙ TQu + O(|u|)uTPu
≤ −(α − δ)uTPu ≤ −α0Φ(u(t)).
Thus, Φ(u(t)) = O(e−α0t).
Here, I have used an a priori estimate of u, namely |u(t)| ≤ δ. This is valid if we choose |u(0)|
small enough. I leave you to fill in the gap.
Project Study the damper of Taipei 101.
9.3 Poincar´e-Bendixson Theorem
We still focus on two-dimensional systems
y0 = f (y), y(0) = y0 (9.10)
where y ∈ R2. As we mentioned, our goal is to characteristized the whole orbital structure. We have seen the basic solutions are the equilibria. The second class are the orbits connecting these equilibria. In particular, we introduce the separatrices and the homoclinic orbits. We have seen in the damped pendulum that solutions enclosed in separatrices go to a sink time asymptotically. In this section, we shall see the case that the solution may go to an periodic solution. In other words, the solution goes to another separatrix. The van de Pol oscillator and the predator-prey system are two important examples.
We first introduce some basic notions. We denote by φ(t, y0) the solution to the problem (9.10).
The orbit γ+(y) = {φ(t, y)|t ≥ 0} is the positive orbit through y. Similarly, γ−(y) = {φ(t, y)|t ≤ 0} and γ(y) = {φ(t, y)| − ∞ < t < ∞} are the negative orbit and the orbit through y. If φ(T, y) = y and φ(t, y) 6= y for all 0 < t < T , we say {φ(t, y)|0 ≤ t < T } a periodic orbit. A point p is called an ω (resp. α) point of y if there exists a sequence {tn}, tn→ ∞ (resp. −∞ ) such
9.3. POINCAR ´E-BENDIXSON THEOREM 177 that p = limn→∞φ(tn, y). The collection of all ω (resp. α) limit point of y is called the ω (resp.
α) limit set of y and is denoted by ω(y) (resp. α(y)). One can show that ω(y) = \
t≥0
[
s≥t
φ(s, y)
Thus, ω(y) represents where the positive γ+(y) ends up. A set S is called positive (resp. negative) invariant under φ if φ(t, S) ⊂ S for all t ≥ 0 (resp. t ≤ 0). A set S is called invariant if S is both positive invariant and negative invariant. It is easy to see that equilibria and periodic orbits are invariant set. The closure of an invariant set is invariant. Further, we have the theorem.
Theorem 9.7. ω(y) and α(y) are invariant.
Proof. The proof is based on the continuous dependence of the initial data. Suppose p ∈ ω. Thus, there exists tn → ∞ such that p = limn→∞φ(tn, y). Consider two solutions: φ(s, p) and φ(s + tn, y) = φ(s, φ(tn, y)), for any s > 0. The initial data are closed to each other when n is enough.
Thus, by the continuous dependence of the initial data, we get φ(s, p) is closed to φ(s + tn, y).
Theorem 9.8 (Poincar´e-Bendixson). If γ+(y) is contained in a bounded closed subset in R2 and ω(y) 6= φ and does not contain any critical points (i.e. where f (y) = 0), then ω(y) is a periodic orbit.
Chapter 10
Numerical Methods for Ordinary Differential Equations
10.1 Two simple schemes
We solve the initial value problem
y0= f (t, y), y(0) = y0. (10.1)
Numerical method is to approximate the solution y(·) by yn ∼ y(tn), where t0 = 0 < t1 < · · · tn are the discretized time steps. For simplicity, we take uniform step size h. We define tk = kh. We want to find a procedure to construct yn+1from the knowledge of yn. By integrating the ODE from tnto tn+1, we get
y(tn+1) = y(tn) + Z tn+1
tn
f (t, y(t)) dt
So the strategy is to approximate the integral by numerical integral hFh(tn, yn).
Below, we give two popular methods 1. Forward Euler method
yn+1 = yn+ hf (tn, yn)
2. Second-order Runge-Kutta method (RK2)
y1 = yn+ hf (tn, yn), yn+1 = yn+1
2h(f (tn, yn) + f (tn+1, y1))
= 1
2(y1+ (yn+ hf (tn+1, y1)) 179
10.2 Truncation error and orders of accuracy
In the forward Euler method, we can plug a true solution y(t) into the finite difference equation, by Taylor expansion, we get
y(tn+1) = y(tn) + hf (tn, y(tn)) + (h) (10.2) where the error term (h) is obtained by
(h) = h Z tn+1
tn
f (t, y(t)) dt − hf (tn, y(tn)) = h Z tn+1
tn
(f (t, y(t)) − f (tn, y(tn))) dt = O(h2).
The error term (h) is called the truncation error. You may view the forward Euler method is a rect-angle method for numerical integration forRtn+1
The error term (h) is called the truncation error. You may view the forward Euler method is a rect-angle method for numerical integration forRtn+1