A Criterion for Hurwitz Polynomials and its Applications

Published Online February 2011 in MECS

Let f (to) = a0to" + a1^-1 +... + an-to + an    (1)

is a polynomial of degree n with real coefficients, then it is said to be a stable or a Hurwitz polynomial if and only all its roots lie in the open left half of the complex plane. i.e. f ( to ) = 0 ^ Re to <  0 .

The well-known Routh-Hurwitz Theorem provides a powerful tool to check if a real polynomial is stable or not.

The Routh-Hurwitz Criterion. (See [1].) The polynomial (1) is stable if, and only if, the following inequalities holds: A 1 >  0, A 2 >  0,..., A _n >  0 , where

	a 1	a 0	0	0	0	0	... 0
	a 3	a 2	a 1	a 0	0	0	... 0
A k =	a 5	a 4	a 3	a 2	a 1	a 0	... 0	. (2)
	a 7	a 6	a 5	a 4	a 3	a 2	... 0
	M	M	M	M	M	M	OM
	a 2 k -1	a 2 k -2	a 2 k -3	a 2 k -4	a 2 k -5	a 2 k -6	...    ak

Supported by the National Science Foundation of Ningbo City (No.2009A610078), and the research foundation of Ningbo University (No.xkl09042, xkl09043).

usually, Д k ( k = 1,2,..., n ) named as the Hurwitz determinant. Whereas the inequalities A _k >  0 ( k = 1,2,..., n ) can be simplified when all of the coefficients of polynomial (1) are positive. The relative result was deduced by Li nard and Chipart in 1914.

The Li nard-Chipart Criterion. (See [2].)Consider a real polynomial f ( X ) = xⁿ + C 1 x^{n 1} + C 2 x^{n 2} + ... + C_n , it is stable if and only if the one of following cases holds:

(1)Cn > 0,Cn-2 > 0,...; A! > 0, A3 > 0,...

(2)Cn > 0,Cn-2 > 0,...; A2 > 0, A4 > 0,...

(3)Cn > 0,Cn-1 > 0,...; A, > 0, A3 > 0,...

(4)Cn > 0,Cn- > 0,...; A2 > 0, A4 > 0,...

Polynomial stability problems of various types arise in a number of problems in mathematics and engineering. We refer to [1, Chapter 15] for deep surveys on the classical stability theory and [3,4,5,6,7] of recent results.

In this paper, we present a new criterion for determining the Hurwitz polynomial.

The rest of this paper is organized in the following ways with the next section stating the main result. And also, we present some computation examples; deduce the explicit criterion for stability of real polynomials with degree no more than four, which is expressed by polynomials in the coefficient of the given polynomial.

In section III, two problems are discussed as the application of our criterion. The first problem is relative to the stability of the control system. In classical theory of control systems stability analysis, determining the stability of a given discrete time system can be converted to the stability of the corresponding continuous time system through a bilinear transformation. We introduced a necessary and sufficient condition to determine the stability of the characteristic polynomial of the continuous time control systems in subsection A. In subsection B, we discussed a pathological case of the bilinear transformation used widely, and brought forward an alternative one. The second problem, which is relative to the inertia of Bezout matrix, was discussed in subsection C. An algorithm named as INER was brought forward.

• II.    A C RITERION FOR H URWITZ POLYNOMIALS

A. Criterrion

Without loss of generality, supposed that the real polynomial f ( o ) = o ” + aO ” - 1 +... + a n is monic. We construct f ( o ) = f ( x + iy ) = f ( x , y ) + i • f 2 ( x , y ) , and get two real polynomials f l ( x , y ) , f ₂ ( x , y ) .

Let f l ( x , y ) and f 2 ( x , y ) be expressed respectively as follows:

f l ( x , y ) = « 0 y k + a y k - 1 +... + a ,           ⁽³⁾

f 2 ( x , y ) = P 0 y m + P l y m - 1 + ... + P m ,         (4)

Where k , m <  ” , k + m = 2 ” -1 and a i , P j are all the real polynomial with variable x . It is not difficult to verify the changes of a ₀, P ₀ according to ” .

Table 1 the changes of a ₀ , P ₀ according to ”

