Stability analysis for sets of solutions

Nonlinear Analysis, Theory, Methods&Applications, Vol.24, No. 6, pp. 801-823, 1995 Copyright© 1995ElsevierScienceLtd Printedin Great Britain.All rightsreserved 0362-546X/95 $9.50+ .00

Pergamon

0362-546X(94)00146-4

STABILITY ANALYSIS FOR SETS OF SOLUTIONS E. A. GALPERINt and V. V. RUMYANTSEV~§ tD6partement de math6matiques et d'informatique, Universit6 du Qu6bec h Montr6al, C.P. 8888, Succ. A, Montr6al, Qu6bec H3C 3P8, Canada; and ~Computing Center, Russian Academy of Sciences, 40, Vavilova Street, 117967, Moscow, Russia (Received 1 August 1993; received for publication 4 July 1994) Key words and phrases: Stability, generalized perturbation equation, contraction analysis.

1. I N T R O D U C T I O N

Given a nonlinear system dx d---t = f ( x , t),

x c R",

t _> to,

(1.1)

under standard conditions ensuring the existence, uniqueness and extendability of solutions in some region A c_ R" of initial data, the Lyapunov stability theory [1] can be applied to investigate stability of a certain particular solution of interest that corresponds to given initial data x°(t) = X(Xo, to, t),

x°(to) = Xo ~ Ao c_ A.

(1.2)

Stability of solutions (1.2) is studied with the use of the perturbation equation which is obtained from (1.1), (1.2) by the transformation (1.3)

x(t) = x°(t) + w(t).

Substituting (1.3) into (1.1) and assuming the function f in (1.1) to be analytic with respect to x, one can use the expansion dx°(t)

- -

dt

+

dw

-d-f

= f ( x ° + w, t) = f ( x °, t) +

7f(x °, t ) w + g ( w , t)

(1.4)

yielding, after cancellation of the first terms, the perturbation equation dw = A ( t ) w + g(w, t), dt

--

g(O, t) - O,

t >_ to.

(1.5)

Here A ( t ) is the Jacobian matrix o f f ( x , t), (1.1), calculated on the solution x°(t), (1.2), and g(w, t) are higher order terms with all partial derivatives calculated on the same solution (1.2). According to (1.3), the unperturbed motion x°(t) of (1.2) corresponds to the trivial solution w(t) =- 0 of the perturbation equation (1.5). This allows us to substitute the problem of stability of the motion x°(t), (1.2), of the nominal equation (1.1) by the problem of stability of trivial § During this research, this author was a visiting scholar at the Universit6 du Qu6bec h Montr6al from the Russian Academy of Sciences under the Canada-USSR (CIS) Academic and Scientific Exchange Program. 801

802

E . A . GALPERIN and V. V. RUMYANTSEV

solution w = 0 of the perturbation equation (1.5). This approach has led to the powerful and elegant methods that constitute the classical stability theory, see, e.g. [1-13] and further references therein. Consideration of perturbation equation (1.5) with all its comfort of using linear approximation d w / d t = A ( t ) w and then, if necessary, successive higher order terms (critical cases) has, however, some specific qualities. First, if a particular solution x°(t), (1.2), is not given as an explicit function of Xo, to, t (i.e. as a formula), then perturbation equation (1.5) cannot be determined. Second, if the solution (1.2) and, thus, the perturbation equation (1.5) are known, then the results of stability analysis on that basis are applicable to that particular solution only.

2. G E N E R A L I Z E D

PERTURBATION

EQUATION

Consider a certain solution x°(t) -- x(x 0, t o , t) of (1.1) corresponding to some fixed initial data x0 e A0 and some other solution x*(t) -- x(x., to, t) of the same equation (1.1) corresponding to another initial data x . e Ao. The deviation w(t) = x*(t) - x°(t),

(2.1)

W(to) = x , - Xo = Wo,

is governed by the equation dw

-~ = f(x

0

+ w, t) - f ( x °, t) = q(w, x °, t),

t _> to

(2.2)

