Document [original]

Robust stability of differential-algebraic equations

Nguyen Huu Du and Vu Hoang Linh∗and Volker Mehrmann†

Abstract This paper presents a survey of recent results on the robust stability analysis and the distance

to instability for linear time-invariant and time-varying differential-algebraic equations (DAEs). Different

stability concepts such as exponential and asymptotic stability are studied and their robustness is analyzed

under general as well as restricted sets of real or complex perturbations. Formulas for the distances are

presented whenever these are available and the continuity of the distances in terms of the data is discussed.

Some open problems and challenges are indicated.

Keywords: differential-algebraic equation, robust stability, stability radius, spectrum, spectral set, re-

stricted perturbation.

AMS(MOS) subject classification: 93B35, 93D09, 34A09, 34D10

1 Introduction

In many areas of science and engineering one uses mathematical models to simulate, control or optimize a

system or process. These mathematical models, however, are typically inexact or contain uncertainties and

thus, the following question is of major importance.

How robust is a specific property of a given system described by differential or difference equations under

perturbations to the data?

Here, we say that a certain property of a system is robust if it is preserved when an arbitrary (but sufficiently

small) perturbation affects the system. An important quantity in this respect is then the distance (measured

by an appropriate metric) between the nominal system and the closest perturbed system that does not

possess the mentioned property, this is typically called the radius of the system property.

Nguyen Huu Du

Faculty of Mathematics, Mechanics and Informatics, Vietnam National University, 334, Nguyen Trai Str., Thanh Xuan, Hanoi,

Vietnam e-mail: [email protected]

Vu Hoang Linh

Faculty of Mathematics, Mechanics and Informatics, Vietnam National University, 334, Nguyen Trai Str., Thanh Xuan, Hanoi,

Vietnam e-mail: [email protected]

Volker Mehrmann

Institut f¨

ur Mathematik, MA 4-5, Technische Universit¨

at Berlin, D-10623 Berlin, Fed. Rep. Germany e-mail: mehrmann@

math.tu-berlin.de

∗Supported by NAFOSTED Grant 101.01-2011.14.

†Supported by Deutsche Forschungsgemeinschaft, through project A02 within Collaborative Research Center 910 Control

of self-organizing nonlinear systems

2 Nguyen Huu Du and Vu Hoang Linh and Volker Mehrmann

In this paper, we deal with robustness and distance problems for differential-algebraic equations (DAEs),

with a focus on robust stability and stability radii. Systems of DAEs, which are also called descriptor

systems in the control literature, are a very convenient modeling concept in various real-life applications

such as mechanical multibody systems, electrical circuit simulation, chemical reactions, semi-discretized

partial differential equations, and in general for automatically generated coupled systems, see [7, 33, 41,

56, 60, 76, 77] and the references therein.

DAEs are generalizations of ordinary differential equations (ODEs) in that certain algebraic equations

constrain the dynamical behavior. Since the dynamics of DAEs is constrained to a set which often is only

given implicitly, many theoretical and numerical difficulties arise, which may lead to a sensitive behavior

of the solution of DAEs to perturbation in the data. The difficulties are characterized by fundamental

notions for DAEs such as regularity, index, solution subspace, or hidden constraints, which do not arise

for ODEs. These properties may be easily lost when the data are subject to arbitrarily small perturbations.

As a consequence, usually restrictions to the allowed perturbations have to be made, leading to robustness

questions for DAEs that are very different from those for ODEs.

This paper surveys robustness results for linear DAEs with time-invariant or time-varying coefficients

of the form

E˙x=Ax +f,(1)

on the half-line I= [0,∞), together with an initial condition

x(t0) = x0,t0∈I.(2)

Here we assume that E,A∈C(I,Kn×n), and f∈C(I,Kn)are sufficiently smooth. We use the notation

C(I,Kn×n)to denote the space of continuous functions from Ito Kn×n, where K=Ror K=C.

Linear systems of the form (1) arise directly in many applications and via linearization around solution

trajectories [16]. They describe the local behavior in the neighborhood of a solution for general implicit

nonlinear system of DAEs

F(t,x,˙x) = 0,(3)

the constant coefficient case arising in the case of linearization around stationary solutions.

Recall the following classical stability concepts for ordinary differential equations

˙x=f(t,x),t∈I,(4)

with initial condition (2), see e.g. [51].

Definition 1. A solution x:t7→ x(t;t0,x0)of (4) with initial condition (2) is called

1. stable if for every ε>0 there exists δ>0 such that

a. the initial value problem (4) with initial condition x(t0) = ˆx0is solvable on Ifor all ˆx0∈Knwith

kˆx0−x0k<δ;

b. the solution x(t;t0,ˆx0)satisfies kx(t;t0,ˆx0)−x(t;t0,x0)k<εon I.

2. asymptotically stable if it is stable and there exists ρ>0 such that

a. the initial value problem (4) with initial condition x(t0) = ˆx0is solvable on Ifor all ˆx0∈Knwith

kˆx0−x0k<ρ;

b. the solution x(t;t0,ˆx0)satisfies limt→∞kx(t;t0,ˆx0)−x(t;t0,x0)k=0.

3. exponentially stable if it is stable and exponentially attractive, i.e., if there exist δ>0, L>0, and γ>0

such that

a. the initial value problem (4) with initial condition x(t0) = ˆx0is solvable on Ifor all ˆx0∈Knwith

kˆx0−x0k<δ;

b. the solution satisfies the estimate kx(t;t0,ˆx0)−x(t;t0,x0)k<Le−γ(t−t0)on I.

If δdoes not depend on t0, then we say the solution is uniformly (exponentially) stable.

Robust stability of differential-algebraic equations 3

Note that one can transform the ODE (4) in such a way that a given solution x(t;t0,x0)is mapped to the

trivial solution by simply shifting the arguments. When studying the stability of a selected solution, one

may therefore assume without loss of generality that the selected solution is the trivial solution, and also

that t0=0.

One can immediately extend Definition 1 verbatim to DAEs. However, one has to be careful with the

initial conditions and the inhomogeneities, since they are restricted due to the algebraic constraints in the

system. This is, in particular, true if one considers the robustness of the stability concepts under perturba-

tions to the system.

The following examples give an illustration for the possible difficulties in the robustness of the stability

concepts for DAEs under small perturbations.

Example 1. Consider the homogeneous linear time-invariant DAE

1 0

0 0 ˙x1

˙x2=0 1

1 0 x1

x2,(5)

which can be written as ˙x1=x2,0=x1and has only the trivial solution x1=x2=0.

If we perturb (5) by a small εas

1 0

0 0 ˙x1

˙x2=0 1

1εx1

x2,(6)

then solving the second equation of (6) for x2and substituting into the first equation, we obtain

˙x1=−(1/ε)x1.(7)

Clearly, if ε<0, then the perturbed DAE (6) is unstable. If ε>0, then the system is asymptotically stable,

but it qualitatively differs from the solution of the original system (5). For an arbitrarily prescribed initial

value x1(0)6=0, the initial value problem for (7) has a unique solution. Furthermore, the value of x2(0)is

not required and is uniquely determined by x1(0). In fact, this small perturbation has changed the index of

the DAE (5), see Definition 3 below.

If we add an inhomogeneity to these DAEs, then more essential differences appear.

In Example 1 the perturbation is affecting only the coefficient A. The situation is even more complicated if

we allow perturbations in the coefficient of ˙x.

Example 2. Consider the well-known singularly perturbed system

In10

0εIn2˙x1

˙x2=A11 A12

A21 A22 x1

x2,(8)

where Ai j,i,j∈ {1,2}, are constant matrices of appropriate sizes, and ε>0 is a small parameter. Let us

assume that A22 is invertible. If εis set to 0, then the leading matrix becomes singular, i.e., we have a DAE,

and we can solve the second equation for x2, and obtain the so-called underlying ODE

˙x1= (A11 −A12A−1

22 A21)x1.

It is well known that for sufficiently small ε, the (asymptotic) stability of (8) depends not only on the

stability of the so-called slow sub-system associated with the underlying ODE, but also on that of the so

called fast sub-system ˙x2=A22x2, see [26, 29, 55], associated with the algebraic equation.

In this example, the rank of the leading matrix is changed, when εmoves from zero to a nonzero value.

In the case ε=0, the initial condition must be consistent to ensure the existence of a solution, but obviously

this is not required in the case of a nonzero ε. The difficulties increase if A22 is singular and/or the leading

matrix involves a singular perturbation of a more general structure.

The presented examples already indicate some of the possible difficulties which will be discussed in this

paper. We present the analysis of robust exponential and asymptotic stability for linear DAEs with time-

invariant or time-varying coefficients. This is a relatively young topic, starting with the work in [15, 75],

4 Nguyen Huu Du and Vu Hoang Linh and Volker Mehrmann

which generalized results on the distance to instability and the concept of stability radius for ODEs in [46]