n       2	3                4	5
a 0        -1	- a ^ - 3 x      1	u , + 5 x
P 0       a y + 2 x      -1               - U y -		4 x      1
n       6	7                8
a 0        1	- U i - 7 x      1
в)       U i + 6 x      -1               - U i -		8 x     …
Then we compute the resultant of		f l ( x , y ) and
f 2 ( x , y ) with	respect to y , denoted by res ( f l , f 2, y )
marked (5).
	a 0 a 1 ••• a k
	a 0 a 1 • •• a k
	a 0 a , •• • a k
res ( fv f ^ y ) =		(5)
	P 0 P l ••• P m
	P 0 P l ••• P m
	P 0 P l ••• P m	(2 ” - 1) x (2 ” - 1)

The following lemma gives the relation between the real roots of res ( f , f ₂, y ) and the x - coordinates of the solutions of PS ( Polynomials System ) :

{fl (x, y ) = 0,f2 (x, y ) = 0} (6)

Lemma 1. (See [8].) In general, res ( f , f ₂, y ) is a polynomial in x whose roots are the x - coordinates of the solutions of PS . According to the classical result from Algebra, it is well known:

Lemma 2. (See [9].) res(f, f2,y) = 0if and only if fl (x, y) and f 2 (x, y) have a common factor in ° [ x, y ], lc(P,y) and lc(f2,y) by GLCx (fl, f2) . Obviously, GLCx (fl, f 2) is a polynomial in x, too.

Theorem 1. Denote the number of distinct real roots ^of res ( f l , f 2 , y ) and GLC x ( f, f 2 ) by ц and m , respectively, and the numbers of distinct x - coordinates of the solutions of PS by у , then у = ц - m . Especially, if GLC x ( f l , f , ) = l, then we have у = ц .

Proof .  It is obvious according to Lemma 1 and

Lemma 2.

Denote the set of distinct x - coordinates of the solutions of PS by Zero x ( PS ) , the distinct roots of res ( f l , f 2 , y ) by Zero ( res ( f l , f 2 , y ) ) , and the distinct roots o^f GLC x ( f , f 2 ) by Zero ( GLC x ( f l , f 2 ) ) .

In general, Zero x ( PS ) c Zero ( res ( f l , f ₂, y ) ) .From Theorem 1 it is easy to see that:

Corollary 1. Zerox (PS)uZero(GLCx (fl, f 2 )) = Zero (res (fl,f2, y)).

Denote the distinct real roots of res ( f l , f ₂, y ) and the distinct negative real roots of res ( f l , f ₂, y ) by l and m , respectively. Then we deduce the main result in following:

Theorem 2. f ( o ) e H [ f ] if and only if l = m ^ 0 .

Proof . Firstly, l = m ^ 0 implies all the real roots of res ( f l , f ₂, y ) are negative. Secondly, according to table l, GLC x ( f l , f ₂ ) = l - Combining corollary l, we deduce that Zero x ( PS ) = Zero ( res ( f l , f ₂, y ) ) . Thus, all of the roots of polynomial f ( o ) have negative real part. That is to say f ( o ) e H [ f ] .

• B. Explicit condition of stability for any real polynomial with degree no more than four

In this subsection, we will deduce the explicit condition of stability for any real polynomial with a degree no more than four. According to Theorem 2, it is important to analyze the number of distinct real roots and distinct negative real roots of a real polynomial. Let’s introduce the relative results from [10], [11] and [12].

Definition 1. Given a polynomial with real coefficients n f (x ) = £ aix”-i( a 0 ^ 0),          (7)

i = 0

The following ( 2 ” + l ) x ( 2 ” +1 ) matrix is called the discrimination matrix of and denoted by Discr ( f ) .

^or a o = в о = 0.

Let lc ( f , y ) and lc ( f ₂, y ) , which are all the polynomials with real coefficients in x , be the leading coefficient of the polynomial f l ( x , y ) and f ₂ ( x , y ) in y . And also we denote the greatest common factor of

^a 0 al    a 2   L an na0 (n -1)a1 L an1 a 0     al    L an-1 an

0     na о L 2 a n - 2 a n - i

a о    a l     a 2    ^L a n

0 naо (n -1)ai L an-i a0    al    a2 L an )

