Multidimensional or D systems (Bose (1984)) are systems of functional equations whose unknown functions depend on independent variables. The stabilizability and stabilization problems are fundamental issues in the study of multidimensional systems in control theory. Nowadays, the problem is well-understood in the case of 1D systems whereas progress for D systems with are rather slow. One approach for handling stabilizability or stabilization issues in systems theory is the fractional representation approach (Vidyasagar (2011)) in which a plant is represented by its transfer matrix where . This transfer matrix admits a left factorization (also called fractional representation of ), where the matrices satisfying and have entries in the integral domain of structurally stable rational fractions, namely the ring of rational functions in which have no poles in the closed unit polydisc of defined by:
Introducing the matrix , it is known (see Quadrat (2003b, a)) that the multidimensional system given by the transfer matrix is then internally stabilizable if and only if the -module is a projective -module of rank , where the closure of in is defined by: This projectivity condition is in turn equivalent to the fact that the reduced minors of the matrix do not have common zeros in (see also Lin (1998)). In other terms, if we denote by the reduced minors of , i.e., the minors of divided by their gcd, by the polynomial ideal generated by the ’s, and by the associated algebraic variety, then the system is internally stabilizable if and only if .
The first contribution of the present paper is to provide an effective algorithm for testing the stabilizability condition for the class of D systems for which the ideal is zero-dimensional, i.e., the ’s have only a finite number of common complex zeros (i.e., consists of a finite number of complex points). Note that this class includes the class of 2D systems. Our main idea is to take advantage of the univariate representation for zero-dimensional ideals (Canny (1988); Becker and Wörmann (1996); Alonso et al. (1996); Rouillier (1999)). This concept, which can be traced back to Kronecker (1882), yields a one-to-one correspondance between the elements of and the zeros of a univariate polynomial . Numerical techniques can thus be applied to compute certified numerical approximations of the roots of and then of those of .
In the case of a stabilizable plant, the next step consists in computing a stabilizing controller which can be achieved by computing a stable (i.e., devoid from zeros in ) polynomial (see Lin (1988)). The polydisc Nullstellensatz, proved by Bridges, Mines, Richman and Schuster (see Bridges et al. (2004)), shows that the existence of a stable polynomial is equivalent to . Several proofs of this result have been investigated in the literature, mainly for the case where is a zero-dimensional ideal (see Raman and Liu (1986); Lin (1988); Bisiacco et al. (1986); Guiver and Bose (1995) and Xu et al. (1994) for instance). Nevertheless none of them is effective in the sense that it provides an algorithm for computing using calculations that can be performed in an exact way by a computer. Indeed, starting from a set of polynomials with rational coefficients (), these algorithms are built on spectral factorization, i.e., factorization of polynomials in into stable and instable factors. For irreducible polynomials in , this factorization requires the explicit computation of the complex roots of the polynomials, which can be done only approximately. This leads to approximate (stable) polynomials that do not belong to the polynomial ideal. As a consequence, these algorithms are able to solve the aforementioned problem only for few simple systems (see Sections 4 and 5 for details).
Our second contribution is to provide an effective algorithm for computing a stable polynomial for the class of systems for which is a zero-dimensional ideal. Our symbolic-numeric method roughly follows the lines of that proposed in Xu et al. (1994) but once again we take advantage of the univariate representation of zero-dimensional ideals (Rouillier (1999)) to control the numeric precision required to achieve our goal.
The paper is organized as follows. In Section 2, we recall some classical computer algebra results on the complex zeros of polynomials and polynomial systems. We also introduce the univariate representation of zero-dimensional ideals which will be our main tool in what follows. In Section 3, we provide an effective stabilizability test, i.e., an algorithm for testing whether a zero-dimensional ideal intersects the closed unit polydisc. In Section 4, we provide an effective polydisc Nullstellensatz namely a symbolic-numeric method for computing, if it exists, a stable polynomial in a zero-dimensional polynomial ideal. Finally, in Section 5, we illustrate our methods on one example and show some running times of prototype implementations.
2 Preliminaries on algebraic systems
In this section, we introduce some notations and we recall some classical material about the computation of certified numerical approximations of the complex zeros of polynomials and polynomial systems.
The bit-size of an integer is the number of bits in its representation and for a rational number (resp., a polynomial with rational coefficients) the term bit-size refers to the maximum bit-size of its numerator and denominator (resp., of its coefficients). For a complex number , we denote by (resp., ) its real (resp., imaginary) part. If , we write if both and . For such that , we shall consider the axes-parallel open box or box for short and its width is defined by . We also introduce the non-negative real number . The box is said to be of rational endpoints if have rational real and imaginary parts, i.e, . Finally, the box is called isolating for a given polynomial if it contains exactly one complex zero of .
The following result concerns the isolation of the complex zeros of a univariate polynomial. We refer, for instance, to Sagraloff and Yap (2009) for more details. Let be a squarefree polynomial of degree . Then, for all , one can compute disjoint axes-parallel open boxes , with rational endpoints such that each contains exactly one complex root of and satisfies .
In the algorithms given in Sections 3 and 4 below, we shall use a routine called Isolate which takes as input a univariate polynomial , a box , and a precision and computes isolating boxes with rational endpoints for the complex roots of that lie inside the given box and such that . If (resp., ) is not specified in the input, we consider all complex roots in (resp., the boxes are computed up to a sufficient precision for isolation).
Let us now recall a standard property about width expansion through interval arithmetic in polynomial evaluation. Here we consider exact interval arithmetic, that is, the arithmetic operations on the interval endpoints are considered exact (see Alefeld and Herzberger (2012)). If is a bivariate polynomial of two real variables and and a box, we denote by the interval that results from the evaluation of the polynomial at the box using interval arithmetic. [Cheng et al. (2010), Lemma 8] Let be a box with rational endpoints satisfying and let be a bivariate polynomial of two real variables and of degree with coefficients of bit-size . Then, can be evaluated at the box by interval arithmetic into an interval of width at most . In particular, a direct consequence of Lemma 2 is that if , then we have .
We now consider a set of polynomials in . We denote by the ideal generated by the ’s and by
the complex variety of their common zeros. In the sequel, we shall always assume that the ideal under consideration is a zero-dimensional ideal, that is, that the ’s have only a finite number of common complex zeros, i.e., consists of a finite number of complex points. Methods for computing certified numerical approximations for the elements of usually proceed in two steps. First, a formal representation (as for instance, a Gröbner basis, a triangular decomposition or a univariate representation) of the set is computed. Then, this formal representation is used to compute, more or less easily, numerical approximations of the elements of . The convenient representation of the variety of zero-dimensional ideals that we shall use in the sequel is the so-called univariate representation introduced by Rouillier in Rouillier (1999). With the previous notation and assumptions, a univariate representation of is the datum of a linear form , with as well as univariate polynomials such that the following two applications
provide a one-to-one correspondence between the elements of and the zeros of . The univariate representation then defines a bijection between the zeros of and those of that preserves the multiplicities and the real zeros. One of its main advantage is that it permits to study the variety , i.e., the set of common zeros of the polynomials ’s of variables through the roots of the univariate polynomial . In particular, if we are interested in computing certified numerical approximations of the complex elements of , a first step consists, in isolating the complex roots of the polynomial (using for instance the algorithm in Sagraloff and Yap (2009)). Then, numerical approximations of the coordinates of can be obtained by evaluating at the obtained numerical boxes. Lemmas 2 and 2 above show that we can thus obtain a given precision for the elements of .
The computation of a univariate representation of can be done by pre-computing a Gröbner basis of the ideal Cox et al. (1992)
, and then performing linear algebra calculations in the quotient vector spacewhich admits a basis composed of the monomials that are irreducible modulo the Gröbner basis . For more details, see Rouillier (1999). For the specific case of bivariate algebraic systems, it should be stressed that very practically efficient algorithms exist for computing such a representation: see Bouzidi (2014). To achieve efficiency, these algorithms replace the “costly” computation of a Gröbner basis of by resultants and subresultants computations, see Basu et al. (2006).
3 An effective stabilizability test
Using the fractional representation approach to multidimensional systems, in Section 1, we have seen that a plant is stabilizable if and only if for a certain polynomial ideal that can be explicitly described in terms of the transfer matrix of the plant.
Given polynomials , the purpose of this section is to provide an effective algorithm to decide whether or not where is a zero-dimensional ideal. To achieve this, we shall use a symbolic-numeric approach. We start by computing a univariate representation of . As explained in Section 2, such a representation allows to describe formally the elements of as
where . In what follows, the degree of is denoted by , and those of are then smaller than (see Rouillier (1999)).
Using a univariate representation, one can compute a set of hypercubes in isolating the elements of . Each coordinate is represented by a box in obtained from the intervals containing its real and imaginary parts. Moreover, from Lemmas 2 and 2, these hypercubes can be refined up to an arbitrary precision. We shall now consider the intersection between those hypercubes and the closed unit polydisc of defined by:
Below, for any , we shall denote by the bivariate polynomial , where (resp., ) is the real (resp., complex) part of the polynomial resulting from .
From the definition of , one can see that the situation is easier when does not contain elements with for some . Indeed, we have: With the previous notations, let us consider such that, for all , . Let be an isolating box for the root of corresponding to in the univariate representation of . Then, there exists such that if , then, for all , the interval does not contain zero. Let (resp., ) denote an upper bound on the degree (resp., bit-size) of the polynomials , . For , the real (resp., imaginary) part of is a bivariate polynomial in and of degree (resp., bit-size) bounded by (resp., ). Consequently, the bivariate polynomial has degree and bit-size respectively bounded by and respectively. Now, for , let and let be an isolating box for the root of corresponding to in the univariate representation of (1), and such that . From Lemma 2, if we refine so that , where , then, we have that for all , the interval satisfies so that, by definition of , it does not contain zero. Therefore, if does not contain elements with one coordinate in the unit circle, one can easily test the stabilizability condition . Indeed, with the previous notations, we isolate the roots of inside boxes . Then, for , we refine until, for all , the interval does not contain zero. If one of the intervals is included in , we proceed to the next box , otherwise we have found an element in so that the system is certainly not stabilizable. After having investigated all the boxes , we can then conclude about the stabilizability of the system.
We shall now consider the case where contains (at least) one element having some coordinates on the unit circle. In this case, we cannot proceed numerically as before since if , then, using the above notations, we cannot fulfill the condition that the interval does not contain zero. To guarantee the termination of the algorithm, we shall then have to compute, for each variable , the number of elements satisfying . Let be a zero-dimensional ideal and the associated algebraic variety. Then, for all , one can compute the non-negative integer . Let be a univariate representation of and . Computing the resultant of the polynomials and with respect to the variable we get a univariate polynomial that can be written , where the multiplicity corresponds to . Then, using the classical Bistritz test (see Bistritz (2002)), one can compute the number of complex roots counted with multiplicity of that lie on the unit circle and obtain the non-negative integer .
An alternative to the Bistritz test, which we use in practice, consists in applying to the polynomial , the Möbius transform , which maps the complex unit circle to the real line . The number of complex roots of on the unit circle is then given as the number of real roots of the gcd of two polynomials and , where (resp., ) is the real (resp., complex) part of the numerator of the rational fraction .
Using Lemma 3, we can test the stabilizability condition as follows. We start with the variable . We refine the isolating boxes for the roots of until exactly intervals contain zero. We throw away the boxes ’s such that the interval is included in and we proceed similarly with the next variable . If at some point we have thrown away all the boxes ’s, then the system is stabilizable. Otherwise the boxes which remain at the end of the process lead to elements of so that the system is not stabilizable.
We summarize our symbolic-numeric method for testing stabilizability in the following IsStabilizable algorithm.
4 An effective stabilization algorithm
When a multidimensional system is stabilizable, one is then interested in computing a stabilizing controller. Within the fractional representation approach, this problem reduces to the following task (see Section 1 or Lin (1988)): given a polynomial ideal satisfying , compute a polynomial such that belongs to the ideal and is a stable polynomial, i.e., is devoid from zeros in .
For an ideal , the polydisc Nullstellensatz (see Bridges et al. (2004)) asserts that the existence of a stable polynomial is equivalent to . [Polydisc Nullstellensatz] Let us consider a polynomial ideal such that . Then, there exist polynomials in such that:
Several proofs of Theorem 4 have been invesigated in the literature. Nevertheless none of them is effective in the sense that it provides an algorithm for computing and the cofactors ’s using calculations that can be performed in an exact way by a computer. In Xu et al. (1994), the authors study D systems (i.e., ) for which the ideal under consideration is zero-dimensional. The idea of their method for computing a stable polynomial is to compute univariate elimination polynomials and with respect to each variable and and to factorize them into a stable and an unstable factor, i.e., for , , where the roots of (resp., ) are outside (resp., inside) the closed unit disc , and then, to compute the stable polynomial as . However, this approach presents a major drawback with respect to the effectiveness aspect. Indeed, when the elimination polynomial (resp., ) is an irreducible polynomial in (resp. ), its stable factor (resp., ) could not be computed exactly since it will have coefficients in , and thus, only an approximation of this polynomial can be obtained. As a consequence, the polynomial will not belong to the ideal .
In the sequel, we present a symbolic-numeric algorithm for computing and the ’s that follows roughly the approach of Xu et al. (1994) while we provide a way for tackling the effectiveness issue. Our main ingredient is the univariate representation of zero-dimensional ideals which allows us to compute and refine approximate factorizations over of the elimination polynomials ’s.
Let be a zero-dimensional ideal such that . For simplicity reasons, in what follows, we further assume that the ideal is radical, i.e., . The elements of are given by a univariate representation
where , for , , and . Since is a radical ideal, the polynomial is a squarefree polynomial and for distincts . Moreover, from Definition 2, if we introduce the polynomial ideal , then we have . In particular, if , then .
Let us first explain how we can compute approximations of the stable polynomials appearing in the method of Xu et al. (1994) sketched above. Using Lemma 2, we can compute a set of boxes with rational endpoints, isolating the distinct complex roots of . Then, according to the stabilizability condition , for all , the box can be refined so that there exists satisfying . We then set to the midpoint of the refined box and we add the factor to the polynomial . We finally obtain a set of stable univariate polynomials , such that .
Let us now introduce the polynomial . By construction has rational coefficients and it vanishes on , where the polynomial ideal is defined by with . Hence, according to the classical Nullstellensatz theorem (Cox et al. (1992)), belongs to the ideal so that there exist polynomials such that . Moreover can be explicitly computed as the quotient of the Euclidean division of by in .
We shall now show that if we refine enough the boxes ’s isolating the roots ’s of , then the stable polynomial that we are seeking for can be obtained from the polynomials and constructed as explained above. For 111small enough so that the previous process can be applied., we denote , , , and the objects constructed by the previous process where the roots of are isolated up to precision (i.e., , for all ). Using the previous notations, the main result of this section can be stated as follows: The polynomial belongs to the ideal . Moreover, there exists such that the polynomial is a stable polynomial. The proof of Theorem 4, given below, requires the following lemma. For , the polynomial has coefficients bounded by , where is a positive real number that does not depend on . Let and , with , denote the expansion of the polynomials and on the monomial basis. By the standard Vieta’s formulas, for all , we have:
By assumption, for all , so that
Now, we can write , where the ’s denote the symmetric functions associated to and, in particular, . Consequently, since , we get:
On the other hand, the polynomial can be computed as the quotient of the Euclidean division of by . Formally, these polynomials can be considered as polynomials in , where are considered as new indeterminates so that their quotient denoted by can be computed independently from . The coefficients of are thus bounded by a certain positive real number that does not depend on (use, for instance, Mignotte’s bound Mignotte (1989)). Now, for , since , we have . Thus if we denote , then the coefficients of the evaluation of for particular values of in are bounded by . Finally, we have proved that the coefficients of are bounded by , which ends the proof. We are now in position to give a proof of Theorem 4. With the previous notations, for all , we have so that vanishes on , which implies . Let us now prove that we can choose so that , viewed as , is stable, i.e., , . According to Lemma 4, for , we have , where (resp., ) denotes a positive real number (resp., the degree of the polynomial ) that does not depend on . On the other hand, the polynomial can be written as for some not necessarily distinct. Thus, for all , we have . Now, if we denote by , the minimum distance between the elements of and , then we have , which yields:
Finally, for sufficiently small , we have so that:
which ends the proof. The following StablePolynomial algorithm summarizes our method for computing a stable polynomial in a zero-dimensional ideal satisfying . The routine IsStable is used to test if a polynomial