and [84].

The outline of the paper is as follows. In Section 2 we summarize recent results on the robust stability

and stability radii for linear time-invariant DAEs. In Section 3 we study the robust stability of linear time-

varying DAEs. Stability radii, their dependence on the data, and the robustness of stability spectra are

analyzed. Some discussions and topics for future research close the paper.

2 Robust stability of linear time-invariant DAEs

In this section we study homogeneous linear time-invariant DAEs of the form

E˙x=Ax,t∈I,(9)

where Eand Aare given constant matrices in Kn×n,K=Cor R.

Definition 2. A matrix pair (E,A),E,A∈Kn×nis called regular if there exists λ∈Csuch that the deter-

minant of (λE−A), denoted by det(λE−A), is different from zero. Otherwise, if det(λE−A) = 0 for all

λ∈C, then we say that (E,A)is singular.

If (E,A)is regular, then a complex number λis called a (generalized finite) eigenvalue of (E,A)if

det(λE−A) = 0. The set of all (finite) eigenvalues of (E,A)is called the (finite) spectrum of the pencil

(E,A)and denoted by σ(E,A). If Eis singular and the pair is regular, then we say that (E,A)has the

eigenvalue ∞.

In the following we only consider regular pairs (E,A). Such pairs can be transformed to Weierstraß-

Kronecker canonical form, see [7, 35, 37], i.e., there exist nonsingular matrices W,T∈Cn×nsuch that

E=WIr0

0NT−1,A=WJ0

0In−rT−1,(10)

where Ir,In−rare identity matrices of indicated size, J∈Cr×r, and N∈C(n−r)×(n−r)are matrices in Jordan

canonical form and Nis nilpotent. If Eis invertible, then r=n, i.e., the second diagonal block does not

occur.

Definition 3. Consider a regular pair (E,A)with E,A∈Kn×nin Weierstraß-Kronecker form (10). If r<n

and Nhas nilpotency index ν∈ {1,2,...}, i.e., Nν=0,Ni6=0 for i=1,2,...,ν−1, then νis called the

index of the pair (E,A)and the associated DAE (9) and we write ind(E,A) = ν. If r=nthen the DAE has

index ν=0.

The finite spectrum σ(E,A)is given by the spectrum of σ(J)in the Weierstraß-Kronecker form (10)

and it is easy to verify that for the degree of the characteristic polynomial degdet(λE−A) = rankE=r

holds if and only if ind(E,A)≤1.

For the DAE (9) with initial condition (2), there always exists a projection matrix P∈Kn×rso that the

projected initial condition P(x(t0)−x0) = 0 is consistent, i.e., the DAE (9) with this initial condition has

a unique solution, see [7, 37]. Shifting again the solution to the trivial solution x=0, for a DAE of the

form (9) with regular pair (E,A),E,A∈Kn×nwe say that this solution is exponentially stable if there exist

δ>0, L>0, and γ>0 such that the initial value problem

E˙x=Ax,P(x(0)−ˆx0) = 0,

is solvable on Ifor all ˆx0∈Knwith kˆx0k<δ, and the solution satisfies the estimate kx(t;t0,ˆx0)k<

Le−γ(t−t0)on I. If the trivial solution is exponentially stable, then we say that (9) is exponentially stable.

We remark that for linear time-invariant systems, the concept of exponential stability is equivalent to that

of asymptotic stability and hence we do not have to distinguish these concepts, and discuss only asymptotic

stability as is usually done in the literature.

Using the canonical form (10), one can easily verify the following statement, see e.g., [15].

Robust stability of differential-algebraic equations 5

Proposition 1. Consider a DAE of the form (9) with regular pair (E,A), E,A∈Kn×n. System (9) is asymp-

totically stable if and only if the pair (E,A)is asymptotically stable, i.e., the finite spectrum satisfies

σ(E,A)⊂C−, where C−denotes the open left-half complex plane.

After introducing the basic notation, in the next subsection we discuss the stability radius of a DAE.

2.1 Stability radii for lime time-invariant DAEs

In this section we study the behavior of the finite spectrum of a regular pair (E,A)under structured pertur-

bations in the matrices E,A. Suppose that the system (9) is asymptotically stable and consider a perturbed

system

(E+B1∆1C1)˙x= (A+B2∆2C2)x,(11)

where ∆i∈Kmi×qi(i=1,2) are perturbations and Bi∈Kn×mi,Ci∈Kqi×nare given matrix pairs that restrict

the structure of the perturbations. The matrix pair (B1∆1C1,B2∆2C2)is called a structured perturbation.

For simplicity, let us consider the case that the restricting matrices satisfy C1=C2=C(alternatively we

could consider the other simplifying case that B1=B2). Set

∆=∆1

∆2,B=B1B2,

and introduce m=m1+m2and q=q1=q2. Then we consider the set of destabilizing perturbations

VK(E,A;B,C) = ∆∈Km×q,(11)is singular or not asymptotically stable

and have the following definition.

Definition 4. The structured stability radius of (9) subject to structured perturbations as in (11) is defined

rK(E,A;B,C) = inf{k∆k,∆∈VK(E,A;B,C)},

where k.kis a matrix norm induced by a vector norm. Depending on K=Cor K=R, we talk about the

complex or the real stability radius, respectively.

Obviously, we have the estimate

rC(E,A;B,C)≤rR(E,A;B,C).

To obtain a computable formula for the complex stability radius, one introduces the matrix functions

G1(s) = −sC(sE −A)−1B1,G2(S) = C(sE −A)−1B2,G(s) = G1(s)G2(s),

for s∈C. Denoting by iRthe imaginary axis of the complex plane, the following result is analogous to that

for linear time-invariant ODEs of [47].

Theorem 1. Suppose that the matrix pencil (E,A)is regular and asymptotically stable. Then with respect

to any matrix norm induced by a vector norm, the complex stability radius of (9) has the representation

rC(E,A;B,C) = sup

s∈iR

kG(s)k−1

.(12)

Proof. The proof can be obtained by using the same techniques as in [27, 28, 32].

Note that in [27, 28] the complex structured stability radius is considered with respect to perturbations

either in Eor in A, i.e., either B1=0 or B2=0. In [32], formula (12) has been proven for perturbations in

both E,A. Note further that these papers discuss both the continuous-time and the discrete-time case. Fur-

thermore, it has been shown in [32] that it is always possible to construct a rank 1 destabilizing perturbation

∆, with a norm that approximates the value of rCwithin an arbitrarily prescribed accuracy.

6 Nguyen Huu Du and Vu Hoang Linh and Volker Mehrmann

Remark 1. The concept of stability radius can be extended to more general sets. Suppose that all the eigen-

values of the unperturbed matrix pencil lie in a prescribed open subset Cgof the complex plane. Then we

want to determine the largest perturbations that the system can tolerate so that its spectrum remains in Cg.

In the asymptotic stability analysis of differential equations, the open subset Cgis chosen to be C−. As in

the other cases, it is trivial to obtain a formula of a Cg-stability radius analogously to (12). In fact, we sim-

ply replace iRby the boundary set of Cb=C\Cg. As a consequence of the definition, the strict positivity

of a Cg-stability radius with a relevant subset Cgimplies the continuity of the spectrum with respect to the

data.

Unlike for the complex stability radius, a general formula for the real stability radius measured by an

arbitrary matrix norm is not available. However, if we consider as vector norm the Euclidean norm, then

a computable formula has been obtained in [74]. This formula is based on the notion of real/complex

structured singular values, which for a given M∈Kp×m, are defined by

µK(M) = inf{σ1(∆),∆∈Km×p,and det(I−∆M) = 0}−1,

respectively, depending on K=Cor K=R. Here σ1(∆)denotes the largest singular value, see [36], of

the matrix ∆.

Clearly, if Mis real, then the complex and the real structured singular values coincide. While for the

complex structured singular values, the formula µC(M) = σ1(M)follows trivially, for the real case the

formula is more sophisticated.

Proposition 2. [74] The real structured singular value of M ∈Km×qis given by

µR(M) = inf

γ∈(0,1]σ2ReM−γImM

γImMReM,

where σ2(A)denotes the second largest singular value of A.

Using the definition of the stability radius with respect to the Euclidean norm, we thus have that

rK(E,A;B,C) = inf{σ1(∆),∆∈VK(E,A;B,C)}

=infs∈iRinf{σ1(∆),∆∈Km×qand det(s(E+B1∆1C)−(A+B2∆2C)) = 0}

=infs∈iRinf{σ1(∆),∆∈Km×qand det(I+(sE −A)−1(sB1∆1C−B2∆2C)) = 0}

=infs∈iRinf{σ1(∆),∆∈Km×qand det(I−∆G(s)) = 0}

={sups∈iRµK(G(s))}−1,

and hence we obtain the following theorem.