Definition 3. By { d 1 , d 2 ,..., d _{2 n + 1} } denote the principal minor sequence of Discr ( f ) , we call the 2 n - tuple { d 1 d ₂, d ₂ d 3 ,..., d ₂ n d ₂ n _{+ 1} } the negative root discrimination sequence, and denote it by n . _r . d ( f ) •

Definition 4. We call {sign (D1), sign (D2),..., sign (Dn )} the sign list of a given sequence {D1 (f), D2 (f),...Dn (f)}, where sign (•) is a sign function.

Definition 5. Given a sign list   [s1,s2,...,sn] , we construct its revised sign list, [t1,12,..., tn], as follows:

• 1)    if [ s i , s _{i + 1},..., s i ₊ j ] is a section of given list, where

s i ^ 0, s i _{+ 1} = • • • = s i ₊ j 1 = 0, s i ₊ j ^ 0, then we replace the

subsection [si+1, si+2, • • •, si+j. 1 ] by the first terms of [-s, - s, s, s, - s,- s, s, si., • • •].

• 2)    Otherwise, let t_k = s_k , i.e. no changes for the others.

Lemma 3 . Given a polynomial f ( x ) with real coefficients, if the numbers of the sign changes of the revised sign of { D 1 ( f ) , D ₂ ( f ) ,... D n ( f ) } is д , then the number of the pairs of distinct conjugate imaginary roots of f ( x ) equals to д . Furthermore, if the number of nonvanishing members of revised sign list is l , then, the number of the distinct real root of f ( x ) equals l - 2 д .

Lemma 4 . Let { d 1 , d ₂ ,..., d _{2 n + 1} } be the principal minor sequence of Discr ( f ) , the discrimination matrix of polynomial f ( x ) with f ( 0 ) ^ 0 • Denote the number of sign changes and the number of non-vanishing members of revised sign of sequence { d 1 d ₂, d ₂ d 3 ,..., d ₂ n d ₂ n _{+ 1} } by д and 2 m , respectively. Then, the number of distinct negative roots of f ( x ) equals to д - 2 m .

According to the Lemma 3 and Lemma 4, we can give the explicit condition to describe the number of the distinct real roots and the negative roots.

Let us consider the following polynomial with real coefficients. All the computations are made with the Computer Algebraic System MAPLE14.

f (to) = to + a1to+a 2

• f . ( to ) = to + a 1 to ² + a₂to + a ₃

• f . ( to ) = to + a 1 to ³ + a ₂to + a₃to + a 4

Using the method discussed above, it is not difficult to establish the explicit criterion for stability, which is expressed by polynomials in the coefficient of the given polynomial.

• 1)    f ( to ) e H [ f ] if and only if a 1 > 0, a ₂ > 0 .

• 2)    f , ( to ) e H [ f ] if and only if a 1 > 0, a 1 a ₂ >  a ₃, a ₃ > 0 .

• 3)    f . ( to ) e H [ f ] if and only if a 1 > 0, a 1 a ₃ + a ₂ ² - 4 a 4 > 0, a 1 a ₂ a ₃ - a 1 ² a ₄ - a ₃ ² > 0, a ₄ > 0 .

• III.    A PPLICATIONS

• A.    The stability of control system

The stability analysis is one of the most important aspects in the design of control systems. A discrete time system is called stable if all the zeros of the system characteristic polynomial lie in the unit circle. There are some methods to check the stability of a given discrete time system [13][14]. However, in some cases, instead of direct analysis on the discrete time system, it is more convenient to perform stability analysis on an equivalent continuous time system. The stability of a continuous time system is determined by the root locations of the system characteristic polynomial with respect to the imaginary axis, the system being stable if and only if all roots lie in the open left half complex plane. The characteristic polynomial of a discrete time system is called Schur stability if the corresponding discrete system is stable, and the same as Hurwitz stability of a continuous system characteristic polynomial. There exist also many methods to check the stability of a given continuous time system [1][2][15]. Using a bilinear transformation, the determination of Schur stability can be converted to the determination of Hurwitz stability of an equivalent polynomial. As a matter of convenience, we marked f ( z ) e S [ f ] if f ( z ) being Schur stability, and g ( w ) e H [ g ] if g ( to ) being Hurwitz stability.