that can be obtained from (1. l) with substitution x* = x ° + w, same as (1.3). In contrast with equations (1.4), (1.5), the composite function q in (2.2) contains an unknown solution x°(t) as its argument. On the other hand, q(0, x °, t) - 0 for all x°(t), t, thus, w(t) - 0 is the solution of (2.2) for any x°(t). It means that trivial solution w - 0 can be put in correspondence to any particular solution x°(t), serving, therefore, the whole region A 0 of possible initial data. If f ( ' ) in (1.1) is analytic with respect to x, then q(-) of (2.2) is analytic with respect to w, yielding the generalized perturbation equation dw

d---t = A ( x ° ( t ) ' t)w + g(w, x°(t), t),

g(0, x°(t), t) - 0,

t > to.

(2.3)

For some particular x°(t), it is, of course, identical to (1.5) with corresponding A(t), g(w, t), where we use the same notation A, g for different functions. However, without fixing x°(t), it represents a bundle of equations given on a continuum of different particular solutions. With this meaning, we shall drop sometimes the indication of a particular solution, writing simply dw m

dt

A ( x , t)w + g(w, x, t),

t _ to,

(2.4)

with the understanding that (2.4) is a corresponding perturbation equation for every solution x(t) of (1.1). It means that the form (2.4) is conserved while the terms are different for different

x(t).

Stability analysis 3. S T A B I L I T Y

ANALYSIS

803

FOR SETS OF SOLUTIONS

In stability analysis, deviations w(t) are studied in a neighborhood H of the origin and one is interested to determine whether or not for every r / > 0 there exists fi(r/) > 0 such that if IIw011 ~ '~('7),

(3.1)

then [[w(t)[[ < q,

for all t > to,

(3.2)

where 11-[[ is the Euclidian norm. If the answer to this question is in the affirmative, then the motion w(t) - 0 is called stable, otherwise, unstable. It means that if there exists r/0 > 0 such that, whatever small ~ > 0 may be, there is a moment t . > to at which [1w(t.)[[ = r/0, then the motion is unstable. If at some moments t* > to, perturbations grow to a fraction of the magnitude of a nominal coordinate, [wj(tD[ -- oglxj(tD[, 09 = const _> 1, 1 _< j _< n, then the motion is unstable. A stable motion with the additional property lim I]w(t)ll -- 0

(3.3)

t~eo

is called asymptotically stable. These are the classical definitions of stability given by Lyapunov [1]. With the notation (2.1), it refers, of course, to the stability of the solution x°(t). In this research we do not fix the initial condition x°(to) = Xo ~ Ao, considering instead a collection of nominal solutions [x°(t)] = x°([x0], to, t) corresponding to a set {Xo] c A o of initial conditions; the notation Ix0] may mean a finite collection or a set, a continuum. To study and solve the problem by Lyapunov's second (direct) method, C~-functions V(w, t), W(w), W~(w), W*(w) are considered that vanish if w = 0, V(0, t) = W(0) = W~(0) = W*(0) -- 0,

t _ t o,

(3.4)

and have some additional properties. Recall the basic theorems of Lyapunov's second method. THEOREM 3.1 (Lyapunov, [1]). If there exists a function V(w, t) satisfying the conditions

(a)

V(w, t) >_ W(w) > O,

w ~ H,

dV OV --~ (z.z) = 0----[+ V V. q(w, x °, t) <_ O,

(b)

w # 0,

w e H,

(3.5)

t _> to;

(3.6)

t >_ to,

then the solution w(t) =- 0 is stable. THEOREM 3.2 (Lyapunov, [1]). If there is a function V(w, t) satisfying condition (a) and (c) (2,2)

(d)

=--+ Ot

VV.q(w,x°,t)<.-W~(w)_ V(w, t),

w e H,

then the solution w(t) ~ 0 is asymptotically stable.

weH, t ___ to,

w~0,

t>-to;

(3.7) (3.8)

804

E.A. GALPERIN and V. V. RUMYANTSEV

TnEOREM 3.3 (Chetaev, [2]). If there exists a function V(w, t) satisfying the conditions: (e)

the set 8t~ = {w ~ HI V(w, t) > O, t >_ tol N IIIwll < ,I, ,7 > ol ~ ®

(3.9)

is nonempty for all t >_ to and any small r / > 0;

V(w, t) is bounded within 8t~;

(f) (g)

dV OV -d-i (2.2) = O--td-

V V"

(3.10)

t w e 8~,

q(w, x °, t) > 0,

w ~ 0,

(3.11)

meaning that d V/dt is positive definite in 8t~, in the sense that for every small e > 0 there is ), > 0 such that if V(w, t) >_ e, then dV (2.2)

~-~ ~)

for all t _> t 0,

(3.12)

and the motion w(t) - 0 is unstable. Geometrically, condition (3.11) together with (3.12) mean that if w ( t ) ~ 8t~ is uniformly separated from the boundary OSt~ for all t _> to, then d V / d t _> ), > 0 is uniformly separated from zero for all t _> to, see [2, Section 13]. As distinct from (3.5), a function V(w, t) in (3.9) need not be positive definite. Of course, stability, asymptotical stability or instability o f the solution w(t) - 0 implied by theorems 3.1-3.3 means the same property of all nominal solutions [x°(t)l for which (3.6), or (3.7), or (3.11)-(3.12), respectively, are fulfilled. To demonstrate how those theorems work with solutions x°(t) included in (3.6), (3.7), (3.11), consider an example from [3, Sections 4, 44].

Example 3.1. Consider the equation dx - dt

-

x(a2

-

x2)'

O( > 0,

t > t 0.

(3.13)

According to (2.1), (2.2), we have fv=(x+

w)[a 2 - ( x +

3xw 2 - w 3,

w)2 ] - x ( o ~ 2 - x 2 ) = ( a 2 - 3 x E ) w -

t>_t o

(3.14)

which is the generalized perturbation equation (2.4) in our case of (3.13). With (3.14) we can do the standard stability analysis for (3.13) as follows. Equation (3.13) has three stationary solutions x~ -= 0, x2.3 = _+a. Substituting those solutions in (3.14), we immediately obtain instability for x~ -- 0 and asymptotic stability for x2,3 = _+or, all by the first approximation in (3.14). These results can also be established by considering the Lyapunov function V -- ½w2 which has the following derivative on trajectories of (3.14)

f, = diE_ lag

= w2(o/2 -

3X 2 --

3xw

-

w2).

(3.15)

(3.14)

For Xl ~- 0, we have from (3.15) that f" > 0 if ot2 - w 2 > 0, asserting instability and yielding domain of repulsion w e (-c~, c¢) with respect to nominal solution Xl(t) - 0.

Stability analysis

805

For x2,3 = _+o~, we have f r o m (3.15)

V2,3 = w2( -2Ot2:1: 3 o t w -

w2),

(3.16)

asserting asymptotic stability of both solutions for small w. To find domain of attraction for xz = or, we take the upper sign in (3.16) and solve the inequality w E + 3o~w + 2or2 > 0, yielding w > -ct or w < -2c~, which in coordinates tOx corresponds to x > 0 or x < -o~ since in this case w = x - x2 = x - o~. However, in the region x < - a there is another attractor, namely, x3 = -or; hence, domain of attraction for x 2 - - ~ is x e (0, oo). For x3 = -o~, the same arguments with the lower sign in (3.16) yield domain of attraction x e (-oo, 0); details are left to the reader. We see that generalized perturbation equation can be used for all known solutions of the nominal equation. Moreover, it can be used for stability analysis of solutions that cannot be expressed as explicit integrals and for which one cannot write the specific perturbation equation (1.5) corresponding to a particular solution x°(t), see (1.2)-(1.3), n o t given as a f o r m u l a . In such cases, the generalized perturbation equation represents a new and important tool for stability analysis.

4. C O N T R A C T I V E AND REPULSIVE DOMAINS Consider x ° in (2.2), (2.3) and x in (2.4) not as a particular solution, but as a parameter. Then inequalities (3.6), (3.7), (3.11) become characteristics of a domain (simply connected open set) E = ~D × (t', t' + T),

~D c_ A c_ R n,

t' > t o fixed,

T > 0,

(4.1)

where ~3 m a y vary with t e (t', t' + T). With x, t considered as independent variables, the left-hand side of (3.6), (3.7), (3.11) becomes a function F: R n × R n × R ~ R of three arguments OV F(w, x, t) = - ~ + V V . q(w, x, t),

(4.2)

which coincides with the total derivative lk = d V / d t of a chosen function V(w, t) on trajectories x(t) of the nominal equation (1.1). Consideration of such functions (4.2) and domains (4.1) is motivated by the need to evaluate the rate of attraction of perturbed motions to a nominal solution of (1.1) within a finite time interval, and for all nominal trajectories passing through domain E of (4.1). For processes evolving in a finite space-time region, such information may be useful irrespective of stability properties on [to, oo). In such considerations, perturbations w do not have to be small. Definition 4.1. I f for a chosen V(w, t) satisfying (3.5) on an interval (t', t' ÷ T), the condition (3.6) or (3.7) holds for (x, t) e E, then domain E is called neutral or contractive, respectively. Definition 4.2. I f for a chosen V(w, t) satisfying (3.9), (3.10) on an interval (t', t' + T), the condition (3.11) holds for (x, t ) e E , then domain E contains a repulsive sector ~t; such domain, E, is called repulsive.

806

E . A . GALPERIN and V. V. RUMYANTSEV

The statement that a certain domain E is contractive, neutral or repulsive means that there is a function V(w, t) mentioned in definitions 4.1, 4.2 which renders the corresponding property of E. The availability of such a function defines the corresponding domains. For example, if V(w, t) satisfies (3.5), (3.7) for all t _> to, then our domain becomes a contractive band E = ~ x [to, oo) with one sole Lyapunov function which is the classical case.

Remark 4.1. The names contractive or repulsive domain relating to the (x, t)-space should not be confused with the names domain o f attraction or repulsion relating to the w-space, as in example 3.1 above. To illustrate the geometry corresponding to definitions 4.1, 4.2, we can use the standard argument o f the Lyapunov stability theory [1, 2]. Consider, for example, a neutral domain E1 = $31 x [to, tl). For a given r / > 0, let 71 =

inf W1; Ilwll= n

due to (3.5), 7~ > 0.

(4.3)

Since ~(w, to) does not depend on t, so due to (3.4) and to the continuity of II1, there is fi > 0 such that for IIw011 ~ ~ we have Vl(wo, to) < 7~. Choosing such initial conditions and due to the relation V1 - Vl(W0, to) =

f'

~ dt,

121 -< 0 as of (3.6),

t e [to, G],

(4.4)

to

we obtain that w(t) is such that the following conditions are satisfied W~ _< ~(w, t) _ ~(w0, to) < Yl,

t • [to, t~l,

(4.5)

implying IIw(t)ll < '1 for t • [to, t~]. It means that, over a neutral domain, perturbations within a ball IIwll < ,7 where (3.4)-(3.6) are satisfied cannot escape this ball whatever (x, t) • E~ = ~ × [to, tt). If t~ = oo, stability follows. If we have strict inequality ~ < 0 in (4.4), compare with (3.7), then domain El is contractive. If t~ = oo and we use the additional condition (3.8), then asymptotic stability follows by the standard argument [1, 2]. However, if we consider two adjacent domains with different functions I/1, I/2 (with one common function it would be one single domain), then neutrality or contractivity of the union does not follow from the same property for component domains. Indeed, continuing the argument (4.3)-(4.5) f o r E 2 ~)2 X [tl, t2), we denote/71 IJW(tl)ll. Clearly, r / > rh > 0 since, otherwise, the value w(tO = 0 of the solution w(t) ~ 0 would contradict the uniqueness of a solution emanating from the point (t~, 0) due to the existence of the trivial solution w(t) =- O. Let =

=

Y2=

inf W2. Ilwll= 71

(4.6)

Since V2(w, tl) does not depend on t so due to (3.4) and to the continuity of V2, there is '~2 > 0 such that for Ilwlll we have ~ ( w l , tl) < Y2. However, wl = w(tO = W(Wo, to, tO comes from El and cannot be chosen so as IIw~ll ~ ~2 for appropriate c~2 > 0. Hence, to continue the argument and to assure that finite or countable union of adjacent neutral (contractive) domains be also neutral (contractive), we have to impose the following condition.

Stability analysis

807

Consistency condition. A sequence of adjacent or overlapping neutral (contractive) domains E l , E2 . . . . with functions V~, V2. . . . . acting on [to, fi), [h, t2) . . . . . and satisfying (3.4)-(3.6), or (3.7) for contractive domains, is called consistent if the functions V~, Vz. . . . are such that, with the initial condition Ilwoll -< for a given r / > 0 in (4.3), we have Vz(w~, fi) < Y2 for W l = W(Wo,to, q) and any x e : D 2 ( q ) , then V3(w2,tz)
Proof. There is a sequence of functions V~, Vz. . . . acting on [to, q), [q, t2) . . . . and satisfying (3.4)-(3.6) that corresponds to a cover by a finite or countable chain of consistent neutral domains. If the chain is finite, we prove the theorem after a number of repetitions of the above argument (4.3)-(4.5) since the last tk = ~ . If the chain is countable, then t, --. oo, thus, for every t e [to, ~ ) there is a subsegment to which it belongs, yielding IIw(t)ll < ~/ for all t>-to. • A solution which is entirely in a contractive band may not be asymptotically stable though its stability follows from theorem 4.1 since (3.6) is implied by (3.7). If the chain is finite and for the last function Vg(w, t) acting on [tt,, oo) the condition (3.8) is satisfied, then asymptotic stability follows from the classical Lyapunov theorem [l]. For a countable chain of consistent contractive domains, consider a sequence of corresponding functions

Vi(w, t), t ~ [ti-l, ti),

t i ~ oo,

i = 1,2 . . . . .

(4.7)

each acting over corresponding domain Ei of finite time length Ati = ti - t~_ 1 -> r > 0. Functions (4.7) may be regarded as components of a piecewise continuous function V(w, t) acting on [to, ~ ) , which components should satisfy the consistency condition stated above. Now, condition (3.8) can be extended onto the sequence (4.7) as follows. From (3.5), (3.8) we have

W*(w) >- V(w, t) >- W(w) > O,

w ~ H,

w # O,

(4.8)

where V(w, t) represents V/(w, t) over each [ti-1, ti) of (4.7). Since W*(w) -* 0 as IIwl[ ~ o, so for appropriate r / > q . > 0 the surface V(w, t) = y is enclosed in the spherical ring r/_> Ilwll -> r/.,

(4.9)

808

E . A . G A L P E R I N and V. V. R U M Y A N T S E V

provided that r / > y > t/, and the ring (4.9) is in the region H. Indeed, it is sufficient to take such ~/, r/. that the sphere Ilwll = ,1 is circumscribed around W*(w) = r/1 _< r/, and the sphere Ilwll = ,1, is inscribed in W(w) = ~/2 -> t/,, T/I > /'/2" Since r/1 = W*(w) --, 0 as Ilwll --' o, we can take r / ~ 0. Vice versa, if (3.8) holds, then for any spherical ring (4.9) in the region H , by virtue of (4.8), (3.8), there exist functions of (4.7) acting over this ring (we say in such case that ring (4.9) is covered by consistent contractive domains). Take a decreasing sequence r / = /~1 > "'" > /~k > /~k+l > . . . . lim q~ = 0, and consider rings R k = {w ~ H I r& >_ Ilwll --- ,7 ÷11, k = 1,2 . . . . . Consider all functions V/(w, t) from (4.7) acting over the ring Rk. By (3.7), every ~ < 0 which means that there exists ~ l ( w ) such that over the segment of definition of E(w, t) we have definite negative and bounded from zero total derivatives - V i ( w , x , t)

>

~/I(w) > 0,

weH,

w#O,

(x,t) eEi=

~i×

[ti-i,ti).

(4.10)

Let Yik =

inf

(4.11)

w / l ( w ) _> Yk > O.

w~R k

The uniform bound Yk > 0 exists for all V/acting over Rk since otherwise ~/l(w) would not be separated from zero within closed Rg not containing zero, in contradiction with definition of a positive definite function. Now, integrating the piecewise continuous function V(w, t) with components (4.7) along a trajectory (or a part thereof) lying entirely within Rk, we obtain by (4.10), (4.11) V-

V(w,, t,) =

I"dt <_ - ~ t,

yik(ti - ti_l) <_ - y k ( t - t,),

(4.12)

i

where the sum covers all components V~(w, t) acting over Rk and t . is the starting time of a perturbed trajectory. From (3.5), (4.12), we get 0 < V(w, t) <_ V(w,, t,) - yk(t - t,), meaning that there is only finite time t - t, within Rk. Since the band is contractive, approaching zero, so that for t > t. + Tk we k = 1, 2 . . . . . there are Vi from (4.7) that act This proves the following theorem,

Yk > O,

(4.13)

_< Tk < o0 during which a trajectory can stay the perturbed trajectory w(t) will leave R k , have [Iw(t)ll < ek÷l. By (3.8), for any ring R~, over that ring, hence limllw(t)H = lim JTk = o. t~ k-~

THEOREM 4.2. If a contractive band is such that for any r / > 0 there is N(v/) such that for all i _ N(r/) functions V/(w, t) of (4.7) satisfy the condition r/_> V/(w, t) > 0, w E H, w # 0, t ~ [ti_l, ti), ti ~ oo, then every solution passing entirely within such a band is asymptotically stable. • R e m a r k 4.2. The above arguments resemble the analysis based on property (A) or (B) in [4, Sections 4, 5], under which there exists a Lyapunov function V(w, t) acting on [to, oo) with sign definite derivative that renders asymptotic stability of certain nominal solution x°(t). However, it may be difficult to find such a function and, if found, it serves one particular solution only. Functions (4.7) may be easier to construct, and they serve all solutions passing through corresponding domains E i. If considered as components of one function V(w, t), this

Stability analysis

809

function, though generally discontinuous, renders, under certain conditions, the same conclusions about stability or asymptotic stability as a classical Lyapunov function. Remark 4.3. In contrast and similarity with vector Lyapunov functions introduced, e.g. in [11], that create a space mosaic based on the idea that each subsequent function (all acting on [to, oo)) covers a manifold (or a part thereof) where preceding functions are inconclusive (e.g. where lk = 0), the functions (4.7) correspond to a time-space mosaic of consistent domains which domains, if forming a band extending over [to, oo), may deliver the same stability properties as a conventional Lyapunov function. 5. R A T E O F C O N T R A C T I O N

ESTIMATES

Functions V~of (4.7) corresponding to a chain of consistent contractive domains can be used to obtain quantitative results concerning the measure of contraction within every domain Ei. For simplicity, suppose that V~ = V~(w), independent of time, since otherwise we would have to use a pair of functions W/*(w) and W/(w) of (4.8) for each V/(w, t). Dropping index i, we consider a generic domain E of (4.1), its closure/~ and the function F(w, x, t) of (4.2) defined over H × / ~ (in our case F = VV. q(w, x, t)). Along a trajectory x(t) passing through E, we have d V / d t = (/" = F(w, x(t), t). Suppose that the range of initial perturbations w(t') = w' for a trajectory x(t) ~ ~ ( t ) , t E [t', t' + T], is given: IIw'll -- ,7'. By definition of a contractive domain, we have F(w, x, t) < 0 for all (x, t) e / ~ -- if) × [t', t' + T] (not only for x(t) along a trajectory). If a nominal solution x°(t) is known, then we have V(w) - V(w') =

F(w(t), x°(t), t) dt.

(5.1)

t'

The use of (5.1) is, however, impractical since it requires the integration of the perturbation equation (1.5) together with nominal equation (1.1). Suppose, first, that V(w) = I[wll z and denote IIw(t)ll = r/(t). Then, from (5.1), we have the following estimates

~lztt)

-

r/, z =

IIw(t)ll z

IIw'il 2

-

max Ilwll

t' ,l(t)

=

F(w, x°(t), t) dt

,l(t)

max

F(w, x°(t), t) dt

-- Ilwll <- ~ '

_< (t - t')

max ~(t) ~ Ilwll <

F(w, x, t) < 0,

0 < t - t' < T,

(5.2)

,1'

(x,t) E E

where all maxima exist due to continuity of functions and compactness of respective closed domains. One takes that estimate from (5.2) which is easier to compute. In general, a sphere Vo(w) = cllwl[2 may not yield Fo(w,x, t) < 0 for all (x, t) E/~. I f / ? is contractive, then there is another function V(w) that yields F(w, x, t) < 0 for all (x, t) ~/~. For this general V(w), estimates (5.2) also hold if I[wll is replaced by V(w) everywhere in (5.2). However, maximization over a level surface V(w) = rtZ(t) or a slice t/2(t) _< V(w) <_ ,7'z for a

810

E . A . GALPERIN and V. V. RUMYANTSEV

general function V(w) is complicated, and computations are dependent on changing functions V/(w), cf. (4.7). To simplify the calculations, we wish to retain minimization over spheres or spherical slices as in (5.2) but apply it to nonspherical functions V/(w). This leads to the ~ * set construction [14]. Define v+=

(5.3)

max V(w) II~ll= ~'

(5.4)

~+ = Iwl v(w) < v+l v-=

min V(w)

(5.5)

~ - = [wl v(w) < v - I ,

(5.6)

where rt' > 0 and r/, 0 < r / < r/', are given. Clearly, v- < v + due to contractive V(w) and r/ < t/'. According to the result in [14, p. 191], for any solution x(t) passing through the region /~ = ~ ( t ) × [t', t' + T], deviations w(t) will shrink from I[wll = ~', at t = c , to I[wll -< ~ < ~' at t = t' + TO or earlier, if I: +

F(w, x, t) <_

--

v-

,

TO

TO <_ T.

(5.7)

Let (01

=

min

min

F(w, x, t)

(5.8)

max

max

F(w, x, t) < 0.

(5.9)

(x,t)~ ~< Ilwll~ '

¢P2 =

(x,t)~£: ~< Ilwll<~'

Then, obviously, (Pl -< (o2 and we have, regarding TO as the actual time of contraction, ~'+

--

-

>

1~21

V-

~'+

-

TO>--

--

It-

(5.10)

I~,1

Note that if r / ~ 0 in (5.5), then v- ~ 0 and (& ~ 0 uniformly with respect to (x, t) e E (due to V(0) = 0 and q(0, x, t) - 0). Parameters r/, TO are linked; if for a chosen r / i n (5.5), time period TO > T by (5.10), then such contraction cannot be achieved within/~ and r/should be increased. The speed of contraction or repulsion at different points (x, t) e E can be estimated by the following relative bounds

qq(x, t)

=

1

min F ( w , x , t ) , ~ Ilwlr=

1

~,~(x, t) = ~ , wmax , = ~ F(w, x, t),

(x,t) e E

(5.11)

(x, t) e E,

(5.12)

where qq _< ~/-/2< 0 for contraction and 0 < ~//2 for repulsion. It should be noted that all min and max in this section are global which in a nonconvex situation requires application of the full global optimization methods, see, e.g. [15].

Stability analysis

811

Example 5.1. Application o f estimates is a p p a r e n t in the simple example 3.1, where we have, according to (3.15), (4.2) for V = l w 2 F(w,x, t) = w2(ot 2 - 3x 2 - 3 x w - w2).

(5.13)

F r o m (5.11), (5.12), we have for w = _+r/ ~q(x) = ~2 _ 3x 2 _ 3xr/ - r/2,

r / > 0,

(5.14)

~-/2(X) = Ot2 -- 3x 2 + 3xr/ - r/2,

r/ > 0.

(5.15)

T o a p p l y estimates (5.7) to (5.10), let us take a b a n d [Xl,X2] × [to, oo), x2 > xl > r/' > 0. Then we have (,01 =

t/'2(Ot 2 --

3x 2 - 3XEr/' - r/'2)

(fiE = /72( 0/2 -- 3x 2 - 3Xlr/ -

/12) (

(5.16) O.

(5.17)

In order that (fi2 < 0 for all t / > O, we require that a 2 - 3x 2 _< O, i.e. x~ __ a / v ~ yielding d o m a i n o f contraction in the u p p e r half plane as [ct/v~, oo) different f r o m d o m a i n o f attraction (0, oo) in the same half plane obtained in example 3.1. F o r V = { w 2, we have, according to (5.3)-(5.6), V+ = ~ r1 / p2, ~')+ = {W I W 2 < ~ , 2 j ,

v - = ~ r1 / , 2

(5.18)

~')- = iW I W 2 < /~2j,

(5.19)

and inequality (5.10) yields (r/'/r/) 2 -

1

1 -

(r//r/') 2

la 2 _ 3x12 _ 3xlr / _ r/El _> 2TO _> [Ot2 _ 3X~ -- 3XEr/' -- r/'21 "

(5.20)

With r/, r/' small in c o m p a r i s o n with x~, x2, we can drop respective terms in (5.20), writing the inequality with small uncertain ~ > 0, ~2 > 0 as follows (r/'/r/) 2 - 1 3x 2 _ tx2

1 - (r/It/') 2 ~1 >-- 2To _> 3x 2 _ cx2

(2,

x2 > xl > et/x/3.

(5.21)

T h e factor 2 in (5.20), (5.21) is implied by the choice o f Euclidean n o r m in V = c[[wll 2, c > 0, in our example V = cw 2. The last inequality is convenient for practical estimation o f the m e a s u r e o f contraction existing in various locations o f the (x, O-plane. Let, for example, x~ = ct, x2 = 2a; then the b a n d E = [ct, 2a] × [to, ~ ) contains whole trajectories for which decrease o f perturbations by half ( t / ' / r / = 2) is attained in time period To estimated as follows 0.75 3 3 0.034 cz2 - 4ot2 > T o > - - - 8 8 a 2 - ~ 2 a2

~2>0,

(5.22)

where we d r o p p e d correction t e r m ~ o f (5.21). F o r ct _> 1, we have fast contraction in E, and for a _< 0.1 the contraction is relatively slow. If r/, r/' are c o m p a r a b l e with x~, x2 (say, at 25%), then (5.20) should be used. For example, if v / ' = 0.5o~, calculation by (5.20) yields 0.54 > otzTO > 0.026 where equality is lifted due to r o u n d - o f f errors.

812

E . A . G A L P E R I N and V. V. R U M Y A N T S E V

One m a y notice that the last estimate in (5.2) coincides with the left estimate in (5.10). H o w e v e r , two middle estimates in (5.2) would give m u c h better values for transient period To than (5.10), (5.21); this is due to the k n o w n solution x°(t) in (5.2). T o illustrate h o w it works, we use F(w, x, t) f r o m (5.13) in the first two estimates of (5.2) which in our one dimensional example coincide and yield the equality since in b o t h integrals m a x F = F(r/(t), x°(t), t). We have r/2(t ) _ r/,2 =

f t r/2(a 2 -

3x 2 - 3xr/ - r/E)dt,

r/(t) > 0.

(5.23)

¢,

This is equivalent to differential equation 217 dr/ = r/a[c~2 - 3x2(t) - 3r/x(t) - r/El dt,

t >- t', r/(t') = r/'.

(5.24)

For n o m i n a l solution x(t) = 0, equation (5.24) is easily integrated, yielding 1 , r/2(od2 - r/t2) t - t' = ~-~ m r/,2(o~2 _ r/E),

t' t >--

~/'. >-- t 0, r/(t') =

(5.25)

It provides exact time To = t - t' for a given decrease o f perturbations. T o c o m p a r e it with (5.22), let us take o~ = 1, r/' = 4, r / = 2, so as r / ' / r / = 2 as in (5.22), and there is no stationary solution x(t) -- const o f (3.13) in the b a n d [2, 4] x [to, oo). According to (5.25), we obtain To = In

4(1 - 16) 16(1 - 4)

- In 1.25 = 0.223,

which is within the bounds o f (5.22). For n o m i n a l solution x(t) = o~, equation (5.24) becomes 2 r / d r / = r/2[-2ot2 - 3otr/ - r/2] dt,

t >_ t' >_ to, r/(t') = r/'.

(5.26)

With ot = 1, r/' = 0.2, r/ -- 0.1, meaning the same 50°7o p e r t u r b a t i o n decrease as in (5.22), we can integrate (5.26) to obtain exact period To for such decrease. I f we a p p r o x i m a t e (5.26) as d r / = -t~2r/dt, then for those data we get To = ln(r/'/r/) = In 2 = 0.69, within the b o u n d s o f (5.22). I f we use a m o r e precise a p p r o x i m a t i o n , d r o p p i n g only r/2 in the bracket o f (5.26), we obtain the integral 1 2or 2 + 3ctr/ t = ~-~ In + const, r/

(5.27)

which for ot = 1, r/' = 0.2, r / = 0.1 yields •

To = m

r/'(2 + 3r/)

r/(2 + 3r/')

= In 1.77 = 0.57,

(5.28)

also in the b o u n d s o f (5.22). In general case o f multidimensional systems, integral estimates with m a x o p e r a t i o n in (5.2) provide a measure o f p e r t u r b a t i o n decrease in terms o f contracting surfaces v = cllwll2 = r/2(t), without integration o f the vector p e r t u r b a t i o n equation for a given system. Estimates (5.7), (5.10), and the last estimate o f (5.2) eliminate also the integration o f the n o m i n a l equations o f a system at the expense o f less precision which m a y still be acceptable.

Stability analysis

813

Speed of contraction bounds (5.11), (5.12) evaluate the quality of different nominal solutions relative to perturbation behavior not only asymptotically (stability) if E is a band, but also instantaneously and in a bounded space-time region E. With the notions of contractive, neutral or repulsive domains and bands, one can make global analysis of perturbations behavior. To illustrate the situation, consider again (5.13). In order to have F(w, x, t) < 0, we require W2 +

3xw + 3x 2 - a z > 0.

(5.29)

If discriminant 4or 2 -- 3X2 < 0 , that is, Ix] > 2a/v~, then in those regions we have I2 < 0. It means that in the regions x _> 1.1548a > 2a/x/3 or x _< -1.1548ot < -2or/x/3, every nominal solution, while passing in one of those regions, attracts its every perturbed solution (within the same region) with certain speed of contraction. This contraction is global in the sense that perturbations may be of any magnitude if the perturbed motion remains in the same region. Analysis of small perturbations is much simpler. Considering the term ot2 - 3x 2 in (3.15) or (5.29), one can conclude that every solution passing in the region Ix[ -< 0.57ot < a / v ~ is repulsive for small w(t), i.e. repels every perturbed solution near nominal one, and every solution passing in the region x _> 0.58a > or/x/3 or x _< - 0 . 5 8 a < - c t / v ~ is attractive for small w(t), i.e. attracts close neighboring solutions passing within the same region. Combining this with the fact that x(t) ~ ~ sign x 0 as t ~ oo, we conclude that every solution of (3.13) except for x - 0, is ultimately attractive within its half plane in the sense that, given r / > 0, x o # 0 and x . ;~ Xo, X.Xo > 0, there is a moment t,(r/, Xo, x , ) _> t o such that Iw(t)l = Ix(x,, to, t) - X(Xo, to, t)l < ~/

for all t _> t,

(5.30)

and limlw(t)l = 0. t~oO

Ultimate attractivity resembles the notions of the asymptotic stability uniform with respect to coordinates (see [4, property (D) in Section 10, definition 10.1]) or with respect to time t o and coordinates (see [4, property (B) in Section 5, definition 5.1]). A noticeable difference is that we consider any nominal solution x°(t), not fixing certain particular solution (in notations of [4], a special case of property (D) follows from (5.30), if we identify in (5.30) Xo = 0, to = 0, x , = Xo, thus, t,(r/, 0, x0) = T(r/, 0), w(t) = X(Xo, O, t) since x(0, 0, t) = 0 is unperturbed particular solution in [4]). Furthermore, in ultimate attractivity w0 is not confined to regions H o or G~ as in definitions 10.1 and 5.1 of [4]. Due to continuity of solutions with respect to initial data, there exists ~ = ~(r/, x0) > 0 such that Iw(t)l < r/also within [to, t.] if Ix. - x01 < ~, yielding asymptotic stability in the sense (3.1)-(3.2)-(3.3). In ultimate attractivity (5.30), such ¢~ > 0 is missing so that initial deviation Ix. - X o l may be arbitrarily large within the corresponding half plane. Nevertheless, the solution x*(t) is attracted to x°(t) and enters its r/-neighborhood (5.30) at some distant moment t. > to. Such ultimate attractivity under large deviations may be of interest in population dynamics, environmental sciences, and celestial mechanics. In this example, we see that the band E 0 = ( - a / v ~ , ct/x/3) x (to, o0) is repulsive. However, this band contains only one solution, namely, Xl = 0, which is, thus, unstable. Every other solution originating in this band stays in it for a finite time interval. For small perturbations, the bands corresponding to domains (-oo, -a/x/3), (ct/x/3, o0) on the x-axis are contractive and every solution with x o in one of those domains stays in that domain for all t _> t o, being, thus,

814

E . A . G A L P E R I N and V. V. RUMYANTSEV

asymptotically stable in the classical sense (3.1)-(3.2)-(3.3). In contrast, solutions started with 0 < Ixol < ~/v~ enter one of contractive bands only at a moment t,(Ix01), being, thus, ultimately attractive on [to, oo) and asymptotically stable in the classical sense on [t., oo).

6. U S E O F B U N D L E S

OF FIRST INTEGRALS

Chetaev's method of construction of Lyapunov functions in the form of bundles of first integrals [2] (see also [5, Section 10] and further references therein) can be used with the generalized perturbation equation, that is, for stability analysis of sets of solutions. We demonstrate it by the following classical example. Consider the Euler case in the motion of a rigid body around its fixed center of mass without external forces. Equations of such motion are usually written in the form

dp

A - z 7 + ( C - B ) q r = O,

(6.1)

B ~t + (A - C)rp = O,

(6.2)

cUr d t + (B

(6.3)

(ll

- A ) p q = O,

where p, q, r are projections of the vector of angular velocity on coordinate axes taken as principal axes of the ellipsoid of inertia, and A, B, C are principal moments of inertia of the rigid body. Suppose that p°(t),q°(t), r ° ( t ) is some particular solution of (6.1)-(6.3). Substituting p = p0 + ~, q = q0 + r/, r = r ° + ( into (6.1)-(6.3), eliminating terms that are cancelled by virtue of nominal equations (6.1)-(6.3) and dropping the superscript, we obtain the generalized perturbation equations A ~ = (B - C)(rrl + q ( + riO,

(6.4)

Bil = ( C - A ) ( p (

(6.5)

+ r~ + (~),

C ( = (A - B ) ( q ~ + prl + ~rl).

(6.6)

Here upper dot denotes time derivative, and p(t), q(t), r(t) are fixed particular solutions of (6.1)-(6.3) defined by certain initial conditions p(to) = Po, q(to) = qo, r(to) = ro. By inspection, one can see that equations (6.1)-(6.3) have the following first integrals T 1 = A p 2 + B q 2 -k C r 2 =

const,

T2 = A 2 p 2 -k B 2 q 2 q- C 2 r 2 = c o n s t ,

(6.7)

(6.8)

where constants T1, T2 are defined by the vector (Po, qo, ro) of initial velocity. Integrals (6.7)-(6.8) express the laws of conservation of energy and angular momentum, respectively. To investigate stability of nominal solutions with respect to p, q, r, we consider several cases.

Stability analysis

815

Case 1. A = B = C. In this case all solutions are stationary, p = Po, q = qo, r = ro, and all are

stable. Case 2. p = q = r - O. Equations (6.4)-(6.6) coincide with (6.1)-(6.3). Therefore, integrals

T~, T2 with ~, r/, ( instead of p, q, r are also first integrals of perturbed motions. Being positive definite, they can be used as Lyapunov functions to conclude about stability of this trivial solution (at rest). and po2 + q o 2 + r o 2 > 0 . In this case, and taking into account (6.1)-(6.3), generalized perturbation equations (6.4)-(6.6) have the following first integrals C a s e 3. A ; ~ B ; ~ C ; ~ A ,

T~' = A ( p + ~)z + B ( q + r/)2 + C(r + ()2 = const, 7.2. = A 2 ( p + ~)a + B2(q + r/)a + C2(r + O 2 = const,

(6.9) (6.10)

where constants TI*, T2*are defined by initial data Po, qo, ro and initial perturbations Go, r/o, (o. Since T~*, T2* do not vanish at ~ = r / = ( = 0, they cannot be taken as Lyapunov functions. Consider the function V = (T~* - T1)2 + (T2* - T2)2.

(6.11)

This function is nonnegative, V _> 0; vanishes if ~ = r / = ( = 0, and its total derivative on trajectories of perturbed motions (6.4)-(6.6) of the system (6.1)-(6.3) is zero, I)" - 0, since V is a bundle of integrals. If V were positive definite, one would conclude about stability of all motions. Unfortunately, this is not the case. If ~ , r / , ( a r e not all zero, [~[ + [r/] + I([ > 0 , then V = 0 if and only if Tl*= Tl and T2* = T2. To find the manifold on which V = 0, we can write, by virtue of (6.7)-(6.10) T~* - T1 = A ( 2 p ~ + ~2) + B(2qrl + r/2) + C ( 2 r ( + (2) = 0,

(6.12)

T2. _ T2 = A2(2p(

(6.13)

+

~2) q_ B2(2qr/ + r/2) + C2(2r( + (2) = 0.

Denoting parentheses in (6.12), (6.13) as x , y , z, we obtain for the case A ;~ B ;~ C ;~ A the integral-invariant manifold in the (r/(-space x

y -

BC(B - C)

z -

CA(C - A)

-

AB(A

,~(t).

(6.14)

- B)

Physically, it means that ~, q, ( satisfying (6.14) do not affect the energy nor the angular momentum of the body. From conservation property of integrals at the left-hand side of (6.12), (6.13), it follows that perturbed trajectories either lie entirely on the manifold (6.14) or do not intersect it at all. If for nominal motions p ( t ) , q(t), r(t) there are no perturbed trajectories that lie on the manifold (6.14), then those motions are stable by Lyapunov's theorem on stability [1, Section 16] with the function V of (6.11) which is positive definite if (6.14) does not contain perturbed trajectories. On the other hand, if, for certain nominal motions, the manifold (6.14) contains perturbed trajectories which deviate from the region (3.2), then those nominal motions are unstable.

816

E . A . G A L P E R I N and V. V. R U M Y A N T S E V

N o m i n a l solution p ( t ) , q(t), r(t) will have a perturbed trajectory lying on the manifold, if initial values Po, qo, r0, and ~o, r/o, Co are such that (6.14) is satisfied. Then the perturbed m o t i o n will belong to the m a n i f o l d for all t _> to which can be verified directly as follows. Taking derivative by virtue of (6.1), (6.4), we have

dx dt

1 d2 A - ABC . . . . 2 dt 2 ( B - C)

A C

B

= qrg + (p + O ( m + q( + ~ 0 .

(6.15)

Replacing in (6.15) x by y = 2 q r / + /72, then by z = 2r( + (2, and taking derivatives by virtue o f (6.2), (6.5), then (6.3), (6.6), we get exactly the same expressions as (6.15), proving that (6.14) holds for all t _> t 0. Let us study the existence and b e h a v i o r o f perturbed trajectories lying on the manifold. Assume that A > B > C, in which case (6.14) contains perturbed trajectories if and only if x > 0, orx=y=z=2

y < 0,

z > 0,

=0atisolatedmomentst x < 0,

i>t

y > 0,

o,j=

z < 0,

thus, ~.(t) > 0,

(6.16)

1,2 . . . . . or so that, ~.(t) < 0.

(6.17)

Case 3.1. Let Po = qo = 0, ro # 0. Then, according to (6.1)-(6.3) the unique solution is p ° ( t ) = q°(t) - O, r°(t) = ro = const. For this solution, we have x = ~2 > 0, y = r/2 > 0, contradicting (6.16) and (6.17). This means that there are no perturbed trajectories on the m a n i f o l d (6.14). If p0 # 0, qo = ro = 0, then p ° ( t ) = Po = const, q°(t) = r°(t) = O, and we have y = r/2 > 0, z = (2 > 0, again with no trajectories on (6.14) due to contradictions with (6.16) and (6.17). W e conclude that constant rotations a r o u n d the greatest (r°(t) = ro = const) or smallest (P°(t) = Po -- const) axes of the ellipsoid of inertia are stable. This case was p r o v e n stable in [2, Section 15] by the construction o f a special integral, as given in what follows. Case 3.2. Let Po = ro h a v e x = ~2 > 0, z = is satisfied and (6.14) In this case, ;t(t) >

= 0, qo # 0. In this case, p ° ( t ) = r°(t) = 0, q°(t) = qo = const, and we (2 > 0, y = 2qor/ + r/2 < 0 for small r/such that q o r / < 0, so that (6.16) contains perturbed trajectories. 0 and f r o m (6.15) we have 1 d3~ = A B C - ; 7 = ~ ( q o + r/) > 0 z

ut

(6.18)

for small t/, if ~(qo > O. It m e a n s that I(], It/], ](] will grow in the sectors where ~(qo > O, qoV/< 0 and n o m i n a l solution p = r - O, q = qo = const is repulsive in those sectors for perturbed trajectories lying on (6.14).

Stability analysis

817

To prove instability in this case, let, for example, qo > 0, ~o > 0, r/o < 0, Go > 0 and take 0 < e < qo. Then, if I~/ol -< e, we have from (6.18)

ABC~t

= 2~((q o + r/) _> 2(o(o(q o - e) > 0

(6.19)

uniformly for all ~ _> ~0, ( -> (0, 0 > r/_> - e . This means that )~(t) > 0 will grow for t _> to until such t = t* for which r/(t*) = - q o . If initial perturbations ~2 + r/2 + ~2 = 62, then, however small 6 z may be, ~(t), ((t) are increasing and ~/(t) < 0 is decreasing on the interval [to, t*), so that at t = t* there will be ~2(t*) + r/2(t *) + (2(t*) > q2, implying instability. We conclude that constant rotation around the middle (q°(t) = qo = const) axis of the ellipsoid of inertia is unstable. This case was proven unstable in [2, Section 15] by the construction of a special function, as in what follows.

Case 3.3. Assume that two o f p o , qo, ro are nonzero. From (6.1)-(6.3) with different A, B, C, it follows that at most one o f p ( t ) , q(t), r(t) can vanish and only at an isolated moment t _> to. Indeed, if two of p, q, r are zero at some moment, then we return to cases 3.1, 3.2; and if one, say, p(t) - 0 on some interval (t, t + 6), 6 > 0, t >_ t 0, then on that interval also d p / d t =- O, thus, either q = 0 or r = 0, due to (6.1), and we again return to cases 3.1, 3.2. Following Appell [16, Sections 387-392, pp. 162-186], let Tl

T2 fi3, 7, =

so that T~ = ~ p 2 , T2 = ~D2p2

(6.20)

and suppose first that if) is different from A, B, C. From (6.7) written as A p 2 + Bq 2 + C r 2 = fi)/t 2 = const, we see that p is a projection of the vector to of angular velocity and :D has the dimension of the moment of inertia, so that the constancy of ~Dp2 represents the energy conservation law. Denote, as in [16] f2

2 ~D(~D - C )

=/~ ~ - C )

g2

>0,

=/J

2 ~)(A - ~D)

~-A

B) > 0 ,

(6.21)

where A > fi3 > C follows from the integrals (6.7), (6.8) by eliminating p or r which yields also the relations

p2 = B ( B A(A r2

_ B(A

C)

-

C) (f2 _ q2),

(6.22)

(g2 _ q2).

(6.23)

B)

C(A - C) In order that p, r be real, it should be

qE(t ) ___min(f2, g2).

(6.24)

818

E.A. GALPERIN and V. V. RLIMYANTSEV

By s y m m e t r y , we can assume that g2 > f 2 which is the case if B > ~ , see (6.21). In this case, r(t) o f (6.23) does not vanish and conserves the sign o f r o = r(to). Substituting (6.22), (6.23) with p r o p e r signs o f roots into (6.2) yields the solution ff))(B -A--ffC

r = n(t - to) = U~ ( A

q = f s n r,

I

C)

(t - to),

(6.25)

where to represents the m o m e n t when q(to) = 0 while increasing which corresponds to f > 0 in (6.25). W i t h o u t loss o f generality, assume that ro > 0. Then, while q(t) is increasing f r o m zero to f > 0, see (6.24), the value r(t) is decreasing, see (6.23), f r o m

[B(A B) -

m a x r = r ( t o ) = ro = g-~J C--~

g > O,

C)'

(6.26)

to

/s(A

m i n r = p = ~g2 - j ~ - ~

-

n)

C) > 0.

(6.27)

Let t 1 be the first m o m e n t when q = f . Then on (to, tl) with q(t) increasing, we have d q / d t > 0 and, with r(t) > 0, it follows f r o m (6.2) that p ( t ) < 0 is increasing f r o m B(B-

p(to) = Po = - f . ~ A ( A

C)

C) - m i n p < 0

(6.28)

to p(tl) = 0. This yields d p / d t > 0, d r / d t < 0 on (to, tl) in compliance with (6.1)-(6.3). On ( q , t2) the value q(t) is decreasing f r o m q(tl) = f to q(t2) = O, d q / d t < 0, and, with r(t) > 0, we have p ( t ) > 0 and d p / d t > 0, the value p ( t ) increasing f r o m p(tl) = 0 to p(t2) = m a x p = - P o > 0. The value r(t) is increasing on ( q , tz) f r o m r(tl) = p = m i n r > 0 to r(tz) = m a x r = r o, thus dr~dr > 0, in compliance with (6.1)-(6.3). The reader can continue this sign investigation for intervals (tz, t3) and (t 3 , t4). T h e period is given in [16, p. 167] T=

t4 - t o = -

r/

0 X/(1 - S2)(1 - k2s 2)

f2 with n f r o m (6.25) and k 2 - g2"

(6.29)

After this period, at t = t4, we have p(t4) = Po < 0, q(t4) = q(to) = qo = 0, r(t4) = ro > 0, and the process is repeated. If r o < 0, then certain signs a b o v e are reversed but p , q, r are passing t h r o u g h the same extreme values defined by (6.21)-(6.25). Consider a perturbed trajectory lying on the m a n i f o l d (6.14). Since qo = 0, this will be the case if ~o, r/o, Co are such that

2polo + B C ( B - C)

2ro;o + CA(C - A)

A B ( A - B)

- ;~(to) < o.

(6.30)

Let ~g + r/2 + (2 = ~2 with ~ > 0 arbitrarily small, and take certain (o > 0, go > 0, (o < 0 in agreement with ).(to) < 0, (6.30), since Po < 0, r o > 0.

Stability analysis

819

W h e n t > to increases, we have according to (6.14) 2p~ + ~2 BC(B -

C)

2r( + (2

2qt/ + t/2 -

CA(C

- A)

-

AB(A

- B)

=

2(t).

(6.31)

At t = tl we have P(tl) = O, q ( t l ) = f > O, r(t 0 = p > O, so perturbations ~l, r/l, (1 satisfy the relations ~12

2ft/l + r/2

BC(B - C)

CA(C

- A)

2p~, + (z AB(A

- ).(tl) > 0,

- B)

(6.32)

since ~2 > 0, B - C > 0. Note that ~, r/, ( cannot vanish together (this being the trivial solution o f the generalized p e r t u r b a t i o n equations), so if it were ~ = 0, then either t/1 = - 2 f or (1 = - 2 p and instability would follow. The function ).(t) is continuous, so there is t, e (to, tl) such that 2 ( t . ) = 0. D e n o m i n a t o r s o f (6.31) being nonzero, the perturbations ~., r/., ( . at t = t. satisfy the relation ~ , ( 2 p , + ~,) = r/,(2q, + r/,) = ( , ( 2 r , + ( , ) = 0.

(6.33)

Since ~,, t/,, ( , cannot vanish together, one of them, in magnitude, attains double the m a g n i t u d e o f the corresponding nominal value for any small ~, i.e. for arbitrarily small initial perturbations. This means instability o f all nominal solutions in this case. • We have assumed above that B > 33 > C and r o > 0. Cases A > 33 > B a n d / o r r o < 0 present the same situation o f instability with modified sign conditions where a p p r o p r i a t e . Different particular cases are obtained if 33 equals one of m o m e n t s of inertia: A, B or C [16, Section 389, pp. 169-171]. If if) = C, then, eliminating r z f r o m (6.7)-(6.8) and using notations (6.20), we have ApZ(A

- C ) + B q 2 ( B - C) =/2233(33 - C) = 0,

(6.34)

so t h a t p = q = 0 (this follows f r o m (6.21), (6.22) in the limit if 33 - C ~ 0, thus, f 2 ~ 0). If 33 = A, then, similar to (6.34), by eliminating p2 f r o m (6.7)-(6.8), or in the limit as A - 33 --, 0, g2 ~ 0, in (6.21), (6.23), we obtain q = r - 0. Both cases comprise our case 3.1 a b o v e with stable nominal motions r = r o = / 2 = const, or p = Po = / 2 = const. • S i n g u l a r c a s e 3.3b. If 33 = B, then, eliminating

q2 f r o m

(6.7)-(6.8) and using notations (6.20),

we have Ap2(A

- B ) - C r 2 ( B - C) =/225)(33 - B) = 0.

(6.35)

In this case P0 ;~ 0 and ro ~ 0 (otherwise, it would be case 3.2 above). F r o m (6.35) it follows /

p ( t ) _ P o _ + [ C ( B - C ) _ b, r(t) ro - "~ A ( A B)

t >_ t o .

(6.36)

With 33 = B, we h a v e f 2 = g2 =/2z, see (6.21), and f r o m (6.22), (6.23) it follows that q2
ds dt

-

n(1 - s 2) s i g n b = +n(1 - s2), -

n =/2/(A "~1

- B)(B - C) A C '

(6.37)

820

E.A. GALPERIN and V. V. RUMYANTSEV

where n is the same as in (6.25) with ~ = B. Taking the sign plus in (6.37) setting r = n(t - to) and integrating, we obtain 1 - s _ 1 - s o e_2Z' 1 +s o

z > 0,

s(0) = So.

(6.38)

1 +s

As r ~ ao, we have s ~ 1, thus q ~ / ~ and the energy integral (6.7) with T1 = n/./2 implies p ~ 0, r - 0. I f we t o o k the sign minus in (6.37), we would have s ~ - 1 , q ~ -/~, p ~ 0, r ~ 0 as r ---, oo. We see that the vector o f angular velocity (p, q, r) tends as t ---, oo to a limiting position (0, +~t, 0) and convergence is m o n o t o n i c , if qo = ~tSo -> 0 with the choice o f the sign + in (6.37), or qo =/tSo ~ 0 with the choice o f sign minus in (6.37), see f o r m u l a e (24) in [16, p. 171] where So = 0, in which cases p , q, r do not change sign. Let us demonstrate that m o t i o n s in this case are unstable. Consider perturbations Co2 + v/2 + (2 = $2 with $ > 0 arbitrarily small and such that ~o ~ b(o. Then perturbed m o t i o n s (p*(t), q*(t), r*(t)) = (p + C, q + vl, r + (), t >_ to, belong to the regular case 3.3 since P~

Po + Co ~ b, :

thus T2*

II)*

T2

(6.39)

r0 +

This means that, according to the case 3.3 above, the variablep*(t) periodically changes sign. The nominal variable p ( t ) does not change sign, hence, C(t) is o f the magnitude o f p ( t ) at some m o m e n t s t = t i, i = 1 , 2 . . . . . and instability follows. It is easy to conclude instability in terms o f some fixed deviation e > 0 attained at some t = t, whatever small initial perturbations ~o,r/o,(o m a y be. Indeed, the variable q*(t) = q(t) + rl(t) also periodically changes sign and, since the nominal variable q(t) -~ +_~, it follows that at some m o m e n t t -- t . < t o we shall have Ir/(t,)l = p / 2 = e whatever small Co, r/o, (o m a y be, implying instability. We conclude that all motions in the case 3.3 (including ~ = B) are unstable.

Case 4. A ~ B - - C . In this case equations (6.1), (6.4) yield p ( t ) = P o - - c o n s t , and C(t) -- Co -- const so that systems (6.1)-(6.3) and (6.4)-(6.6) are reduced to the following dp O, d--t =

B dq - ~ + (A - B)Por = O,

dr B - ~ + (B - A ) p o q = 0, = O,

Bil = (B - A ) ( p o ( + Cor + CoO,

B~ = (A - B)(Porl + Coq + CoVl)•

(6.40)

(6.41) (6.42) (6.43)

N o w the problem o f stability can be solved by consideration o f the simple positive definite function V = ~2 + r/2 + (2 (6.44)

Stability analysis

821

with its total derivative BF" = B -~d V = 2(B - A)Co(rrl - q~)

(6.45)

on trajectories of (6.42)-(6.43). C a s e 4.1. I f qo = ro = 0, then q ( t ) = r(t) ==- O. In this case, I / = 0 and we conclude that constant rotations a r o u n d the axis of revolution are stable; this is essentially the same as in case 3.1. C a s e 4.2. C o n s t a n t rotation a r o u n d an axis in the equatorial plane o f the ellipsoid o f revolution, B = C, occurs i f p ( t ) = q ( t ) -- O, r(t) = r0 = const where we identify coordinate r with the axis o f rotation. In this case, we have BI? = 2(B - A)Cororl > 0 in the region where r/ has the same sign as ( B - A ) C o r o. This suggests possible instability, though lk is not positive definite in the region V > 0 and so V o f (6.44) does not qualify for use in C h e t a e v ' s instability t h e o r e m , see [2, Section 13]; here in (3.9) we have gt = H N IIIwll < ,t, ,r > 0}, (3.10) is satisfied but (3.11) is not. Let us verify directly that the m o t i o n is, indeed, unstable. In this case, equations (6.42)-(6.43) b e c o m e = O,

(6.46)

Bil = ( B - A)Co(r o + ()

B~ = (A - B)~or/.

(6.47)

Let B - A > 0, r0 > 0, r/0 = 0, C0 = (0 = 8 > 0, 8 small. T h e n V(to) = 282 , and f r o m (6.46), (6.47) we have ~ / > 0, ~ < 0 on (to, t.) with r/(t) increasing, ((t) decreasing until r 0 + ( ( t . ) = 0, that is, ( ( t , ) = - r o , yielding V ( t . ) = C2 + tl 2 + (2, > 82 + r g > r 2 irrespective o f h o w small 8 > 0 m a y be. This implies instability of constant rotations a r o u n d an equatorial axis o f the ellipsoid o f revolution. Note that i f p ( t ) -- 0, q ( t ) = qo, r(t) = ro, then (6.47) would be replaced by B ~ = ( A - B)Co(qo + r/); taking q0 > 0, one gets the same p r o o f o f instability. This is, obviously, the same case. C a s e 4.3. F o r a n o m i n a l m o t i o n p = P0 = const, q(t), r(t), called regular precession, one can use d e c o m p o s e d structure o f the system (6.42)-(6.43) and apply the idea o f L y a p u n o v a b o u t conditional stability for restricted perturbations, see [1, Section 1] or [2, Section 3], to obtain unconditional stability.

F r o m the equality C = C0, we conclude that all motions in case 4 are stable with respect to p , whatever its value p = P0. U n d e r perturbations restricted by C0 = 0, the m o t i o n (Po, q, r) is stable with respect to q, r since f r o m (6.45) we have 1? = 0, so that V(0, r/, () = r/2 + (2 = r/2 + (2 = const for all t ___ to. This means that with C0 ;~ 0 it is possible to reduce this case to the previous one corresponding to the values/50 = Po + Go, (o = 0, since the integrals p = Po, q2 + r 2 = q2 + r g are continuous functions. Thus, all regular precessions in the case o f Euler are stable with respect to p , q, r. This conclusion follows also f r o m the explicit solution o f (6.40)-(6.41) P = P0,

q = q0 cos o~t - r o sin o~t,

since Icoso~t[ _ 1 and Isinoal -< 1.

•

r = q0 sin oct + r o cos a t ,

~ -

A-B B

Po,

822

E.A. GALPERIN and V. V. RUMYANTSEV

For comparison, we briefly reproduce here the analysis given in [2, Section 15], for cases 3.1, 3.2, 4.1. Let Po = qo = 0, ro > 0, which defines the unique solution p = q - 0, r(t) = ro = const. In this case, perturbation equations become A ~ = (B - C)vl(r o + ~),

(6.48)

Bil = (C - A)~(ro + (),

(6.49)

C ( = (A - B)~tl.

(6.50)

If A _ > B > C or A _ < B < C , stability follows from the positive definite integral of (6.48)-(6.50) proposed in [2, Section 15]

V-

A-

___~2

_]_

/72 -+ [A~ 2 + Br/2 + 2Cro(

+

C(2] 2,

(6.51)

where we take sign plus, i f A _> B > C, and sign minus, i f A _ B < C. This establishes stability of constant rotations around the greatest or the smallest axes of ellipsoid of inertia. Instability of rotation around the middle axis (A < C < B) follows by consideration of the function [2, Section 15], F1 = ~v/

(6.52)

with its total derivative ~=(ro+

() ( B - - ~

v/2 + ~------~A C ,2)

(6.53)

on trajectories of perturbed motions (6.48)-(6.50). With fixed ro > 0 and small (, we have > 0, and instability follows by Chetaev's theorem [2] in the sectors ~ r / > 0, cf. (3.9)-(3.11). In case A > C > B, we have ~ < 0 with instability in the sectors ~ r / < 0. • S u m m a r y . The rest p = q = r = 0 and all motions in trivial case A = B = C are stable. The

motion p = q - 0, r(t) = r 0 = const in cases A _< B < C or A ___B > C (i.e. constant rotation around extreme axis C, including circular ellipsoids of inertia) is also stable. From the above analysis, we see that all other motions in the case A ~ B ~ C ~ A are unstable. In the case of circular ellipsoid of inertia (ellipsoid of revolution, A ~ B = C), constant rotation around an equatorial axis is unstable and all other motions are stable.

7. CONCLUSIONS On the basis of the generalized perturbation equation proposed in this paper, stability analysis can be made for whole sets of trajectories without the knowledge of explicit particular solutions. A notion of consistency is introduced to allow for the use of discontinuous V-functions with component functions acting over bounded time-space regions in cases where it is difficult to construct a classical continuous Lyapunov function satisfying corresponding Lyapunov conditions on the whole half-axis [to, oo). With this methodology, computational methods are proposed to estimate the measure of contraction or repulsion over bounded time-space domains and for the analysis of behavior of perturbed motions over given bounded time-space regions. It is demonstrated how to use bundles of first integrals for stability analysis

Stability analysis

823

of sets of solutions. As an application, stability properties are investigated for all motions in the Euler case of the motion of a rigid body around its fixed center of mass without external forces. REFERENCES 1. LYAPUNOV A. M., Probl~me g~n~ral de la stabilit6 du mouvement, Ann. Fac. Sci. Univ. Toulouse, 2 s6rie IX, 203-474 (1908). Photographic reproduction in Annals of Mathematical Studies, No. 17. Princeton University Press, Princeton, N. J. (1947). (First published in Russian by the Kharkov Mathematical Society in 1892.) 2. CHETAEV N. G., Stability of Motion. Pergamon Press, New York (1961). (First published in Russian by Gostekhizdat in 1946.) 3. MALKIN I. G., Theorie der Stabilitat einer Bewegung. Oldenbourg, Munich (1959). English translation by Atomic Energy Commission, Washington, D. C. (First published in Russian by Gostekhizdat in 1952.) 4. KRASOVSKII N. N., Stability of Motion. Stanford University Press, Stanford, CA (1963). (First published in Russian, Moscow, 1959.) 5. RUMYANTSEV V. V. & OZIRANER A. S. Stability and Stabilization of Motion with Respect to a Part of Variables. Nauka, Moscow (1987). ( In Russian.) 6. BELLMAN R., Stability Theory of Differential Equations. Academic Press, New York 0953). 7. CESARI L., Asymptotic Behavior and Stability Problems in Ordinary Differential Equations. Springer, Berlin (1959). 8. LASALLE J. P. & LEFSCHETZ S., Stability by Lyapunov's Direct Method with Applications. Academic Press, New York (1961). 9. LEFSCHETZ S., Stability of Nonlinear Control Systems. Academic Press, New York (1965). 10. LEITMANN G., On stabilizing a linear system with bounded state uncertainty, in Topics in Contemporary Mechanics, CISM monograph No. 20. Springer, Vienna (1974). ll. MATROSOV V. M., On the stability of motion, J. Appl. Math. Mech. 26, 885-895 (1962) (English translation, pp. 1337-1353). 12. MASSERA J. L. & SCHAFFER J. J., Linear Differential Equations and Function Spaces. Academic Press, New York (1966). 13. LAKSHMIKANTHAM V., LEELA S. & MARTYNYUK A. A., Stability Analysis of Nonlinear Systems. Dekker, New York (1989). 14. GALPERIN E. A. & SKOWRONSKI J. M., Geometry of V-functions and the Liapunov stability theory, Nonlinear Analysis 11, 183-197 (1987) 15. GALPERIN E. A. & ZHENG Q., Global Solutions in Optimal Control and Games. NP Research Publ., Montreal (1991). 16. APPELL P., Traitd de Mdcanique Rationnelle, I~dition 3, Tome 2. Gauthier-Villars, Paris (1911).