Theorem 2. [74] Suppose that the matrix pencil (E,A)is regular and asymptotically stable. Then the

structured stability radii of (9), measured in Euclidean norm, are given by

rC(E,A;B,C) = sup

s∈iR

σ1(G(s))−1

,(13)

and

rR(E,A;B,C) = sup

s∈iR

inf

γ∈(0,1]σ2ReG(s)−γImG(s)

γImG(s)ReG(s)−1

,(14)

respectively.

It is important to note that the presented results on the structured stability radii do not reflect the fact that

an eigenvalue at ∞may become finite or conversely a finite eigenvalue may move to ∞, i.e., it may happen

that the index of the pair (E,A)changes or that the pair becomes singular.

An increase of the index is not problematic for homogeneous systems with zero initial conditions, since

the part of the solution associated with the infinite spectrum is always 0. In the case of inhomogeneous

systems or nonzero initial conditions, however, an increase of the index may lead to a loss of solvability

Robust stability of differential-algebraic equations 7

of the equation due to inconsistent initial values or a lack of smoothness for the inhomogeneity. As we

have demonstrated in Examples 1 and 2, this can even happen with infinitesimally small perturbations.

Furthermore, while the stability radii of an ODE are always strictly positive, those of a DAE may be zero.

To see this, considering a pair in canonical form (10), we have

(sE −A)−1=T(sIr−J)−10

0−∑k−1

i=0(sN)iW−1.

It then obvious that if N6=0, then ||G1(s)|| and ||G2(s)|| may tend to ∞as |s| → ∞, which implies that rC=0.

Hence, the perturbations in (11) must be further restricted so that the stability radii are strictly positive.

Partitioning the restriction matrices C,B(after transformation to Weierstraß-Kronecker form) as

CT = [C1C2],W−1B1=B11

B12 ,W−1B2=B21

B22 ,(15)

according to the structure of (10), it is easy to see that if ind(E,A) = 1 then

sup

s∈iR

kG1(s)k<∞if and only if C2B12 =0

and if ind(E,A)>1 then

sup

s∈iR

kG2(s)k<∞if and only if C2B12 =0 and C2B22 =0.

These observations are summarized in the following result.

Proposition 3. Consider a regular pair (E,A)and the associated DAE of the form (9). If ind(E,A) = 1,

then the structured stability radii of (9) are strictly positive if and only if C2B12 =0. If ind(E,A)>1, then

the structured stability radii of (9) are strictly positive if and only if C2B12 =0and C2B22 =0, where the

transformed structure matrices are defined by (15).

Moreover, if C2is of full rank, then rK(E,A,B,C)>0if and only if B12 =0for the case ind(E,A) = 1

and B12 =B22 =0for the case ind(E,A)>1.

It follows that the perturbations have to be further restricted by choosing

B1=WB11

0,(16)

if ind(E,A) = 1 and

B1=WB11

0,B2=WB21

0(17)

if ind(E,A)>1.

Definition 5. [15] A structured perturbation as in (11) is called admissible if it does not alter the nilpotency

structure of the Weierstrass-Kronecker form (10) of (E,A), i.e., the nilpotent matrix Nand the correspond-

ing left invariant subspace associated with the eigenvalue ∞are preserved.

In the case that ind(E,A) = 1, one has the following characterization of admissible perturbations.

Proposition 4. [15] Consider the regular DAE (9) with ind(E,A) = 1, subject to a general unstructured

perturbation

(E+F)˙x= (A+H)x.

Then there exist an orthogonal matrix P and a permutation matrix Q such that

PEQ =E11 E12

0 0 ,PAQ =A11 A12

A21 A22 ,

8 Nguyen Huu Du and Vu Hoang Linh and Volker Mehrmann

where rank[E11,E12] = rankE=r, and rankA22 =n−r. Furthermore, if (F,H)is an admissible perturba-

tion, then

PFQ =F11 F12

0 0 ,PHQ =H11 H12

H21 H22 .

Note that in Proposition 4, the transformation by the matrices P,Qdoes not change the structure, the sta-

bility, and consequently, the stability radii of (9). Note further that Proposition 3 can also be used to char-

acterize admissible perturbations for the case ind(E,A)>1.

After these observations we can introduce the distance to the nearest pair with a different nilpotency

structure

dK(E,A;B,C) = infk∆k,∆∈Km×qand (11) does not preserve the nilpotency structure,

and obtain the following result, see [5].

Theorem 3. Consider a regular DAE with Weierstraß-Kronecker form (10), subject to transformed per-

turbations satisfying (16) and (17). Then the distance to the nearest system with a different nilpotency

structure is given by

dK(E,A;B,C) = C1B11 C2B22 −1

if ind(E,A) = 1and

dK(E,A;B,C) = ||C1B11||−1,

if ind(E,A)>1. Moreover, dC(E,A;B,C) = dR(E,A;B,C).

Proof. The proof is similar to that for the stability radii, see Remark 1. The nilpotency structure of the

perturbed system (11) is preserved if and only if the perturbed matrix

I+B11∆1C1B11∆1C2

B22∆2C2I+B22∆2C2

is nonsingular in the case ind(E,A) = 1 and if I+B11∆1C1is nonsingular in the case ind(E,A)>1. Thus,

in both cases we have a problem of characterizing the distance of a matrix to singularity, which we obtain

by taking Cg=C\{0}. For the equality of the complex and the real stability radii, we note that if the data

are real, then the smallest perturbation that makes a matrix singular can always be chosen to be real [36].

We summarize the results in the following theorem.

Theorem 4. Consider a regular and asymptotically stable DAE (9) with Weierstraß-Kronecker form (10).

Let the perturbation structure satisfy (16) and (17) for index ind(E,A) = 1and ind(E,A)>1, respectively.

If ||∆|| <rK(E,A,B,C), then the perturbed DAE (11) preserves not only the regularity and the asymptotic

stability, but also the nilpotency structure of the DAE (9). In other words, we obtain the structured stability

radii as

rK(E,A;B,C) = inf{k∆k,∆∈V(E,A;B,C)or (11) has different nilpotency structure}.

These structured stability radii are given by the formulas (12), (13), and (14), respectively.

Proof. Assuming (16) and (17), respectively, it is easy to check that

lim

s→∞G1(s) = C1B11,and lim

s→∞G2(s) = C2B22.

Since sups∈iR||G(s)|| ≥||G(∞)||, the statement for the case K=Cfollows directly from Theorem 1 and

Theorem 3.

If all the data are real, then G(∞)is real as well. Since the complex and the real structured singular

values of G(∞)are equal, again by Theorem 3, the same statement holds for the case K=R.

Robust stability of differential-algebraic equations 9

Remark 2. The first results for stability radii for linear time-invariant DAEs of index one were given in [75]

and [15]. In [75] only unstructured perturbations in Awere considered and the formula for the unstructured

complex stability radius measured in Euclidean norm is exactly a special case of (13) with restriction ma-

trices B1=0 and B2=C=I. A more general result was obtained in [15], where the authors considered

admissible perturbations as in Proposition 4 and formulated the complex stability radius using the Frobe-

nius norm which is not a matrix norm induced by a vector norm. However, using the fact that the Frobenius

norm gives an upper bound for the Euclidean norm, the stability radii in Frobenius norm are upper bounds

for those in Euclidean norm. On the other hand, in the proof of Theorem 1, a rank one destabilizing pertur-

bation can be constructed whose (Euclidean) norm approximates the true value of the stability radius with

an arbitrarily small accuracy. Since the Euclidean and the Frobenius norm of a rank one matrix are equal,

the formula of the complex stability radius given in [15] and (13) yield the same value, i.e., the former one

can be considered as a special case of (12) and (13).

Remark 3. In many applications, one does not only have static perturbations as in (11), but also linear time-

varying, nonlinear or even nonlinear dynamic perturbations. In the context of regular DAEs (9) of index

at most one, if perturbations are admitted in Aonly, then it is possible to extend the concept of structured

stability radii with respect to linear time-varying, nonlinear, or nonlinear dynamic perturbations. It can be

shown that all the complex structured stability radii with respect to different classes of perturbations are

equal as it is well known in the ODE case [49], see also Corollary 2 below.

Remark 4. Numerical algorithms for computing the stability radii for ODEs are proposed in a number of

works, e.g., see [4, 6, 12, 34, 38, 39, 40, 42, 45, 71, 72, 73, 81]. Some extensions to DAEs are discussed

in [5]. Since the robust stability of a linear time-invariant system is closely related to the sensitivity of the

spectrum, the characterization and the computation of stability radii is also very closely related to the topics

of spectral value sets [44, 50] and pseudospectra [11, 83].

Remark 5. Robustness questions can also be discussed for other fundamental concepts of control theory

such as controllability, observability, stabilizability, or detectability. These concepts have been extended to

DAEs in many different publications, see e.g., [8, 9, 10, 17, 21, 68, 70]. It is natural to analyze the robustness