Now we can describe our algorithm for determining the Hurwitz stability of the characteristic polynomial of a continuous time system as follows. The correctness of the algorithm is guaranteed by the above discussion. Without loss of generality, we assume that the polynomials in following algorithm have no multiplicative roots. Otherwise, we may use the algorithm [11][16] to determine the multiplicities of real roots for a given polynomial.

Step 1. Input the monic polynomial f ( to ) , and then obtain the polynomial f ( x , y ) and f 2 ( x , y ) by acting the substitution to = x + i • y upon f ( to ) .

Step 2. Compute the resultant res ( f , f 2 , y ) , and then analysis the distinct real roots of res ( f , f 2 , y ) , which denoted by l , by using Lemma 3. Analysis the distinct negative real roots of res ( f₁ , f 2 , y ) , which denoted by m , by using Lemma 4.

Step 3. Deduce the conclusion by the Theorem 2.

Our method has been implemented by using Computer Algebraic System MAPLE14 for numerical and algebraic manipulations. The following example illustrated our approach.

Example 1. Consider the determination of Hurwitz stability of the following continuous time polynomial f (о) = о3 -Ю + 1.                   (9)

It is not difficult to deduce that l = 2 and m = 1 . Thus, according to Theorem 2, the polynomial f ( о ) is not Hurwitz stability. That is to say, the corresponding continuous time system is not stable.

Example 2. Consider the determination of Hurwitz stability of the following polynomial f (о) = 2.16о4 -0.42о3 + 6.58о2 -0.42о+2.16. (10)

l = 2 and m = 0 . According to Theorem 2, f ( о ) is not Hurwitz stable.

Example 3. Consider the determination of Hurwitz stability of the following polynomial f (о) = 1.70396616о8 + 4.49238296о7 +12.29115880 • о6 + 16.93754704о5 + 22.042 343 84о4 + 16.93754704о3 +12.29115880о2 + 4.49238296о+1.70396616.               (11)

l = m = 2 . According to Theorem 2, f ( о ) is Hurwitz stable.

• B.    The pathological case of the bilinear transformation Denote

g ( о ) = f (° -у ) ( о - 1 ) n ,           ⁽¹²⁾

V о - 1 )

Then f ( z ) e S [ f ] if and only if g ( w ) e H [ g ] .

However, there exists a pathological case while using bilinear transformation to determine the Schur stability via the Hurwitz stability.

Considering the following discrete-time polynomial:

Example 4 . f ( z ) = z ³ - 1 . 21 z ² - 0 . 20 z + 0 . 41

Using the bilinear transformation z = ( о + 1 )/( о - 1 ) ,                (13)

We obtain

g ( w ) = f | w ⁺ M ( ^w -1 ) 3 = 0-76 w ² + 5.64 w +1.6 .    ⁽¹⁴⁾

V w -1J

It has two real roots w =-0.2954505652, w2 =-7.125602066.

Thus, g ( w)e H [ g ], furthermore, f (z) e S [ f ] . However, the above polynomial f ( z) has a root at z = 1 , which leads to non-Schur stability. Indeed, there exist much more examples verifying the pathological case while the discrete time polynomial has a root at z = 1 .

Let us introduce an alternative transformation without any question.

Given a polynomial with real coefficients

n f (z) = £aizn-i (i = 0,1,2,...,n;a0 * 0).( i=0

Using the transformation z = (1 + о + о2 )/(1 - о + о2) ,(16)

Then we get the polynomial g0 (о) = g (о)/(1 - о + о2)n,

Where

g ( о ) = b 0 о ^m + Ь ₁ о ^{2 m - 1} +... + b _{2 m} ( b j e ° , j = 0,1,2,...,2 m )

and m <  n .

Let о = X + iy , which leads to

12     |1+ о + о ²

|1- о + о ²

( 1+ x + x ² - y ² ) ² + ( 2 xy + y ) ²

( 1- x + x ² - y ² ) ² + ( 2 xy - y ) ² '

Denote g ( x , y ) as the expression when the numerator minus the denominator of (18), it is easy to simplify that g ( x , y ) = 4 x ( x ² + y ² +1 ) .              ⁽¹⁹⁾

Combining (18) and (19), we obtain the following result

• 1)    | z | < 1 / g ( x , y ) < 0 / x < 0 / о lies in the open left half plane;

