A Sum-of-Squares approach to the Stability and Control of Interconnected Systems using Vector Lyapunov Functions

0.22 MB·English

by  Soumya Kundu

A Sum-of-Squares approach to the Stability and Control of
Interconnected Systems using Vector Lyapunov Functions*
Soumya Kundu1 and Marian Anghel2
Abstract—Stabilityanalysistoolsareessentialtounderstand- Primarily we will focus on an example of a randomly
ing and controlling any engineering system. Recently, sum- generatednetworkofninemodified1 VanderPoloscillators.
5 of-squares (SOS) based methods have been used to compute
Each Van der Pol oscillator can be represented as a two-
1 Lyapunov based estimates for the region-of-attraction (ROA)
state system with state dynamic equations as polynomials
0 ofpolynomialdynamicalsystems. Butforareal-lifelargescale
2 dynamicalsystemthismethodbecomesinapplicablebecauseof of degree three [12]. The network is then decomposed into
growing computational burden. In such a case, it is important many interacting subsystems. Each subsystem parameters
n
todevelopasubsystembasedstabilityanalysisapproachwhich are so chosen that individually each subsystem is stable,
a
isthefocusoftheworkpresentedhere.Aparallelandscalable
J when the disturbances from neighbors are zero. Sum-of-
algorithmisusedtoinferstabilityofaninterconnectedsystem,
2 squares based expanding interior algorithm [13], [14] is
with the help of the subsystem Lyapunov functions. Locally
2 computablecontrol lawsareproposedtoguaranteeasymptotic used to obtain estimate of region of attraction as sub-
stability under a given disturbance. level sets of polynomial Lyapunov functions for each such
]
S subsystem. Finally a sum-of-squares based scalable and
D I. INTRODUCTION parallel algorithm is used to certify stability in the sense
of Lyapunov of the interconnected system by using the
.
h Toward the end of the nineteenth century the Russian subsystem Lyapunov functions computed in the previous
t
a mathematicianA.M.Lyapunov[1]hasintroducedanumber step. A distributed control strategy is proposed that can
m of powerful tools for the stability analysis of nonlinear guarantee asymptotic stability of the interconnected system
[ dynamical systems. His stability results have been later under given disturbances. Following some brief background

generalized by Barbashin, Krasovskii, and LaSalle, while in Sec. II we outline the problem statement in Sec. III.
2
v control system engineers have used Lyapunov’s methods to An algorithmic approach to certifying asymptotic stability
6 design stabilizing feedback controllers — see [2] and refer- is presented in Sec. IV while a distributed control strategy
6 encestherein.Nevertheless,constructinga systemLyapunov is discussed in Sec. V. Sec. VI shows an application of our
2 function is usually a difficult task, which becomes daunting stability analysis and control approach to a network of Van
5
when the size of the dynamical system increases. der Pol oscillators. We conclude the article in Sec. VII.
0
. A more practical approach is to define the Lyapunov
1 II. BASIC CONCEPTS AND BACKGROUND
function of the interconnected system as some function of
0
Beforeformulatingtheproblem,letusbrieflyreviewsome
5 the subsystemLyapunovfunctions.Thereare differentfunc-
ofthekeyconceptsbehindouranalysis.Wewillfirstdiscuss
1 tionalformsfortheLyapunovfunctionoftheinterconnected
: system, such as a scalar Lyapunov function expressed as how the stability of a dynamicalsystem can be analyzed by
v
constructing suitable Lyapunov functions. Then we briefly
i a weighted sum of the subsystem Lyapunov functions, or
X refer to sum-of-square polynomials and a very useful result
applications of vector Lyapunov functions and comparison
r principles[3]–[6].FormulationsusingvectorLyapunovfunc- which helps us in formulatingthe sum-of-squaresproblems.
a
tions [7], [8] are computationally attractive because of their A. Lyapunov Stability Methods
parallel structure and scalability. In [4], it was shown that
Let us consider the dynamical system
if the subsystem Lyapunov functions and the interactions
satisfy certain conditions, then application of comparison x˙(t)= f(x(t)), t≥0, x∈Rn, f(0 )=0 (1)
n×1 n×1
equations [9]–[11] can provide a certificate of exponential
We assume that x = 0 (for simplicity, t will be often
stability ofthe interconnectedsystems. In thisworkweseek n×1
dropped when obvious) is an equilibrium point of the dy-
an algorithmic certification of asymptotic stability via the
namical system2, and f :Rn→Rn is locally Lipschitz. An
vector Lyapunov function approach, where each subsystem
important notion of stability is as follows:
Lyapunovfunctionsareexpressedinsomepolynomialform.
Definition 1 The equilibrium point at origin is called
*This work was supported by the U.S. Department of Energy through
asymptoticallystable if it is stable in the sense of Lyapunov,
theLANL/LDRDProgram.
1SoumyaKunduiswiththeCenterforNonlinearStudiesandInformation
Sciences Group (CCS-3), Los Alamos National Laboratory, Los Alamos, 1WechoosetheoscillatorparametersinsuchawaythattheVanderPol
[email protected] oscillators haveastable equilibrium atorigin.
2Marian Anghel is with the Information Sciences Group (CCS-3), Los 2Note that this is not arestrictive assumption, since byshifting ofstate
AlamosNationalLaboratory,LosAlamos,[email protected] variables, theorigincanalways bemadeanequilibrium point.
and if toolbox SOSTOOLS [18], [19] along with a semidefinite
programming solver such as SeDuMi [20].
∃d˜ >0 s.t. kx(0)k <d˜ =⇒ lim kx(t)k =0. (2)
2 t→+¥ 2 SOStechniquecanbeusedtofindapolynomialLyapunov
function V(x):Rn →R, with V(0 )=0, which satisfies
The Lyapunov stability theorem [1], [15] presents a suf- n×1
the following SOS conditions [13], [18], [19], [21]–[24],
ficient condition of stability through the construction of a
certain positive definite function. V(x)−f (x)∈S ,∀x∈D (7a)
1 n
−V˙(x)−f (x)∈S ,∀x∈D (7b)
2 n
Theorem 1 Theequilbriumpointx=0 ofthedynamical
system in (1) is stable in the sense of Lny×a1punov in D ∈Rn, for some domain D around the origin and positive definite
functionsf (x), f (x).OftenitisconvenienttochooseD:=
if there exists a continuously differentiable positive definite 1 2
function V˜ : D → R (henceforth referred to as Lyapunov {x∈Rn|p(x)<b }, for some positive definite function p(x)
and b >0.
function) such that,
An important result from algebraic geometry called Puti-
V˜(0)=0 (3a) nar’sPositivstellensatztheorem[25],[26]helpsintranslating
V˜(x)>0,∀x∈D\{0 } (3b) the SOS conditions into SOS feasibility problems. Before
n×1
stating the theorem, let us define:
and, −V˙˜(x)≥0,∀x∈D (3c)
Further, if V˜ satisfies −V˙˜(x)>0,∀x∈D\{0 }, then the Definition 3 Given gj ∈ Rn, for j = 1,2,...,m, the
equilibrium point at origin is asymptotically snt×a1ble in D. quadratic module generated by gj’s is M(g1,g2,...,gm):=
s +(cid:229) s g s ,s ∈S ,∀j
0 j j j 0 j n
When there exists such a function V˜(x), the region of
(cid:8)Then the Putin(cid:12)ar’s Positivestel(cid:9)lensatz theorem states,
attraction (ROA) of the stable equilibrium point at origin (cid:12)
can be (conservatively) estimated as
Theorem 2 Let K ={x∈Rn|g (x)≥0,...,g (x)≥0} be
1 m
RA:= x∈D V˜(x)≤g max (4a) a compact set. Suppose there exists u(x)∈Rn such that
where, g max:=(cid:8)argmgax(cid:12) x∈Rn V˜(x)(cid:9)≤g ⊆D (4b) u(x)∈M(g1,g2,...,gm), (8a)
(cid:12)
(cid:8) (cid:12) (cid:9) and, {x∈Rn|u(x)≥0} is compact. (8b)
It can be noted that, without any(cid:12)loss of generality, the
Lyapunov function can be scaled by g max, so that the ROA If p(x) is positive on K , then p(x)∈M(g ,g ,...,g ).
1 2 m
is given by,
Note: often in this work, for the g’s used, the constraints
i
R :={x∈Rn|V(x)≤1} (5a) (8) would be redundant, i.e. the existence of u(x) would be
A
where,V(x)=V˜(x)/g max (5b) guaranteed [26]3.
III. PROBLEM OUTLINE
Henceforth, for simplicity, we would assume, without any
serious loss of generality, that the ROA is estimated to be Let us assume that the dynamical system in (1) is in
sub-level set of V(x)=1. polynomialform,i.e. f is a vectorof n polynomials4. Given
the full dynamical system (1), we can decompose it into
B. Sum-of-Squaresand Putinar’s Positivestellensatz
subsystems
WhiletheTheorem1givesasufficientconditionforstabil-
∀i=1,2,...,m,
ity,itisstillnotatrivialtasktofindasuitablefunctionV(x)
thatsatisfiestheconditionsofstability,evenwhentheorigin x˙i= fi(xi)+gi(x), xi∈Rni, (9a)
is a stable equilbrium point. Relatively recent studies have f(0 )=0 , g(0 )=0 (9b)
i ni×1 ni×1 i n×1 ni×1
exploredhow sum-of-squaresbased optimization techniques x= xT,xT,...,xT T ∈Rn (9c)
can be utilized in finding Lyapunov functions by restricting 1 2 m
m
the search space to sum-of-square polynomials [13], [14], n=(cid:0)(cid:229) n , x ∩x =(cid:1) 0/ (9d)
i i j
[16], [17]. Let us denote R as the set of all polynomialsin
n i=1
x∈Rn. Then,
In the decomposed system description, the f’s denote the
i
isolated subsystem dynamics, and g’s are the interactions
i
Definition 2 Amultivariatepolynomial p(x)∈Rn isasum- from the neighbors [14]. Let us also denote by
of-squares (SOS) if there exist some polynomial functions
h(x),i=1...r such that p(x)=(cid:229) r h2(x), and the set of Ni:={i}∪ j ∃xi,xj, s.t. gij(xi,xj)6=0 , ∀i, (10)
i i=1 i
all such SOS polynomials is denoted by
3Simplicityinfor(cid:8)mul(cid:12)(cid:12)atinganSOSproblemmotivatedu(cid:9)stousePutinar’s
S :={p(x)∈R | p is SOS} (6) version of Positivstellensatz over other, more general, versions of the
n n Positivstellensatz theorem [26].
Givena polynomial p∈R , checkingif itis SOS is a semi- 4Ifthe dynamics is not inpolynomial form, ithas toberecasted into a
n polynomialform,withpossibleadditions ofequalityconstraints [14],[21],
definite problem which can be solved with a MATLAB(cid:13)R [23],[27],acasenotconsidered inthis work.
the set of neighbors (including the subsystem itself). We could be expressed as some sub-level set of that Lyapunov
assumethattheisolatedsubsystemsareindividually(locally) function. While it is very hard to obtain a scalar Lyapunov
stable,andthereexistLyapunovfunctionsforeachoftheiso- function for the full interconnected system, one could use
latedsubsystems.Thegoalistodevelopaframeworkforthe vector Lyapunovfunctionapproachto obtain certification of
stability analysis of the full interconnected system by using stability in a scalable way.
the localsubsystem Lyapunovfunctionsand consideringthe In thispresentwork,we choosenotto imposeanyfurther
neighbor interactions. restriction on the Lyapunov functions V(x) than requiring
i i
Given a decomposition (9), the next step is to find poly- that those are in polynomial forms, and concern ourselves
nomial Lyapunov functionsVi(xi):Rni →R, such that [13], with asymptotic stability. In general, it is difficult to test
a necessary and sufficient condition of asymptotic stability,
such as the one given in (14). Nevertheless, it is possible
∀x ∈D, V(x)−f (x)∈S , (11a)
i i i i i1 i ni to derive sufficient conditions of asymptotic stability under
and, −(cid:209) Vi(xi)Tfi(xi)−f i2(xi)∈S ni (11b) certain scenarios. Let us now present a distributed iterative
where, D :={x ∈Rni|p(x)≤b } procedure which can be used to certify asymptotic stability
i i i i i
in a domain defined by sub-level sets of the subsystem
for some positive definite functions f (x), f (x), p(x)
i1 i i2 i i i Lyapunov functions.
and positive scalars b . Starting from an initial Lyapunov
i
function candidate obtained using (11) and a corresponding A. Algorithmic Test of Asymptotic Stabiltiy
estimateoftheregionofattraction,aniterativeprocesscalled Before proceedingto explainingouralgorithm,let us first
expandinginterioralgorithm,[13],[14],isusedtoiteratively note the following result:
enlarge the estimate of the region of attraction by finding a
better Lyapunov function at each step of the algorithm. At Lemma 1 Suppose, for all i ∈ {1,2,...,m}, there exists
the completion of this iterative step, the stability of each a strictly monotonically decreasing sequence of scalars
isolated subsystem (assuming no interaction) is quantified e k ,k∈{0,1,2,...}, such that
i
by its Lyapunov function V(x), with an estimation of the
boundary of the domain oif iattraction given by R = (cid:8)∀i,(cid:9)k, V˙i(x):=(cid:209) Vi(xi)T(fi(xi)+gi(x))<0, ∀x∈Dik (14a)
A,i
{xi∈Rni|Vi(xi)≤1}. where, Dk:= x∈Rn eik+1≤Vi(xi)≤eik, (14b)
The Lyapunov level-sets can be used to express the i V (x )≤e k ∀j6=i
(cid:26) (cid:12) j j j (cid:27)
strengthofadisturbance.Theequilibriumpointofthesystem (cid:12)
Then the system (1) is a(cid:12)symptotically stable in the domain
iastaordigisintucrboarrnecsepofrnodmsttohitsheeqleuvileilbsreiutmVi(p0oniin×t1,)th=e0s,t∀atie.sIfotfhethree x∈Rn mi=1Vi(xi)≤ei0(cid:12) , if limk→+¥ eik=0,∀i5.
systemwouldmovetosomepointx(0)awayfromtheorigin. (cid:8) Pro(cid:12)oTf: Please refer(cid:9)to Appendix A.
(cid:12)
Thisdisturbedinitialconditionwouldresultinpositivelevel- Using the Lemma 1 we can devise a simple iterative SOS
setsVi(xi(0))=gi0∈(0,1]forsomeorallofthesubsystems. algorithm to certify whether or not a domain D defined by
A necessary and sufficient condition of asymptotic stability
m
can then be translated into the condition D := x∈Rn V(x)≤v , (15)
i i oi
∀i, Vi(xi(0))=gi0 =⇒ ∀i, t→lim+¥ Vi(xi(t))=0 (12) for some scalars vo(i ∈(0,1](cid:12)(cid:12)(cid:12)(cid:12),i\=∀1i, is a regio)n of asymptotic
In the rest of the article, we present SOS algorithms to stability for the system in (1(cid:12)). It is to be noted that, using
test stability conditions and design local (subsystem-level) the Putinar’s Positivstellensatz theorem (Theorem 2), the
control laws to achieve asymptotic stability. condition in (14) essentially translates into equivalent SOS
feasibility conditions
IV. STABILITY UNDER INTERACTIONS
It is assumed that the isolated subsystems in (9) are all ∀i,k, ∃s ik0,s ikj∈S n¯i, s.t. (16)
(locally) asymptotically stable, and there exist subsystem −(cid:209) VT(f +g)−s k V −e k+1 − (cid:229) s k e k−V ∈S
LyapunovfunctionsV(x),∀i.Theestimatedregionofattrac- i i i i0 i i ij j j n¯i
tionoftheinterconnecitedi systemundernointeraction,R0,is (cid:229) (cid:16) (cid:17) j∈Ni (cid:16) (cid:17)
A where, n¯ = n .
givenbythe cross-productofthe regionsof attractionof the i j
isolated subsystems, R , which are defined as sub-unity- j∈Ni
A,i
Thealgorithmicstepsto ascertainasymptoticstabilityare as
level sets of the corresponding (properly scaled) subsystem
outlined below:
Lyapunov functions (as in (4)-(5)), i.e.
1) We initialize e 0=v ,∀i∈{1,2,...,m}, and choose a
RA0:= RA,1×RA,2×···×RA,m (13a) sufficiently smiall e¯ ∈oiR+.
where, RA,i={xi∈Rni|Vi(xi)≤1}, ∀i (13b) 2) At the start of the k-th iteration loop, we assume to
know the scalars e 0,e 1,...,e k ,∀i, and our aim is
In presence of non-zero interactions, the resulting ROA i i i
would be different. If there exists a Lyapunov function for 5Ifthelimit condition does(cid:8)nothold, wecan(cid:9)onlyguarantee stability in
the interconnected system, the ROA for the whole system thesenseofLyapunov[1],[15]
to compute the scalars e k+1, ∀i such that (16) holds. whichcanbecheckedlocally.Ifi∈/Uk,controlisnot
i
Essentiallywewanttosolvetheoptimizationproblem, necessaryanditsetsFk≡0 andproceedstotask3.
i ni×1
If, however, i∈Uk, it proceeds to task 2 to compute
∀i, min e k+1 (17a) a control law.
ss.ikt0.,,s iktjh,eeik+c1onidition in (16) holds. (17b) 2) RIfnii∈→UR,nia, pFoikl(y0nnoi×m1i)a=ls0tantie×-1f,eeisdbcaocmkpcuotnedtrosluclahwthFaikt:
TmIfhini(si1m7is)umsioselvikine+df1eaobsvyiebrlpeethrafeotrrma0n-itgnhegiat0e,rbeaiitksieo.cnt6ionforseaarncyh ifo∈r Th(cid:209) iVsiTpr(cid:16)ofdiu+cegsi+thFeike(cid:17)q(cid:12)(cid:12)(cid:12)uniVvia=leeikn,tVSj≤OeSkj ∀cjo∈nNdi\it{iio}on,<0. (19)
(cid:2) (cid:3)
{1,2,...,m}, we conclude that the system cannot be
guaranteedto beasymptoticallystable in D,andabort −(cid:209) ViT fi+gi+Fik −r ik eik−Vi − (cid:229) s ikj e kj−Vj ∈S n¯i
the iteration. Otherwise we move on to step 3. (cid:16) (cid:17) (cid:16) (cid:17) j∈Ni\{i} (cid:16) (cid:17)
3) If eik−eik+1 ≥e¯,∀i, we continuefrom step 2 for the r ik∈Rn¯i, s ikj∈S n¯i∀j6=i, n¯i= (cid:229) nj (20)
(k+1)-thiterationloop.Otherwisewestoptheiteration j∈Ni
dec(cid:0)idingthat(cid:1)thelimitsofthesequences eik ,∀i,have 3) Finally it performs the search over minimum e k+1, as
i
been attained. Further, if the limits are all zero, we
in (17), with the un-controlled subsystem dynamics
certify asymptotic stable in D. (cid:8) (cid:9)
(f +g) in the feasibility condition (16) is replaced
i i
B. Remarks by the controlled dynamics fi+gi+Fik .
The algorithm presented in Sec. IV-A describes how one To summarize, each subsystem i computes control laws
(cid:0) (cid:1)
can determine asymptotic stability of an interconnectedsys- Fik : Rni → Rni, with Fik(0ni×1) = 0ni×1, during each k-
tem in a domain D defined by the subsystem sub-level sets. th iteration, so that the subsystem dynamics under control
This test can be performed locally, and in a parallel way, at becomes:
each subsystemlevel. The LyapunovfunctionsV’s are to be
i ∀i, ∀k, ∀x∈Dk,
found before the start of the analysis, and communicated to i
f(x)+g(x), i∈/Uk
the neighboringsubsystems. Then duringeach analysis, it is x˙ = i i i (21)
assumed that the neighboring subsystems can communicate i (cid:26) fi(xi)+gi(x)+Fik(xi), i∈Uk
with each other the computed sequences eik in real-time. where Dk were defined in (14b).
With the help of the stored Lyapunov functions, and the i
updated e k oftheneighbors,eachsubsy(cid:8)stem(cid:9)willcontinue A. Remarks
i
the iterative process outlined in Sec. IV-A. Since only the Often it is important to impose certain additional con-
(cid:8) (cid:9)
neighborinformationisrequired,thisalgorithmisreasonably straints on the possible control laws, such as bounds on
scalable with respect to the size of the full interconnected the control effort. Although control bounds can be easily
system. Moreover, the algorithm motivates the design of incorporatedin the SOSformulation,we decide to keepthat
a distributed control strategy that can ascertain asymptotic for future studies. We note, however, that since we apply
stability. controlsFk onlyoncertainsubsystemsi∈Uk,andincertain
i
V. DECENTRALIZED CONTROL domains Dk ⊆D, the control effort would be reasonably
i
bounded.
In this section, we discuss design of a local and minimal
controlstrategysuchthatthe system in (1) is asymptotically VI. RESULTS
stable in a domain D defined in (15). We use the term
Letusdescribethemodeloftheinterconnectedsystemthat
minimaltosuggestthatthecontrolbeappliedonlyincertain
we use here, and two examples to illustrate the applications
regions,andnoteverywhere,in thestate space,while by the
of the stability analysis algorithm and control design.
term local we suggest that the control be computable and
implementable on a subsystem level. A. Model Description
Weenvisionthecontroltobecomputedbyeachsubsystem
WewillconsideranetworkofnineVanderPoloscillators
at each iteration loop. At k-th iteration, ∀k ∈{0,1,2,...},
[12], as shown in Fig. 1. The dynamics of each oscillator,
the i-th subsystem, ∀i∈{1,2,...,m} performsthe following
in presence of neighbor interactions, is represented by
tasks:
1) It identifies if it belongs to the following set ∀j∈{1,2,...,9},
(cid:209) ViT(fi+gi)≥0, x˙j1=xj2 (22a)
Uk:=i∈{1,...,m}(cid:12)(cid:12)(cid:12) Vj≤e kjVi∀=j∈eikN, i\{i}  x˙j2=m jxj2(cid:0)1−x2j1(cid:1)−xj1+xj1k(cid:229)6=jz jkxk2 (22b)
 (cid:12)(cid:12)(cid:12) (18) wcihenertse,mz j’s,oarfethcehoinsteenrarcatnidoonmtelrymfsroamre(c−ho2,se0n)raannddothmelycoferoffim-