of these properties when the control systems are subject to uncertain perturbations, which leads to similar

distance problems as for robust stability. For ODEs, such distance problems are extensively investigated in

a number of works, e. g., see [52, 71] and the references therein. Results on the controllability radius for

linear time-invariant descriptor systems are given in [13, 58, 59, 80].

It is well known that solutions of DAEs are more sensitive to data than those of ODEs. This topic

has been discussed in [67] for a perturbed index two DAE in semi-explicit form and in [24, 82] for general

singularly perturbation problems of DAEs. But a general perturbation theory for linear time-invariant DAEs

is still open. This is partly due to the fact that no complete characterization of the distance to the nearest

singular pencil is available [14].

2.2 Dependence of stability radii on the data

In view of the numerical computation of the stability radii, a natural question is whether the structured

stability radii rK(E,A;B,C)depend continuously on the data E,A,B,C. In the ODE case, i.e., E=Iand

if only Ais perturbed, it was shown in [48] that the complex structured radius depends continuously on

data, but the real one does not. Extending these results to DAEs, it follows that the complex structured

stability radius rC(E,A;B,C)depends continuously on the data, provided that the nilpotency structure is

preserved, i.e., we are restricted to the set of DAEs (9) that have the same nilpotency structure and to the

set of admissible perturbations. In [29, 30] the robust stability of the parameterized DAE system

(E+εF)˙x=Ax,(18)

is considered, where ε>0 is a small parameter and the unperturbed DAE (ε=0) is assumed to be regular,

of index at most one and asymptotically stable. The classical singularly perturbed system (8) is a special

10 Nguyen Huu Du and Vu Hoang Linh and Volker Mehrmann

case of this more general system. If εFbelongs to the class of admissible perturbations characterized by

Proposition 4, then it is easily shown that the complex structured stability radius depends continuously on

the parameter ε. This, however, is not the case when the appearance of εFchanges the index and/or the

number of finite eigenvalues, i.e., the nilpotency structure of (E,A).

Sufficient conditions can be given to ensure that (18) is asymptotically stable for all sufficiently small

and positive ε. Namely, if the unperturbed DAE (ε=0) and the fast subsystem (that is associated with the

algebraic part of the DAE) are simultaneously asymptotically stable, then for all sufficiently small ε, so

is the parameterized system (18). Furthermore, the complex structured stability radius of (18) converges

to the minimum of the stability radius of the reduced system and that of the fast subsystem. This implies

that the stability radius of (18) may be discontinuous at ε=0, when the nilpotency structure is no longer

invariant. As a special case, the result for the robust stability of (8), investigated in [25] by a different

approach, then follows immediately.

Remark 6. In [61], the robust stability of a DAE subject to perturbations of the form E˙x= (A+εH)xis

considered, where E,Aare given as in (9), His a given matrix, and εis a uncertain parameter. Assuming that

the unperturbed system is regular, of index at most one and asymptotically stable, a computable formula

for the maximal stability interval (ε1,ε2)is derived, i.e., the perturbed system retains the index and is

asymptotically stable for all ε∈(ε1,ε2).

A somewhat more general and extensive analysis of stability radii for DAEs is given in [5], considering

both the structure and the spectrum, which are treated separately based on the notion of structured singular

values. This approach makes the characterization of stability radii for higher index DAEs possible as well

as that for the special uncertainty structure of affine perturbations. The latter one also extends the result

in [61] to the case when both coefficient matrices are perturbed with a one-parameter family. In addition,

the work in [5] is partially devoted to robust stability of second order DAEs with applications in electrical

networks.

In [28], complex structured stability radii for the discrete-time analogue of DAEs, i.e., singular differ-

ence equations are analyzed. In particular, the complex structured stability radius of the discretized system

of (9) using the implicit Euler method is shown to converge to the corresponding one of the continuous-time

system as the stepsize tends to zero.

It has been shown in [79] that for the class of positive systems the complex and the real stability radii

coincide and they are easily computable. Attempts to extend this result to DAEs and to other implicit

dynamic equations are presented in [27] and [32], respectively.

3 Robust stability of linear time-varying DAEs

In this section, we investigate the exponential and asymptotic stability and its robustness for linear time-

varying DAEs of the form

E˙x=Ax,t∈I,(19)

with matrix functions E,A∈C(I,Kn×n),K∈ {C,R}.

Analyzing the different stability concepts for (19) is, however, much more complicated than for linear

time-invariant systems. Even if for all t∈I, the finite eigenvalues of (E(t),A(t)) have negative real part,

system (19) may be unstable, as many well-known examples demonstrate for the ODE case, see e.g., [51].

Example 3. For all t∈R

A(t) = cos2(3t)−5sin2(3t)−6cos(3t)sin(3t)+3

−6cos(3t)sin(3t)+ 3 sin2(3t)−5cos2(3t)

has a double eigenvalue at −2 but the solution of ˙x=Ax, with x(0) = c1

0is given by x(t) =

c1etcos(3t)

−c1etsin(3t), which is obviously unstable.

Robust stability of differential-algebraic equations 11

In the following we assume that (19) is of index at most one, in the sense of the tractability index [37]

or that it is strangeness-free in the sense of the strangeness-index [56]. Let us briefly introduce these index

concepts.

We use the following standard notation as in [31, 54]. Let X,Ybe finite dimensional vector spaces. For

every p,1≤p<∞and s,t,t0≤s<t<∞, we denote by Lp(s,t;X)the space of measurable functions

fwith values in Xand norm kfkp:=Rt

skf(ρ)kpdρ1/p<∞and by L∞(s,t;X)the space of measur-

able and essentially bounded functions fwith kfk∞:=ess sup ρ∈[s,t]kf(ρ)k. We also consider the spaces

Lloc

p(t0,∞;X)and Lloc

∞(t0,∞;X), which contain all functions f∈Lp(s,t;X)and f∈L∞(s,t;X)for some

s,t,t0≤s<t<∞, respectively. We, furthermore, use the notation L(Lp(t0,∞;X),Lp(t0,∞;Y)) to denote

the Banach space of linear bounded operators Pfrom Lp(t0,∞;X)to Lp(t0,∞;Y)supplied with the norm

kPk:=supx∈Lp(t0,∞;X),kxk=1kPxkLp(t0,∞;Y).

To introduce the tractability index and the projector chain approach, see e.g., [37, 60], if E,A∈

Lloc

∞(0,∞;Kn×n), then for S:=kerEthere exists an absolutely continuous projector Qonto S, i.e.,

Q∈C(0,∞;Kn×n),Qis differentiable almost everywhere, Q2=Q, and Im Q=Sfor all t∈I. If we

assume in addition that ˙

Q∈Lloc

∞(0,∞;Kn×n), then P=I−Qis a projector along Sand system (19) can be

rewritten in the form

dt (Px) = b

Ax,(20)

where b

A:=A+E˙

P∈Lloc

∞(0,∞;Kn×n). Setting G:=E−b

AQ, the DAE (2.1) is said to be tractable of index

one, if G(t)is invertible for almost every t∈[0,∞)and G−1∈Lloc

∞(0,∞;Kn×n).

Multiplying both sides of (20) by PG−1,QG−1,we obtain

dt (Px) = ( d

dt P+PG−1b

A)Px,

Qx =QG−1b

APx,

which decomposes the DAE into a differential part and an algebraic part. With z=Px, the dynamics of the

system is given by the inherent ODE

˙z= ( ˙

P+PG−1b

A)z(21)

of (19). Let Φ0(t,s)denote the fundamental solution associated with the inherent ODE (21), i.e., the matrix

function satisfying

dt Φ0(t,s) = ( ˙

P+PG−1b

A)Φ0(t,s),Φ0(s,s) = I,

for t>s≥0, then the fundamental solution operator generated by (19) is defined by

dt Φ(t,s) = AΦ(t,s),P(s)(Φ(s,s)−I) = 0

and Φcan be expressed as Φ(t,s) = (I+QG−1b

A)(t)Φ0(t,s)P(s).

The concept of the strangeness index uses the DAE and its derivatives to construct a so-called

strangeness-free DAE with the same solution [56]. In the homogenous case this strangeness-free system

has the form

E˙x=Ax,t∈I,(22)

where E=E1

0,A=A1

A2,with E1∈C(I,Kd×n)and A2∈C(I,K(n−d)×n)such that the matrix ˆ

E(t):=

E1(t)

A2(t)is invertible for all t∈I.

By using a global kinematic equivalence transformation, see [62, Remark 13], (22) can be transformed

to the special form, ˜

E11 ˜

E12

0 0 ˙

˜x=˜

A11 ˜

A12

0˜

A22 ˜x,t∈I.(23)

12 Nguyen Huu Du and Vu Hoang Linh and Volker Mehrmann

so that the underlying ODE is given by E11 ˙x1=A11x1with nonsingular E11.