• ²⁾    | z | = 1 / g ( x , y ) = 0 / x = 0 / о ^{lies in the} imaginary axis;

• 3)    | z | >  1 у g ( x , y ) > 0 У x > 0 У о lies in the open right half plane.

According to the discussion above, we obtain the following theorem.

Theorem 3. The transformation (16) maps the imaginary axis to the unit circle and vice versa, the open left half plane to the inside of unit circle and vice versa, the open right half plane to the outside of unit circle and vice versa.

Theorem 4. f ( z ) e S [ f ] J g ( о ) e H [ g ] , where the definitions of f ( _z ) and g ( о ) are the same as mentioned above.

According to Theorem 4, determining the stability of a given discrete time system can be converted to the stability of the corresponding continuous time system through the transformation (16).

Now, we reconsider example 4. By using (16), we obtain g (о) = 0.76о5 + 5.64о4 + 3.12о3 + 5.64о2 + 0.76о Obviously, g(0) = 0, thus g(о)^ H[g]. According to Theorem 4, f (z )^ S [ f ].

Example 5. Consider the determination of Schur stability of the following discrete time polynomial f (z) = z2 - 0.05z + 1.21.               (20)

Using the transformation (16), we have

g ( го ) = 216 го ⁴ - 0.42 го ³ + 6.58 го ² -

0.42 го + 2.16 6 H [ g ] .              (21)

Thus, f ( z ) б 5 [ f ] .

Example 6. Consider the determination of Schur stability of the following discrete time polynomial f (z) = z4 + 0.338z3 + 0.28006z2 + 0.0800038z +

0.00590236                          (22)

Using the transformation (16), we have g (го) = 1.70396616®8 + 4.49238296®7 +

12.29115880 го ⁶ +16.93754704 го ⁵ +

22.04234384 го ⁴ +16.93754704 го ³ + 12.29115880 ® ² + 4.49238296 го + 1.70396616                         (23)

The above polynomial has four pairs of conjugate imaginary roots which are all own negative imaginary part. As a result, the corresponding continuous time system is stable. Thus, the discrete time system considered is stable, too.

• C.    The inertie of Bezout matrix

Let u (x) and v (x) be two polynomial in integer domain of degree n and m , n > m , respectively, u (x) = £ ux , v (x) = £ vx , Un, vm (x) ^ 0. i=0                       i=0

A n x n Bezout matrix B = (bi j )e Cnxn associated with u (x) and v (x) is defined by u (x) v (y) - u (y) v (x) x - y

n

= J j ^- y j ^{- 1}

i , j = 1

The entries of B can recursively be computed by means of the formula bi,j+1 = bi+1,j + u^j - ji, (25)

Complemented with the initial conditions b i 0 = b n _{+ 1} j. = 0 . This rule clearly shows that the entries of B are polynomials of degree 2 in the coefficients of u ( x ) and

v ( x ) .

Definition 6. (See [17][18].) The inertia of Bezout matrix B is defined by a triple

In ( B ) = ( n ( B ) , v ( B ) , J ( B ) ) , (26) where n ( B ) , v ( B ) and £ ( B ) are, respectively, the numbers of Eigen values of B with positive, negative, and zero real parts.

It is well known that all the classical tools for investigation of the roots of algebraic equations, such as the Sturm theorem and continued-fraction criteria, can be proved purely algebraically by computing the inertia of Hankel and Bezout matrices associated with suitable real polynomial u ( x ) and v ( x ) of degree n and less than n , respectively[19][20]. In [20], the authors solved certain classical zero-location problems with the help of the inertia of Bezout and Hankel matrix.

In [18], the author revealed the relation between the inertia of Bezout matrix associated with u ( x ) and v ( x ) and the signs of leading coefficients of every quotient in quotient sequence generated by using Euclidean scheme to real polynomials u ( x ) and v ( x ) .

In [21] the evaluation of the inertia of real Bezout matrix is obtained by computing a block LDL t factorization where L is a lower triangular matrix with unit diagonal entries and D is a block diagonal matrix. By Sylvester’s law[19], the factorization preserves the inertia, moreover, the diagonal blocks in D are themselves Hankel matrices , so their inertia can be computed by direct inspection using a set of rules given by Iochvidov’s rule[22]. However, it is difficult to control the growth of coefficient s in these algorithms.