6Because Vi’s, fi’s and gi’s are polynomial, if (17) is feasible at 0-th jk
iteration foralli,thenitisalsofeasible forallsubsequentiterations. (−0.2, 0.2). It is to be noted, that the interactions need not
together oscillators {2,3} and {5,6}, as shown below
3 S :{osc 1}; N :{S ,S ,S }
1 1 2 5 7
4 S :{osc 2, osc 3}; N :{S ,S ,S }
2 2 2 1 3 6
S :{osc 4}; N :{S ,S }
3 3 2 4
S :{osc 5, osc 6}; N :{S ,S } (23)
4 4 3 7
5
1 S5:{osc 7}; N5:{S1,S6,S7}
S :{osc 8}; N :{S ,S }
6 6 2 5
S :{osc 9}; N :{S ,S ,S }
7 7 1 4 5
6
9 Then we can write the subsystem dynamics, along with the
neighborinteractions,in theformof(9).Asan example,the
7 8 states and the dynamics of S are shown below,
2
x˙ = f (x )+g (x); x =(x , x , x , x )T,
2 2 2 2 2 21 22 31 32
x
22
0.12x x −x −0.41x 1−x2
rFeiggi.o1n.s oAfantetrtawcotirokno.fnineVanderPoloscillatorsalongwiththeirisolated f2= 21 32 21x32 22 21 ,
(cid:0) (cid:1)
 0.04x31x22−x31−1.44x32 1−x231 
 0 