A matrix function Φ∈C1(I,Rn×d)is called minimal fundamental solution matrix of the strangeness-

free DAE (22) if each of its columns is a solution to (22) and rankΦ(t) = d, for all t∈I.

In the following we assume that the DAE is of index at most one or alternatively strangeness-free. These

conditions are equivalent if the coefficients are sufficiently smooth [69].

The characterization of the different stability concepts for linear variable coefficient ODEs is well es-

tablished via the concepts of Bohl and Lyapunov exponents, [3, 22] and Sacker-Sell spectra [23, 78]. These

concepts were extended from ODEs to DAEs in [18, 62, 63]. In the following, for simplicity, we assume

that (19) possesses piecewise continuous coefficients.

To analyze exponential stability, we introduce first the Bohl exponent.

Definition 6. The Bohl exponent for an index one system of the form (19) with fundamental solution Φis

given by

kB(E,A) = inf{−α∈R; there exists Lα>0

such that for all t≥t0≥0 : kΦ(t,t0)k ≤ Lαe−α(t−t0)o.

It follows that (19) is exponentially stable if and only if its largest Bohl exponent is negative.

Analogously to the ODE case (see [22]), using the fundamental solution operator Φ, it follows that the

Bohl exponent of (19) is finite if and only if sup0≤|t−s|≤1kΦ(t,s)k<∞. Furthermore, if the Bohl exponent

of (19) is finite, then it can be determined by

kB(E,A) = limsup

s,t−s→∞

lnkΦ(t,s)k

t−s.

In [18], various properties of the Bohl exponent, as well as the connection between the exponential stabil-

ity of (19) and the boundedness of solutions to nonhomogeneous DAE with bounded inhomogeneity are

investigated.

For ODEs, the asymptotic stability of solutions can be characterized by the Lyapunov exponents, see

[66]. The extension of the theory of Lyapunov exponents to linear time-varying DAEs has been given

in [19, 20, 62, 63, 64], using either the projector-based tractability index or the derivative array-based

strangeness index approach.

Definition 7. For a given minimal fundamental solution matrix Φof a strangeness-free DAE system of the

form (22), and for 1 ≤i≤d, we introduce

λu

i=limsup

t→∞

tln||Φ(t)ei|| and λ`

i=liminf

t→∞

tln||Φ(t)ei||,

where eidenotes the i-th unit vector and ||·|| denotes the Euclidean norm. The columns of a minimal fun-

damental solution matrix form a normal basis if Σd

i=1λu

iis minimal. The λu

i,i=1,2,...,dbelonging to

a normal basis are called (upper) Lyapunov exponents and the intervals [λ`

i,λu

i],i=1,2,...,d, are called

Lyapunov spectral intervals.

The strangeness-free DAE system (22) then is asymptotically stable if and only if the largest upper Lya-

punov exponent is negative. Note that for linear time-invariant DAEs (9), the Lyapunov exponents are

exactly the real parts of the finite eigenvalues of pencil (E,A).

This brings us to another major difference between linear time-invariant and linear time-varying DAEs.

The following example shows is that an infinitesimally small time-varying perturbation applied to both

coefficient matrices may change the asymptotic stability, even if the perturbation does not change the

index.

Example 4. [63] Consider the system ˙x1=−x1, 0 =x2, which is strangeness-free and asymptotically stable.

For the perturbed DAE

(1+ε2sin(2nt)) ˙x1−εcos(nt)˙x2=−x1,0=−2εsin(nt)x1+x2,(24)

Robust stability of differential-algebraic equations 13

where εis a small parameter and nis a given number, from the second equation of (24), we obtain x2=

2εsinnt x1. Differentiating this expression for x2and inserting the result into the first equation, after some

elementary calculations, we obtain ˙x1= (−1+nε2+nε2cos(2nt))x1.Explicit integration yields x1(t) =

e(−1+nε2)t+ε2sin(2nt)/2x1(0). Clearly, even if εis arbitrarily small (hence the perturbation in the coefficient

matrices is arbitrarily small in the sup-norm), (24) may become unstable if nis sufficiently large.

3.1 Stability radii for linear time-varying DAEs

In this section we discuss the stability radii for linear DAEs with variable coefficients. Formulas for the

stability radii of exponential stability were derived in [18, 31] extending the results for ODEs in [43, 54].

We assume that the DAE is of index at most one and discuss perturbed systems

E˙x= (A+H)x,t∈I,(25)

with a perturbation function H:I→Kn×nas well as structured perturbations of the form

E˙x=Ax +B∆(Cx),t∈I,(26)

where B:I→Kn×mand C∈C:I→Kq×nare given matrix functions, restricting the structure of the

perturbation and ∆:Lp(0,∞;Km)→Lp(0,∞;Kq)is an unknown perturbation operator. The exponential

stability radius is then defined to be the largest bound Rsuch that the exponential stability is preserved for

all perturbations Hor ∆, respectively, of norm strictly less than R.

To obtain formulas for the exponential stability radius, we assume that (19) is exponentially stable, i.e.,

there exist constants L>0, γ>0 such that kΦ0(t,s)P(s)k≤Le−γ(t−s), for t≥s≥0, and that it is robustly

index one in the sense that, supplied with a bounded projection Q, the matrix functions G−1and −QG−1b

are essentially bounded on I.

To extend the tractability index concept to the perturbed system (26), we assume that the perturbation

operator ∆∈L(Lp(0,∞;Kq),Lp(0,∞;Km)) is causal which is defined as follows. For T∈I, the truncation

operator πTat Ton Lp(0,∞;X)is defined via

πT(u):=u,t∈[0,T],

0,t>T.

An operator P∈L(Lp(0,∞;X),Lp(0,∞;Y)) then is said to be causal, if πTPπT=πTPfor every T∈I.

Let the linear operator e

G∈L(Lloc

p(0,∞;Kn),Lloc

p(0,∞;Kn)) be defined via

Gz) = (E−b

AQ)z−B∆(CQz),t∈I.

If for every T>0, the operator e

Grestricted to Lp(0,T;Kn)is invertible and the inverse operator e

G−1is

bounded, then we say the perturbed DAE (26) is of index one (in a generalized sense). Then we can employ

the concept of mild solution.

Definition 8. We say that the initial value problem for the perturbed system (26) with initial condition

P(t0)(x(t0)−x0) = 0,(27)

admits a mild solution if there exists x∈Lloc

p(t0,∞;Kn)satisfying

x(t) = Φ(t,t0)P(t0)x0+Zt

Φ(t,ρ)PG−1B∆([Cx(·)]t0)(ρ)dρ+QG−1B∆([Cx(·)]t0)(t)

for t≥t0,where

[Cx(·)]t0=0,t∈[0,t0),

C(t)x(t),t∈[t0,∞).

14 Nguyen Huu Du and Vu Hoang Linh and Volker Mehrmann

The existence and uniqueness of mild solutions is given by the following result.

Theorem 5. [31] Consider the initial value problem (26)–(27). If (26) is of index at most one, then it

admits a unique mild solution x ∈Lloc

p(t0,∞;Kn)with absolutely continuous z =Px for all t0∈I,x0∈Kn.

Furthermore, for an arbitrary T >0, there exists a constant L1such that pointwise

kP(t)x(t)k≤L1kP(t0)x0kfor all t ∈[t0,T].

With (26), we associate the perturbation operator

Lt0z=CZt

Φ(t,ρ)PG−1B(ρ)z(ρ)dρ+CQG−1Bz,(28)

which is defined for all t≥t0≥0,z∈Lp(0,∞;Km). For exponentially stable robust index one systems this

operator is linear, bounded, and monotonically non-increasing with respect to t, i.e.,

kLt0k≥kLt1kfor all t1≥t0≥0,

and for all t∈Ithe bound

kLtk ≤ L

γPG−1∞||B||∞||C||∞+CQG−1B∞.

holds.

Denoting by x(t;t0,x0)the solution of (26) with initial condition x(t0) = x0, the trivial solution of (26)

is said to be globally Lp−stable, if there exist a constant L2>0 such that

kx(·;t0,x0)kLp(t0,∞;Kn)≤L2kP(t0)x0kKn,(29)

for all t≥t0,x0∈Kn, see [31]. Note that this kind of stability notion is equivalent to the concept of output

stability, see [53] for various stability concepts in the ODE case.

We then have the following definition of stability radii for time-varying DAEs.

Definition 9. If system (19) is exponentially stable and robustly index one, then the complex/real structured

stability radii of (19) subject to linear, dynamic and causal perturbations as in (26), are defined by

rK(E,A;B,C) = infk∆k,the trivial solution of (26) is not globally Lp−stable

or (26) is not of index one,

where K=Cor K=R, respectively.

In [31] the following formulas for the stability radii were derived.

Theorem 6. [31] If system (19) is exponentially stable and robustly index one, then