A fast and fraction-free procedure for computing the inertia of Bezout matrix was presented in [23], and also it derived the new method to determining the numbers of different real roots and different pairs of conjugate complex roots of a polynomial with integer coefficients. Yet it depends on a fast method to determine the signs of the leading coefficients of every quotient in quotient sequence generated by applying Euclid’s algorithm.

In this section, based on the discussion and results in section II, we proposed a new method to compute the inertia of Bezout matrix with polynomials u ( x ) and

v ( x ) .

Algorithm: INER

Input: two polynomials u ( x ) and v ( x ) with integer coefficients.

Output: the inertia of Bezout matrix associated with u ( x ) and v ( x ) .

Step 1. Compute the Bezout matrix B associated with u (x) and v(x) according to the formula (24), then get characteristic polynomial f (го) of it. Using the transformation го = x + iy upon f (го), we can get two bivarate polynomials f1 ( x, y ) and f2 (x, y) denoted by

(3) and (4), respectively. And also, we can deduce the 1c(f,y) , 1c(f.,y) and res(f1,f,,y) mentioned above.

Step 2 . Analysis the numbers of zero solution of

GLCx (f1, f.) andres(fx, f.,y), denoted by n1,oandn2,0, respectively.

Step 3 . Let

GLC, ( xx ) = f ₁ , f ₂

______________ GLC x ( f 1 , f 2 ) ______________ gcd ( GLC x ( f_x , f 2 ) , diff ( GLC x ( f_x , f 2 ) , x ) )

Compute the number of distinct negative roots, denoted by n ₁ , _- , by using Lemma 4 to GLC x ( f 1 , f 2 ) . Applying

Lemma 3, we can obtain the number of distinct roots of

GLCx (f1, f,) , denoted by ц . Then, the number of positive roots is ц - n _, denoted by n .

1, -                            1, +

Step 4. Let res ( f, f 2, У )

res ! , f 2 ⁽ x ) =

gcd ( res ( f l , f , , y ) ) , diff ( res ( f l , f , , y ) , x )

Then, n and n will be obtained, when we apply the

2,+             2, same method of Step 3 to the polynomial res   (x).

f ₁ , f ₂

Step 5. Compute the inertia

^{( n} 2,+ - ⁿ i,+ , ⁿ 2,0 - ⁿ i,0 , ⁿ 2,- - n i,- ) .

Our method has been implemented by using the routines by MAPLE14 for numerical and algebraic manipulations. The following two examples illustrate our approach.

Example 4. (See [24].) Let us consider the 4 x 4 Bezout matrix B associated with the polynomials u (x ) = 3 + 4 x3 + 8 x4,(27)

and v (x ) = 1 + 2 x + 4 x2.(28)

We find that the Bezout matrix associated with u ( x )

and v ( x ) is
	"- 6 - 12 4   8 "
	- 12   4   16 16	.             (29)
	4    16 32 32
	8    16 32 0
And the characteristic polynomials is
f ( о ) = о 4 -30 о3 1848 о 2 - 3008 о +176128 . (30)

Using transformation о = x + iy , we get the following polynomials:

f (x, y) = y4 + (-6x2 + 90x + 1848) y2 + x4 - 30x3 +1848x2 - 3008x +

176128,(31)

and f 2 (x,y) = (-40x + 30)y3 +

(4x3 -90x2 -3696x -3008)y.(32)

It is obviously that n = n = n _ = 0 , and the

1,+         1,0

polynomial res^ _f ( x ) may get as following:

( x ⁴ + 176128 - 3008 x - 1848 x ² - 30 x ³) •

( x ⁶ - 32848 x ² - 2625780 x - 45 x ⁵ +

24345x3 - 249x4 - 5223856) •

It is not difficult to show that n = n = 2 , n = 0 .

2,+       2,

Thus we have that the inertia is ( 2, 0, 2 ) , which the same as the result in [24].

A CKNOWLEDGMENT

This work was supported by a grant from the National Science Foundation of Ningbo City under grant No.2009A610078, the research foundation of Ningbo University under grant No.xkl09042 and No.xkl09043. The authors are deeply grateful to the referee for valuable remarks and suggestions.