(cid:0) (cid:1)
−0.07x x
g2= 012 21  (24)
 0.01x12x31+0.06x42x31+0.1x82x31 
 
 
Any randomly picked initial condition, x(0)∈R0, where
A
R0 is defined in (13), can be mapped into corresponding
A
subsystem Lyapunov function level sets, g 0 =V(x(0)),∀i.
i i i
Then, by choosing v =g 0, we apply the iterative stability
oi i
analysisalgorithmtodeterminewhetherornotthedomainD
in(15)isaregionofasymptoticstability,andifnot,compute
the necessary control by (19).
B. Example: Certifiably Stable without Control
Fig.2. ComparisonofestimatedROAswithtrueROAforoscillator 1.
InFig.3(a),theevolutionsofallthestatestoanasymptot-
ically stable initial condition is shown. The subsystem Lya-
be symmetric, i.e. in general,z 6=z . Additionally,z =0 punovfunctions,showninFig.3(b),monotonicallydecrease
ik ki jk
if oscillator k is not a neighbor of oscillator j. Such choice to zero starting from the initial level sets:
of m ’s ensure that the oscillators themselves are stable,
j
with correspondingregions of attraction as shown in Fig. 1. g10=0.242, g20=0.728, g30=0.184, g40=0.658,
The regions drawn in ‘red’ around each oscillator shows its g 0=0.357, g 0=0.279, g 0=0.283 (25)
5 6 7
true ROA, while the region in ‘blue’ shows an estimate of
the ROA as sub-unity-level set of its polynomial Lyapunov Then setting v =g 0, and choosing e¯ =0.001, we run the
oi i
function,as in (5). In Fig. 2 we comparewith the true ROA iterative stability algorithm which produces the results in
the estimates obtained using a quartic Lyapunov function Table I. Atthe end ofthe 2nditeration,all the e k’sare zero,
i
and a quadratic one. Also a sequence of estimates using whichcertifiesasymptoticstabilityofthefullinterconnected
thequadraticLyapunovfunctionare shownin‘dottedblack’ system.
lines, which show how the ‘expanding interior’ algorithm
iterativelyexpandstheestimateoftheROA.Clearly,thefinal TABLEI
estimateimprovesasthedegreeofthepolynomialLyapunov ITERATIONRESULTSFORACERTIFIABLYSTABLECASE
function increases, but for this work we choose to stick to
quadratic Lyapunov functions. k ek ek ek ek ek ek ek
1 2 3 4 5 6 7
Wedecomposethesysteminto7subsystems7,bygrouping 0 0.242 0.728 0.184 0.658 0.357 0.279 0.283
1 0.024 0.031 0.147 0.188 0.004 0.010 0.127
2 0.000 0.000 0.000 0.000 0.000 0.000 0.000
7This decomposition is arbitrary. [28] presents a method of decomposi-
tion in weakly interacting subsystem, which however requires symmetric
interactions.
0.8 0.8
V
1
0.6 V
0.7 2
V
0.4 V3
0.6 4
0.2 V5
V
0.5 6
0 V
7
x(t)−0.2 V(t)0.4
−0.4
0.3
−0.6
0.2
−0.8
0.1
−1
−1.2 0
0 5 10 15 20 25 30 35 40 0 5 10 15 20 25 30 35 40
t t
(a) (b)
Fig.3. Evolutionofstates andthesubsystemLyapunovfunctions foracertifiably stableinitial condition.
1 1
V
1
0.9 V
2
V
0.5 0.8 V3
4
0.7 V5
V
6
0 0.6 V
7
x(t) V(t)0.5
−0.5 0.4
0.3
−1 0.2
0.1
−1.5 0
0 5 10 15 20 25 0 5 10 15 20 25
t t
(a) (b)
Fig.4. EvolutionofstatesandthesubsystemLyapunovfunctions foranon-certifiably-stable initial condition.
C. Example: Certifiably Stable under Control Consequently decentralized control is activated for sub-
systems S ,S and S . This results in certifiable asymptotic
1 2 5
Let us now present one example where the iterative
stability, with new e k’s shown in Table III where the ‘∗’
algorithm fails to guarantee stability, and control is applied. i
denotes presence of controllers8 which are documented in
Fig.4showsastableinitialcondition,butthealgorithmfails
(27).
tocertifystabilityforit.InFig.4(b),V (t)isseentoincrease
1
initiallybeforestartingto monotonicallydecrease.Whenwe
TABLEIII
applythe stability analysis algorithmto the initial levelsets:
ITERATIONCERTIFIESSTABILITYUNDERCONTROL
g 0=0.953, g 0=0.990, g 0=0.149, g 0=0.479,
1 2 3 4
g50=0.697, g60=0.220, g70=0.103 (26) k e1k e2k e3k e4k e5k e6k e7k
0 0.953 0.990 0.149 0.479 0.697 0.220 0.103
1 0.004∗ 0.004∗ 0.1010 0.0056 0.003∗ 0.0198 0.0380
the iteration fails at the first iteration because the algorithm
2 0.000 0.000 0.000 0.000 0.000 0.000 0.000
cannotfindfeasiblee 1’sfori=1,2,5(asshowninTableII).
i
TABLEII F11=[0,−0.569x11−2.271x12]T,
ITERATIONFAILSTOCERTIFYSTABILITY F1=[0,−1.237x −0.149x ,0,−0.285x −1.368x ]T,
2 22 21 31 32
F1=[0,−0.504x −1.539x ]T (27)
k ek ek ek ek ek ek ek 5 71 72
1 2 3 4 5 6 7
0 0.953 0.990 0.149 0.479 0.697 0.220 0.103
1 × × 0.1010 0.0056 × 0.0198 0.0380 8Inthisexample,wechosetoseeklinearcontrollers whichprovedtobe
sufficient.
0.8 1
V
1
0.6 0.9 V2
V
0.4 0.8 V3
4
0.7 V5
0.2 V
6
0.6 V
0 7
x(t) V(t)0.5
−0.2
0.4
−0.4
0.3
−0.6 0.2
−0.8 0.1
−1 0
0 5 10 15 20 25 30 0 5 10 15 20 25 30
t t
(a) (b)
Fig.5. Evolutionofstates andthesubsystemLyapunovfunctions, underdecentralized control applied atS1,S2 andS5,thatcertifies stability.
It is to be noted that we decide to apply control only on the Hence we can argue that,
dynamics equations of the states x ,x ,x and x . Fig. 5
12 22 32 72 V(t)≤e 0, ∀t≥t ,∀i
shows that under the action of the controllers in (27), all i i 0
the subsystem Lyapunov functions decrease monotonically =⇒∃t1:=maxt1, s.t. V(t)≤e 1,∀t≥t1,∀i (31)
i i i
i
to zero.
Following similar arguments it is easy to show that,
VII. CONCLUSIONS
In this work, we present an algorithmic approach to cer- Vi(t)≤ei0, ∀t≥t0, ∀i
tify asymptotic stability of an interconnected system whose =⇒∀k, ∃tk≥t , s.t. V(t)≤e k, ∀t≥tk, ∀i (32)
0 i i
dynamics can be expressed in polynomial form. We also
propose the design of decentralized control laws when such Finally combining (28) and (32) we observe,
acertificationisnotpossible.Theapproachpresentedhereis
parallel and scalable. Similar method can also be applicable ∀d ∈ 0,minei0 , ∃tK ≥t0, s.t. Vi(t)<d , ∀t≥tK, ∀i
i
to complex real world systems, such as the power system. (cid:18) (cid:21)
While power system dynamics are non-polynomial, if those which concludes the proof, because of (12).
are transformed into polynomial forms, by introducing ad-
ditional equality constraints [14], the methods developed in REFERENCES
this article can be applied. Future work need to address the
[1] A. M. Lyapunov, The General Problem of the Stability of Motion.
issues such as including bounds on the control effort, and Khatkov, Russia:KharkovMath.Soc.,1892.
relaxingtherequirementofmonotonicdecreaseofsubsystem [2] W.M.HaddadandV.Chellaboina,NonlinearDynamicalSystemsand
Control. Princeton, NewJersey:Princeton University Press,2008.
Lyapunov functions along the flow.
[3] D. Siljak, “Stability of large-scale systems under structural pertur-
bations,” Systems, Man and Cybernetics, IEEETransactions on, vol.
APPENDIX SMC-2,no.5,pp.657–663, Nov1972.
A. Proof of Lemma 1 [4] S.Weissenberger,“Stabilityregionsoflarge-scalesystems,”Automat-
ica,vol.9,no.6,pp.653–663, 1973.
We note that since limk→+¥ eik=0,∀i, [5] A.N.Michel, “Onthe status ofstability ofinterconnected systems,”
AutomaticControl,IEEETransactionson,vol.28,no.6,pp.639–653,
1983.
∀d ∈ 0,mine 0 , ∃K, s.t. e k<d ∀k>K,∀i. (28)
i i [6] M. Araki, “Stability of large-scale nonlinear systems quadratic-
i
(cid:18) (cid:21) order theory of composite-system method using m-matrices,” IEEE
Let us assume, without any loss of generality, that Transactions on Automatic Control, vol. 23, no. 2, pp. 129 – 142,
1978.
m [7] R. Bellman, “Vector Lyapunov functions,” Journal ofthe Society for
∃t ≥0, s.t. x(t )∈ x∈Rn V(x)≤e 0 (29) Industrial & Applied Mathematics, Series A: Control, vol. 1, no. 1,
0 0 i i i
( (cid:12) ) pp.32–34,1962.
(cid:12)i\=1 [8] F. N. Bailey, “The application of Lyapunov’s second method to
(cid:12)
Then, (cid:12) interconnectedsystems,”J.SIAMControl,vol.3,pp.443–462,1966.
t (cid:12) [9] R. Conti, “Sulla prolungabilita` delle soluzioni di un sistema di
∀i, V(t)=e 0+ V˙(t )dt , ∀t≥t (30a) equazioni differenziali ordinarie.” Bollettino dell’Unione Matematica
i i Zt0 i 0 Italiana, vol.11,no.4,pp.510–514,1956.
=⇒∃t1<t + e 1−e 0 /r¯1, r¯1:= supV˙(x)<0 (30b) [10] F. Brauer, “Global behavior of solutions of ordinary differential
i 0 i i i i i equations,”JournalofMathematicalAnalysisandApplications,vol.2,
x∈Di1 no.1,pp.145–158,1961.
s.t. V(t)<(cid:0)e 1, ∀t≥(cid:1)t1 (30c) [11] E.F.Beckenbach andR.Bellman,“Inequalities,” Spring-Verlag,New
i i i York/Berlin, 1961.
[12] B.VanderPol,“Onrelaxation-oscillations,” TheLondon,Edinburgh, http://fewcal.kub.nl/sturm/software/sedumi.html.
and Dublin Philosophical Magazine and Journal of Science, vol. 2, [21] A.PapachristodoulouandS.Prajna,“OntheconstructionofLyapunov
no.11,pp.978–992,1926. functionsusingthesumofsquaresdecomposition,”inProceedingsof
[13] Z.W.Jarvis-Wloszek,“Lyapunovbasedanalysisandcontrollersynthe- theIEEEConferenceonDecisionandControl,Dec.2002,pp.3482–
sisforpolynomialsystemsusingsum-of-squaresoptimization,”Ph.D. 3487.
dissertation, University ofCalifornia, Berkeley, CA,2003. [22] Z. J. Wloszek, R. Feeley, W. Tan, K. Sun, and A. Packard, Positive
[14] M. Anghel, F. Milano, and A. Papachristodoulou, “Algorithmic con- Polynomials in Control. Berlin, Heidelberg: Springer-Verlag, 2005,
struction ofLyapunov functions forpowersystemstability analysis,” ch.ControlApplications ofSumofSquaresProgramming,pp.3–22.
CircuitsandSystemsI:RegularPapers,IEEETransactionson,vol.60, [23] A. Papachristodoulou and S. Prajna, Positive Polynomials in Con-
no.9,pp.2533–2546,Sept2013. trol. BerlinHeidelberg: Springer-Verlag, 2005,ch.Analysisofnon-
[15] J.-J.E.Slotine,W.Li,etal.,Appliednonlinearcontrol. Prentice-Hall polynomialsystemsusingthesumofsquaresdecomposition, pp.23–
Englewood Cliffs,NJ,1991,vol.199,no.1. 43.
[16] P. A. Parrilo, “Structured semidefinite programs and semialgebraic [24] ——,“Atutorialonsumofsquarestechniques forsystemsanalysis,”
geometrymethodsinrobustnessandoptimization,”Ph.D.dissertation, inProceedingsofthe2005AmericanControlConference,June2005,
Caltech, Pasadena, CA,2000. pp.2686–2700.
[17] W. Tan, “Nonlinear control analysis and synthesis using sum-of- [25] M. Putinar, “Positive polynomials on compact semi-algebraic sets,”
squares programming,” Ph.D. dissertation, University of California, IndianaUniversityMathematicsJournal,vol.42,no.3,pp.969–984,
Berkeley, CA,2006. 1993.
[18] A. Papachristodoulou, J. Anderson, G. Valmorbida, S. Prajna, [26] J.-B.Lasserre,Moments,positivepolynomials andtheirapplications.
P. Seiler, and P. A. Parrilo, “SOSTOOLS: Sum of squares WorldScientific, 2009,vol.1.
optimization toolbox for MATLAB,” 2013, available from [27] M.Anghel,J.Anderson,andA.Papachristodoulou,“Stabilityanalysis
http://www.eng.ox.ac.uk/control/sostools. of power systems using network decomposition and local gain anal-
[19] S.Prajna, A. Papachristodoulou, P.Seiler, andP.A.Parrilo, Positive ysis,” in BulkPower System Dynamics and Control-IX Optimization,
Polynomials in Control. Berlin, Heidelberg: Springer-Verlag, 2005, SecurityandControloftheEmergingPowerGrid(IREP),2013IREP
ch.SOSTOOLSandItsControl Applications, pp.273–292. Symposium. IEEE,2013,pp.1–7.
[28] J. Anderson and A. Papachristodoulou, “A decomposition technique
[20] J. F. Sturm, “Using SeDuMi 1.02, a MATLAB toolbox for opti-
for nonlinear dynamical system analysis,” IEEE Transactions on
mization over symmetric cones,” Optimization Methods and Soft-
AutomaticControl, vol.57,pp.1516–1521,June2012.
ware, vol. 11-12, pp. 625–653, Dec. 1999, software available at