rK(E,A;B,C) = min{lim

t0→∞

Lt0

−1,(ess sup t∈I

CQG−1B(t)

)−1}.

This implies the following corollary.

Corollary 1. If system (19) is exponentially stable and robustly index one and the data (E,A;B,C)are real,

then

rC(E,A;B,C) = rR(E,A;B,C).

As special case we obtain the formula (12) for the complex stability radius of time-invariant systems.

Corollary 2. Let E,A,B,C be time-invariant, let the system (19) be index one and exponentially stable. If

p=2, i.e., for the space L2of square integrable functions

rC(E,A;B,C) = kL0k−1=sup

ω∈iR

C(ωE−A)−1B

−1

Robust stability of differential-algebraic equations 15

Note that the formula obtained in Corollary 2 is just a special case of that in Theorem 1 when only Ais

subject to perturbations and in this case the complex stability radius with respect to dynamic perturbations

and that to static perturbations are equal. This results generalizes that for the ODE case in [49].

Example 5. [31, Sect. 5.1] Consider a DAE in semi-explicit form

E=Ir0

0 0 ,A=A11 A12

A21 A22 (30)

with analogous partitioning. The index one assumption means that A22(t)is invertible almost everywhere

in I. One gets Q=diag(0,In−r),

G=Ir−A12

0−A22 ,Φ(t,s) = b

Φ(t,s)

−A−1

22 A21 b

Φ(t,s),

where b

Φ(t,s)is the fundamental solution operator of the underlying ordinary differential equation ˙y=

(A11 −A12A−1

22 A21)y, which is assumed to be exponentially stable. The assumption that the system is ro-

bustly index means the essential boundedness of A−1

22 ,A−1

22 A21 and A12A−1

22 . Partitioning the restriction ma-

trices B,Cas B=BT

1BT

2T,C= [C1C2], analogously, we obtain

(Lt0u)(t) = (C1−C2A−1

22 A21)Zt

t0b

Φ(t,ρ)(B1(ρ)−A12A−1

22 B2(ρ))u(ρ)dρ−C2A−1

22 B2u(t),

and by Theorem 6 we have

rK(E,A;B,C) = min{lim

t0→∞

Lt0

−1,(ess sup t∈IkC2A−1

22 B2(t)k)−1}.

In summary, we have seen in this section that in the index one case the exponential stability radii in the

linear time-varying case are similar to the ones in the linear time-invariant case.

A similar analysis can be performed in principle for the asymptotic stability, by studying the Lyapunov

exponents under perturbations, see the next section. Unfortunately, however, in contrast to the Bohl expo-

nents, the Lyapunov exponents themselves may be very sensitive to small changes in the data.

3.2 Dependence of stability radii on the data

The robustness analysis of the Bohl exponent was extended from ODEs in [43] to DAEs in [18]. In the

context of deriving numerical methods for the numerical computation of Bohl and Lyapunov exponents

for linear time-varying DAEs, see [62, 63, 65] the robustness of these exponents under perturbations was

studied and the concept of admissible perturbations was to extended to the variable coefficient case. In this

section we summarize these results and study how the stability radii depend on the data.

To see that the Bohl exponent is robust under sufficiently small index one preserving perturbations, we

consider the perturbed equation

dt (Px) = (b

A+H)x(31)

and assume that

sup

t∈I

kH(t)k<(sup

t∈I

kQG−1(t)k)−1.(32)

Note that this condition implies the inequality supt∈IkQG−1H(t)k<1, which is essential in the analysis.

Theorem 7. [18] Suppose that system (19) is robustly index one and satisfies (32) and suppose further that

for any ε>0there exists δ=δ(ε)>0such that

16 Nguyen Huu Du and Vu Hoang Linh and Volker Mehrmann

limsup

s,t−s→∞

t−sZt

kPG−1H(τ)kdτ<δ,

then

kB(E,A+H)<kB(E,A)+ε.

As a consequence of Theorem 7, the exponential stability is preserved under all sufficiently small pertur-

bations H.

Remark 7. The robustness analysis of the Bohl exponent is extendable to the case of general perturbations

arising in both coefficients of (19). It is clear that additional assumptions on the perturbation structure

and/or the smoothness of the admissible perturbations are necessary in this case. The same can be said for

the analysis of Bohl exponents for general higher-index DAEs and for nonlinear perturbations, see [2].

We now discuss how the structured stability radius of (25) depends on the perturbation Hand the restric-

tion matrices B,C. To this end, we first establish the asymptotic behavior of the norm of the input-output

operator defined in (28).

Theorem 8. [18] Suppose that system (19 is exponentially stable, robustly index one and satisfies (32) and

suppose, in addition, that the perturbation function H satisfies

lim

t→∞kH(t)k=0.

Then the operator Ltdefined in (28) satisfies

lim

t→∞kLtk=lim

t→∞ke

Ltk.

By invoking Theorem 6, we obtain sufficient conditions for a sequence of perturbations Hkunder which

the structured stability radius of the perturbed systems converges to that of the unperturbed system.

Theorem 9. [18] Suppose that system (19) is exponentially stable, robustly index one and satisfies

(32). Let {Hk(·)}∞

k=1be a sequence of measurable matrix functions and suppose that supt∈IkHk(t)k<

(supt∈IkQG−1(t)k)−1for all k =1,2,... and lim

k→∞supt∈IkHk(t)k=0. Then

lim

k→∞rK(E,A+Hk;B,C) = rK(E,A;B,C).

Theorem 9 implies that the stability radius for the system (19) depends continuously on the entries of the

coefficient matrix function Aand as a consequence of Theorem 2, we get the following corollary.

Corollary 3. Let E,A,B,C be constant matrices, let system (19) be of index one and exponentially stable.

Furthermore, assume that the sequence of time-varying perturbation {Hk}∞

k=1fulfills the conditions of

Theorem 9. Then, for the Euclidean norm, one has

lim

k→∞rC(E,A+Hk;B,C) = ( sup

w∈iR

kC(wE −A)−1Bk)−1.

To illustrate this result, in [18, Example 5.13], a numerical example is given which shows that the complex

stability radius of a linear time-invariant system under time-varying perturbations is computable by that for

the corresponding time-invariant DAE.

Finally, one also obtains that the structured stability radius of (19) depends continuously on the restric-

tion matrices Band Cas in the time-invariant case.

Theorem 10. [18] Suppose that system (19) is exponentially stable, robustly index one and satisfies (32).

Let Bkand Ckbe two sequences of measurable and essentially bounded matrix functions satisfying

lim

k→∞ess sup t∈I||Bk(t)−B(t)|| =0,lim

k→∞ess sup t∈I||Ck(t)−C(t)|| =0

then,

lim

k→∞rK(E,A;Bk,Ck) = rK(E,A;B,C).

Robust stability of differential-algebraic equations 17

We stress once more that, since the dynamics of DAEs is constrained and the index one property should

be preserved, only weaker results hold for the continuity of the stability radius and more restrictive as-

sumptions are required than in the ODE case. Furthermore, for time-varying DAEs, we have to restrict the

analysis to perturbations in Aonly, see (25) and (26), simply because the study of perturbations associated

with the leading term Eis still mainly an open problem.

In order to study the robustness of Lyapunov exponents, we consider the specially perturbed system

[E+F]˙x= [A+H]x,t∈I,(33)

where F=F1

0and H=H1

H2, and where F1and H1,H2are assumed to have the same order of smooth-

ness as E1and A1,A2, respectively. Perturbations of this structure are called admissible, generalizing the

concept for the constant coefficient DAEs studied in [15].

The DAE (22) is said to be robustly strangeness-free if it stays strangeness-free under all sufficiently

small admissible perturbations and it is easy to see that this property holds if and only if the matrix function

Eis boundedly invertible.

If we assume that (22) is already given in the form (23), then the perturbed DAE has the form

(E11 +F11)d

dt ˜x1+(E12 +F12)d

dt ˜x2= (A11 +H11)˜x1+ (A12 +H12)˜x2

0=H21 ˜x1+(A22 +H22)˜x2.

In the following we restrict ourselves to robustly strangeness-free DAE systems under admissible per-

turbations.

Definition 10. The upper Lyapunov exponents λu

1≥... ≥λu

dof (22) are said to be stable if for any ε>0,

there exists δ>0 such that the conditions supt||F(t)|| <δ,supt||H(t)|| <δ, and supt˙

H2(t)<δon the

perturbations imply that the perturbed DAE system (33) is strangeness-free and

|λu

i−γu

i|<ε,for all i=1,2,...,d,

where the γu

iare the ordered upper Lyapunov exponents of (33).

The boundedness condition on ˙

H2, which is obviously satisfied in the time-invariant setting [15], is an extra

condition and it seems to be somehow unusual. However, the DAE (24) shows the necessity.

As in the ODE case, see [1, 23], to have stability of the Lyapunov spectrum, one needs the property of

integral separation, i.e., for the columns of the minimal fundamental solution matrix Φof (22) there exist

constants c1>0 and c2>0 such that

||X(t)ei||

||X(s)ei|| ·||X(s)ei+1||

||X(t)ei+1|| ≥c2ec1(t−s),

for all t,swith t≥s≥0 and i=1,...,d−1. Then we have the following sufficient conditions for the

stability of the upper Lyapunov exponents of (22).

Theorem 11. [62] Consider the DAE (22) in the form (23). Suppose that the matrix ˆ

E is boundedly in-

vertible and that E−1

11 A11, A12A−1

22 and the derivative of A22 are bounded on I. Then, the upper Lyapunov

exponents of (22) are distinct and stable if and only if the system has the integral separation property.

This shows that if perturbations are performed in Eas well, then the perturbation analysis of time-varying

DAEs requires more restrictive conditions than in the time-invariant case. However, for some classes of

structured problems and/or structured perturbation, parts of these conditions can be relaxed.

If the perturbation block H21 disappears, i.e., if Hand Ahave the same block triangular structure, then

for example the restrictive conditions on the derivatives in Definition 10 and Theorem 11 can be omitted.

Similar situations happen in the case that E12 =F12 =0 as discussed in [62, Section 3.2] and the case of

perturbations in Aonly as in (25).

18 Nguyen Huu Du and Vu Hoang Linh and Volker Mehrmann

In [62] another stability concept, the Sacker-Sell spectrum has been extended to linear DAEs with vari-

able coefficients. It is shown also shown that unlike the Lyapunov spectrum, the Sacker-Sell spectrum is

robust in the sense that it is stable without requiring integral separation.

This means that for general strangeness-free time-varying systems, the exponential stability is robust,

but the asymptotic stability is not. However, this remark does not apply to time-invariant systems, for which

the two stability notions are equivalent.

Remark 8. The robustness analysis for linear DAEs of index higher than one and under general perturba-

tions is essentially an open problem. The same is true for the distances to the other important control

properties such as controllability and observability. The robustness of these concepts for linear DAEs with

variable coefficients presents a major challenge.

4 Discussion

In this paper we have surveyed recent results on the robustness of asymptotic and exponential stability for

linear time-invariant and time-varying DAEs. We have analyzed robust stability and its distance measures,

the real or complex structured stability radius and presented formulas and various properties of the stability

radii. We have seen that the robustness analysis for DAEs is much more complicated than that for ODEs. In

general, results already known for ODEs now hold for DAEs only under extra assumptions, mainly restrict-

ing the set of admissible perturbations. DAE aspects also give rise to new robustness and distance problems.

While for time-invariant DAEs, most of the robustness and distance problems are well understood, many

problems for time-varying DAEs are still open. These and robustness analysis for general nonlinear and/or

high-index DAEs and time-delay DAEs are interesting and challenging topics for future work.

References

1. L. Ya. Adrianova. Introduction to linear systems of differential equations, Trans. Math. Monographs, Vol. 146, AMS,

Providence, RI, 1995.

2. T. Berger. Bohl exponent for time-varying linear differential-algebraic equations. Int. Journal of Control, 2012 (to

appear).

3. P. Bohl. ¨

Uber differentialungleichungen. J. f. d. Reine und Angew. Math., 144:284–313, 1913.

4. S. Boyd and V. Balakrishnan. A Regularity Result for the Singular Values of a Transfer Matrix and a Quadratically

Convergent Algorithm for Computing its L∞-norm. Systems Control Letters, 15:1-7, 1990.

5. M. Bracke. On stability radii of parametrized linear differential-algebraic systems. Ph.D. Thesis, University of Kaiser-

slautern, 2000.

6. N.A. Bruinsma and M. Steinbuch. A fast algorithm to compute the H∞of a transfer function matrix. Systems Control

Letters, 14:287–293, 1990.

7. K. E. Brenan, S. L. Campbell, and L. R. Petzold. Numerical Solution of Initial-Value Problems in Differential Algebraic

Equations. SIAM Publications, Philadelphia, PA, 2nd edition, 1996.

8. A. Bunse–Gerstner, R. Byers, V. Mehrmann and N.K. Nichols, Feedback Design for Regularizing Descriptor Systems.

Linear Algebra Appl. 299:119–151, 1999.

9. A. Bunse–Gerstner, V. Mehrmann and N.K. Nichols, Regularization of descriptor dystems by derivative and proportional

state feedback. SIAM J. Matrix Anal. Appl., 13:46–67, 1992.

10. A. Bunse–Gerstner, V. Mehrmann and N.K. Nichols, Regularization of descriptor dystems by output feedback. IEEE

Trans. Automatic Control, 39:1742–1748, 1994.

11. J. Burke, A.S. Lewis, and M.L. Overton. Optimization and pseudospectra, with applications to robust stability. SIAM J.

Matrix Anal. Appl., 25:80–104, 2003.

12. R. Byers. A Bisection Method for Measuring the Distance of a Stable Matrix to the Unstable Matrices. SIAM J. Sci. and

Stat. Comput., 9:875–881, 1988.

13. R. Byers. The descriptor controllability radius. Systems and Networks: Mathematical Theory and Application, Proceed-

ings of MTNS ’93, Akademie Verlag, Berlin, 85–88, 1994.

14. R. Byers, C. He, and V. Mehrmann. Where is the nearest non-regular pencil? Lin. Alg. Appl., 285:81–105, 1998.

15. R. Byers and N.K. Nichols. On the stability radius of a generalized state-space system. Lin. Alg. Appl., 188-189:113–134,

1993.

16. S. L. Campbell. Linearization of DAE’s along trajectories. Z. Angew. Math. Phys., 46:70–84, 1995.

Robust stability of differential-algebraic equations 19

17. S. L. Campbell, N. K. Nichols, and W. J. Terrell. Duality, observability, and controllability for linear time-varying

descriptor systems. Circ. Syst. Signal Process., 10:455–470, 1991.

18. C.J. Chyan, N.H. Du, and V.H. Linh. On data-dependence of exponential stability and the stability radii for linear time-

varying differential-algebraic systems. J. Diff. Equations, 245:2078–2102, 2008.

19. N.D. Cong and H. Nam. Lyapunov’s inequality for linear differential algebraic equation. Acta Math. Vietnam, 28:73–88,

2003.

20. N.D. Cong and H. Nam. Lyapunov regularity of linear differential algebraic equations of index 1. Acta Math. Vietnam,

29:1–21, 2004.

21. L. Dai. Singular Control Systems. Lecture Notes in Control and Information Science, Springer-Verlag, 1989.

22. J.L. Daleckii and M.G. Krein. Stability of solutions of differential equations in Banach spaces. American Mathematical

Society, Providence, RI, 1974.

23. L. Dieci and E. S. Van Vleck. Lyapunov spectral intervals: theory and computation, SIAM J. Numer. Anal., 40:516–542,

2002.

24. M. Dmitriev and G. Kurina. Singular perturbation in control problems. Automation and Remote Control, 67:1–43, 2006.

25. V. Dragan. The asymptotic behavior of the stability radius for a singularly perturbed linear system. Int. J. Robust and

Nonlin. Control, 8:817–829, 1998.

26. V. Dragan and A. Halanay Stabilization of Linear Systems. Birkh¨

auser, Boston, 1999.

27. N.H. Du. Stability radii of differential-algebraic equations with structured perturbations. Systems Control Letters, 57:546–

553, 2008.

28. N.H. Du, D.T. Lien, and V.H. Linh. Complex stability radii for implicit discrete-time systems. Vietnam J. Math., 31:475–

488, 2003.

29. N.H. Du and V.H. Linh Implicit-system approach to the robust stability for a class of singularly perturbed linear systems.

Systems Control Letters, 54:33–41, 2005.

30. N.H. Du and V.H. Linh. Robust stability of implicit linear systems containing a small parameter in the leading term. IMA

J. Math. Cont. Inf., 23:67–84, 2006.

31. N.H. Du and V.H. Linh. Stability radii for linear time-varying differential-algebraic equations with respect to dynamic

perturbations. J. Diff. Equations, 230:579–599, 2006.

32. N.H. Du, D.D. Thuan, and N.C. Liem. Stability radius of implicit dynamic equations with constant coefficients on time

scales. Systems Control Letters, 60:596–603, 2011.

33. E. Eich-Soellner and C. F¨

uhrer. Numerical Methods in Multibody Systems. Teubner Verlag, Stuttgart, Germany, 1998.

34. M.A. Freitag and A. Spence. A Newton-based method for the calculation of the distance to instability. Linear Alg. Appl.,

435(12): 3189-3205, 2011.

35. F.R. Gantmacher, Theory of Matrices. Vol. 2, Chelsea, New York, 1959.

36. G. H. Golub and C. F. Van Loan. Matrix Computations. The Johns Hopkins University Press, Baltimore, MD, 3rd edition,

1996.

37. E. Griepentrog and R. M¨

arz. Differential-Algebraic Equations and their Numerical Treatment. Teubner Verlag, Leipzig,

Germany, 1986.

38. M. Gu, E. Mengi, M. L. Overton, J. Xia and J. Zhu, Fast Methods for Estimating the Distance to Uncontrollability SIAM

J. Matrix Anal. Appl. 28:477–502, 2006.

39. M. G¨

urb¨

uzbalaban and M. L. Overton, Some Regularity Results for the Pseudospectral Abscissa and Pseudospectral

Radius of a Matrix. SIAM J. Optimization 22:281–285, 2012.

40. N. Guglielmi and M. L. Overton, Fast Algorithms for the Approximation of the Pseudospectral Abscissa and Pseudospec-

tral Radius of a Matrix SIAM J. Matrix Anal. Appl. 32:1166–1192, 2011.

41. E. Hairer and G. Wanner. Solving Ordinary Differential Equations II: Stiff and Differential-Algebraic Problems. Springer-

Verlag, Berlin, Germany, 2nd edition, 1996.

42. C. He, G.A. Watson. An algorithm for computing the distance to instability. SIAM J. Matrix Anal. Appl., 20:101–116,

1999.

43. D. Hinrichsen, A. Ilchmann, A.J. Pritchard. Robustness of stability of time-varying linear systems. J. Differential

Equations, 82:219–250, 1989.

44. D. Hinrichsen and B. Kelb. Spectral value sets: a graphical tool for robustness analysis. Systems Control Letters, 21(2),

127–136, 1993.

45. D. Hinrichsen, B. Kelb, and A. Linnemann. An algorithm for the computation of the complex stability radius. Automatica,

25:771–775, 1989.

46. D. Hinrichsen and A.J. Pritchard, Stability radii of linear systems. Systems Control Letters, 7:1–10, 1986.

47. D. Hinrichsen and A.J. Pritchard, Stability radii for structured perturbations and the algebraic Riccati equations. Systems

Control Letters, 8:105–113, 1986.

48. D. Hinrichsen and A.J. Pritchard A note on some difference between real and complex stability radii. Systems Control

Letters, 14:401–408, 1990.

49. D. Hinrichsen, A.J. Pritchard. Destabilization by output feedback. Differential Integral Equations, 5:357–386, 1992.

50. D. Hinrichsen and A.J. Pritchard. On spectral variations under bounded real matrix perturbations Numer. Math., 60:509–

524, 1992.

51. D. Hinrichsen and A. J. Pritchard. Mathematical Systems Theory I. Modelling, State Space Analysis, Stability and

Robustness. Springer-Verlag, New York, NY, 2005.

52. G. Hu. Robustness measures for linear time-invariant time-delay systems. Ph.D. Thesis, University of Toronto, 2001.

20 Nguyen Huu Du and Vu Hoang Linh and Volker Mehrmann

53. A. Ilchmann, I.M.Y. Mareels. On stability radii of slowly time-varying systems. In: Advances in mathematical system

theory, Boston, Birkh¨

auser, 55–75, 2001.

54. B. Jacob. A formula for the stability radius of time-varying systems. J. Differential Equations, 142:167–187, 1998.

55. P. Kokotovic, H.K Khalil, and J. O’Reilly Singular Perturbation Method in Control: Analysis and Design. Academic,

New York, 1986.

56. P. Kunkel and V. Mehrmann. Differential-Algebraic Equations. Analysis and Numerical Solution. EMS Publishing

House, Z¨

urich, Switzerland, 2006.

57. P. Kunkel, V. Mehrmann, W. Rath, and J. Weickert. A new software package for linear differential–algebraic equations.

SIAM J. Sci. Comput., 18:115–138, 1997.

58. S. Lam. Real robustness radii and performance limitations of LTI control systems. Ph.D. Thesis, University of Toronto,

2011.

59. S. Lam and E.J. Davison. Real controllability radius of high-order, descriptor, and time-delay LTI systems. in: Proceed-

ings of the 18th IFAC World Congress, Milano (Italy), August 28 - September 2, 2011, 6 pages, 2011.

60. R. Lamour, R. M¨

arz, and C. Tischendorf. Differential-Algebraic Equations: A Projector Based Analysis. Springer, 2013

(to appear).

61. L. Lee, C-H. Fang, J-G. Hsieh, Exact unidirectional perturbation bounds for robustness of uncertain generalized state-

space systems: continuous-time cases. Automatica J. IFAC, 33:1923–1927, 1997.

62. V.H. Linh and V. Mehrmann. Lyapunov, Bohl and Sacker-Sell spectral intervals for differential-algebraic equations. J.

Dynamics Diff. Eq., 21:153–194, 2009.

63. V.H. Linh and V. Mehrmann. Approximation of spectral intervals and associated leading directions for linear differential-

algebraic systems via smooth singular value decompositions. SIAM J. Numer. Anal., 49:1810–1835, 2011.

64. V.H. Linh and V. Mehrmann. Spectral analysis for linear differential-algebraic systems. In Proceedings of the 8th AIMS

International Conference on Dynamical Systems, Differential Equations and Applications, Dresden, May 24-28, 2010,

DCDS Supplement, 991–1000, 2011.

65. V.H. Linh, V. Mehrmann, and E. Van Vleck. QR methods and error analysis for computing lyapunov and sacker-sell

spectral intervals for linear differential-algebraic equations. Adv. Comput. Math., 35:281-322, 2011.

66. A. M. Lyapunov. The general problem of the stability of motion. Translated by A. T. Fuller from E. Davaux’s French

translation (1907) of the 1892 Russian original. Internat. J. Control, pages 521–790, 1992.

67. R.M.M. Mattheij and P.M.E.J. Wijckmans. Sensitivity of solutions of linear DAE to perturbations of the system matrices.

Numer. Alg., 19:159–171, 1998.

68. V. Mehrmann, The Autonomous Linear Quadratic Control Problem: Theory and Numerical Solution. LECTURE NOTES

IN CONTROL AND INFORMATION SCIENCES, Vol. 163, Springer Verlag, Heidelberg, Juli 1991.

69. V. Mehrmann, Index concepts for differential-algebraic equations. PREPRINT 2012–03, Institut f¨

ur Mathematik, TU

Berlin, 2012. url: http://www.math.tu-berlin.de/preprints/.

70. V. Mehrmann and T. Stykel, Descriptor systems: A general mathematical framework for modelling, simulation and

control. at-Automatisierungstechnik, 54(8):405–415, 2006.

71. E. Mengi. Measures for robust stability and controllability. Ph.D. Thesis, University of New York, 2006.

72. E. Mengi and M. L. Overton, Algorithms for the Computation of the Pseudospectral Radius and the Numerical Radius

of a Matrix IMA J. Numer. Analysis, 25: 48–669, 2005.

73. M. Motscha. Algorithm to compute the complex stability radius. Int. Journal of Control, 48:2417–2428, 1988.

74. L. Qiu, B. Bernhardsson, A. Rantzer, E.J. Davison, P.M. Young, and J.C. Doyle A formula for computation of the real

stability radius. Automatica, 31:879–890, 1995.

75. L. Qiu, and E.J. Davison The stability robustness of generalized eigenvalues. IEEE Transactions on Automatic Control,

37:886–891, 1992.

76. P. J. Rabier and W. C. Rheinboldt. Theoretical and Numerical Analysis of Differential-Algebraic Equations, volume VIII

of Handbook of Num. Analysis. Elsevier Publications, Amsterdam, The Netherlands, 2002.

77. R. Riaza. Differential-algebraic systems. Analytical aspects and circuit applications. World Scientific Publishing Co.

Pte. Ltd., Hackensack, NJ., 2008.

78. R.J. Sacker and G.R. Sell. A spectral theory for linear differential systems. J. Diff. Equations, 27:320–358, 1978.

79. N.K. Son and D. Hinrichsen. Robust stability of positive continuous time systems. Numer. Funct. Anal. Optim., 17:649–

659, 1996.

80. N.K. Son and D.D. Thuan. The structured distance to non-surjectivity and its application to calculating the controllability

radius of descriptor systems. J. Math. Anal. Applic., 388:272–281, 2012.

81. J. Sreedhar, P. Van Dooren, and A. Tits. A fast algorithm to compute the real structured stability radius. Int. Series of

Numer. Math., 121:219–230, 1996.

82. H. Tidefelt. Differential-algebraic equations and matrix-valued singular perturbation. Ph.D. Thesis, Linkping University,

2009.

83. L.N. Trefethen and M. Embree. Spectra and Pseudospectra: The Behavior of Nonnormal Matrices and Operators.

Princeton University Press, 2005.

84. C.F. Van Loan. How near is a stable matrix to an unstable matrix? Contemp. Math., 47:465–477, 1985.