=Paper= {{Paper |id=Vol-2928/paper11 |storemode=property |title=Application of the Gröbner Basis Method for the Study of Nonlinear Control Systems |pdfUrl=https://ceur-ws.org/Vol-2928/paper11.pdf |volume=Vol-2928 |authors=Sergey N. Chukanov }} ==Application of the Gröbner Basis Method for the Study of Nonlinear Control Systems== https://ceur-ws.org/Vol-2928/paper11.pdf
 Application of the Gröbner Basis Method for the Study
              of Nonlinear Control Systems

                                            Sergey N. Chukanov
                                 Sobolev Institute of Mathematics SB RAS,
                                            Novosibirsk, Russia
                                               ch sn@mail.ru




                                                       Abstract
                      The paper considers methods for estimating stability using Lyapunov
                      functions, which are used for nonlinear polynomial control systems.
                      The apparatus of the Gröbner basis method is used to assess the sta-
                      bility of a dynamical system. To apply the method, the canonical rela-
                      tions of the nonlinear system are approximated by polynomials of the
                      components of the state and control vectors. The equilibrium states of
                      a nonlinear polynomial system are determined as solutions of a non-
                      linear system of polynomial equations. An example of determining the
                      equilibrium states of a nonlinear polynomial system using the Gröbner
                      basis method is given. The application of the Gröbner basis method
                      for estimating the attraction domain of a nonlinear dynamic system
                      with respect to the equilibrium point is considered. The coordina-
                      tion of input-output signals of the system based on the construction of
                      Gröbner bases is considered.




Introduction
    The most of the dynamical systems in technology and nature are nonlinear dynamical systems. The canonical
relations of a nonlinear system can be approximated by polynomials of the components of the state and control
vectors. Stability testing using the method of Lyapunov functions is widely applied to nonlinear systems.
   There are several methods in the literature to identify candidates for Lyapunov functions [Krasovsky59]:
    decomposition of the sum of squares [Papach02];
    using the Gröbner basis to select parameters [Forsman91];
    use of homotopy operators for decomposition of the vector field of states of the system [Chukanov12, Edelen85];
    the assumption that the derivative of the Lyapunov function is negative definite, and then obtain by integra-
tion and check the positive definiteness (gradient method).
   Gröbner bases are used to solve problems in the theory of nonlinear systems. Some of the applications of the
Gröbner basis can be named: estimation of equilibrium states of a nonlinear system; finding the critical points
of a given nonlinear system with the Lyapunov function; coordination of input-output signals of the system.
   Gröbner bases facilitate the solution of a system of multidimensional polynomial equations in the same way as
the Gaussian elimination algorithm makes it possible to solve a system of linear algebraic equations. In lexical

Copyright ⃝
          c by the paper’s authors. Use permitted under Creative Commons License Attribution 4.0 International (CC BY 4.0).
In: Sergei S. Goncharov, Yuri G. Evtushenko (eds.): Proceedings of the Workshop on Applied Mathematics and Fundamental
Computer Science 2021 , Omsk, Russia, 24-04-2021, published at http://ceur-ws.org




                                                            1
ordering the Gröbner basis has a triangular structure, reminiscent of the triangular structure in the Gaussian
elimination method.
   The theory of control of dynamic objects can be divided into two subgroups [Khalil02]:
   (1) systems in which the principle of superposition operates, and linear control methods can be used;
   (2) systems in which the superposition principle does not work, and it is necessary to use nonlinear control
methods. To improve the quality of the dynamic object control system, it is necessary to take into account the
nonlinear features of the system.

1    Gröbner bases
      The objects in the theory of Gröbner bases are polynomial ideals and algebraic varieties [Nesic02]. Let
p1 , ..., ps be multidimensional polynomials in variables x1 , ..., xn , whose coefficients lie in the field k (we will
consider the field of real numbers R). The variables x1 , ..., xn are considered ”place markers” in the polynomials:
p1 , . . . , ps ∈ R[x1 , . . . , xn ]. Algebraic variety defined by the polynomials p1 , . . . , ps is the collection of all solutions
in Rn of the system of equations:

                                                         p1 (x1 , ..., xn ) = 0,
                                                         ...                                                                     (1)
                                                         ps (x1 , ..., xn ) = 0.
Formally:

                            V (p1 , ..., ps ) := {(a1 , ..., an ) ∈ Rn : pi (x1 , ..., xn ) = 0, i = 1, ..., s}                  (2)
  The polynomial ideal I, which is generated by p1 , ..., ps , is a set of polynomials obtained by combining these
polynomials by multiplying and adding with other polynomials:
                                                  {                                       }
                                                       ∑s
                            I = ⟨p1 , ..., ps ⟩ := f =      gi pi : gi ∈ R[x1 , ..., xn ]                      (3)
                                                                   i=1

   The polynomials pi , i = 1, . . . , s form the basis of the ideal I. A useful interpretation of the polynomial ideal
I is in terms of the equations (3). Multiplying gi by arbitrary polynomials gi ∈ R[x1 , ..., xn ] and adding them,
we get the consequence from (1):

                                                    f = g1 p1 + . . . + gs ps = 0,

and f ∈ I. Therefore, I = ⟨p1 , . . . , ps ⟩ the ideal contains all ”polynomial consequences” of the equations (3).
   The Gröbner basis method is based on the concept of monomial ordering (a monomial is a polynomial consist-
ing of one term), since it introduces a corresponding extension of the concept of a leading term and a leading coeffi-
cient, familiar for one-dimensional polynomials, to multidimensional polynomials. Let’s consider lexicographic or
lex order [Nesic02]. Let α, β be two n - tuples of integers α = (α1 , . . . , αn ) ∈ Nn , β = (β1 , . . . , βn ) ∈ Nn . n -tuple
α follows β (in lex order), which is denoted as α ≻ β, if in the difference of vectors α−β = (α1 − β1 , . . . , αn − βn )
the leftmost nonzero element is positive. For the polynomial f = x31 x2 x33 + 2x31 x43 using lex order x1 ≻ x2 ≻ x3
results in x31 x2 x33 follows x31 x43 , since the multidegrees of monomials satisfy: (3, 1, 3) ≻ (3, 0, 4). In this order, the
leading coefficient and the leading term are respectively LC(f ) = 1 and LT (f ) = x31 x2 x33 . When using lex of
order x3 ≻ x2 ≻ x1 senior term: LT (f ) = 2x31 x43 , since (4, 0, 3) ≻ (3, 1, 3).
   The ideal I has no unique basis, but for any two different bases ⟨p1 , ..., ps ⟩ and ⟨g1 , ..., gm ⟩ of the ideal I,
the varieties V (p1 , . . . , ps ) and V (g1 , . . . , gm ) are equal; the variety depends only on the ideal generated by its
defining equations. If all polynomials in a given basis of an ideal have a degree lower than the degree of any other
polynomial in an ideal, then this basis is the simplest. For an ideal I and a given monomial order, we denote the set
of leading terms of elements I as LT (I). The ideal generated by elements from LT (I) is denoted by ⟨LT (I)⟩. The
Gröbner basis is formally defined as a set of polynomials g1 , . . . , gm , for which ⟨LT (I)⟩ = ⟨LT (g1 ), . . . , LT (gm )⟩.
When calculating Gröbner bases, a monomial order is specified. We note two properties of Gröbner bases for a
given monomial order:
   1. Each ideal I ⊂ R[x1 , ..., xn ], different from the trivial ⟨0⟩, has a Gröbner basis.
   2. For the ideal I ⊂ R[x1 , ..., xn ], different from the trivial ⟨0⟩, the Gröbner basis of the ideal I can be
calculated using a finite number of algebraic operations.




                                                                    2
  For a given set of polynomials P , there is an algorithm that computes the Gröbner basis for the (ideal
generated by) P in a finite number of steps [Buchberger70]. Buchberger’s algorithm generalizes algorithms:
Gaussian elimination for a system of linear algebraic equations and Euclid’s algorithm for calculating the greatest
common divisor of a set of one-dimensional polynomials. This algorithm was implemented on computers in
symbolic computation programs using Gröbner bases for solving systems of polynomial equations [?, Wolfram03,
Demenkov15].

2     Finding equilibrium states of a nonlinear dynamical system
    The use of the Gröbner basis in finding solutions to a nonlinear system of polynomial equations is similar to
the application of the Gauss method for solving a quadratic system of linear equations. Consider an example
of reducing a nonlinear system of polynomial equations: p1 = x1 − x22 = 0, p2 = x2 + x23 = 0, p3 = x3 − 2x21 =
0, to a triangular form using the Gröbner basis method for lex order: x1 ≻ x2 ≻ x3 . In the WOLFRAM
MATHEMATICA package, the function call
   GroebnerBasis[{p1, p2, p3}, {x1, x2, x3}, {}]
leads to a triangular Gaussian form of polynomial equations:

                                                    x1 − x43 = 0,
                                                    x2 + x23 = 0,
                                                    −x3 + 2x83 = 0,

which allows us to get a solution to this system.
   Consider a nonlinear system without inputs ẋ(t) = f (x(t)); x, f ∈ Rn , t ∈ R, where f (x) = 0 is a vector of
polynomials in x. The equilibrium states for this polynomial system are obtained as solutions of a nonlinear
system of polynomial equations: f (x) = 0.

Example 1
    Equilibrium states of the [Nesic02] polynomial system:

                                                 ẋ1 = x1 + x2 − x23 ,
                                                 ẋ2 = x21 + x2 − x3 ,
                                                 ẋ3 = −x1 + x22 + x3 ,

can be obtained as solutions of system polynomial equations:

                                               p1 := x1 + x2 − x23 = 0,
                                               p2 := x21 + x2 − x3 = 0,
                                               p3 := −x1 + x22 + x3 = 0.
    The Gröbner basis for the ideal (p1 , p2 , p3 ) using lex order: x1 ≻ x2 ≻ x3 , has the form:

                                           g1 := 4x1 − 2x21 − 4x31 + x41 + x61 ,
                                           g2 := −x21 + x41 − 2x2 + 2x21 x2 ,
                                           g3 := −x1 + x21 + x2 + x22 ,
                                           g4 := −x21 − x2 + x3 .

   Algebraic equations gi = 0, i = 1, 2, 3, 4 has the same solutions as pj = 0, j = 1, 2, 3. The polynomial g4
depends only on x3 ; from the algebraic equation g4 (x3 ) = 0, you can determine x3 . If the numerical value of x3
substitute in g3 (x2 , x3 ) = 0, then you can define x2 ; from g2 (x1 , x2 , x3 ) = 0 you can define x1 .
   In the WOLFRAM MATHEMATICA package: Form an ideal of polynomials:
p1 = x1 + x2 − x32 ; p2 = x12 + x2 − x3; p3 = −x1 + x22 + x3.
   Let us define the Gröbner basis:
 grbas = GroebnerBasis[{p1, p2, p3}, {x3, x2, x1}, {}],
 grbas = {4x1 − 2x12 − 4x13 + x14 + x16 , −x12 + x14 − 2x2 + 2x12 x2,
 −x1 + x12 + x2 + x22 , −x12 − x2 + x3}.
   To find the roots of x1 , we define the reduced Gröbner basis:
grbas = GroebnerBasis [{p1, p2, p3} , {x3, x2, x1} , {x3, x2}] ,




                                                            3
        {                                }
grbas = 4x1 − 2x12 − 4x13 + x14 + x16 .
  Let’s execute: Roots[4x1 − 2x12 − 4x13 + x14 + x16 == 0, x1].
  To find the roots of x2 with known x1 , we define the reduced Gröbner basis:
        GroebnerBasis [{p1, p2,}p3} , {x3, x2, x1} , {x3}] ,
grbas = {
grbas = −x1 + x12 + x2 + x22 .
  Let’s execute: Roots[−x1 + x12 + x2 + x22 == 0, x2].
  To find the roots of x3 with known x1 , x2 , execute:
Roots[−x12 − x2 + x3 == 0, x3]. The results are shown in Table 1.

                                                                                                        Table 1.
                                        x1                    x2                     x3
                Solution 1 :         x1 = −1,          x2 = 0.5 − 1.32i,      x3 = 0.5 − 1.32i,
                Solution 2 :          x1 = 0,               x2 = 0,                x3 = 0,
                Solution 3 :          x1 = 1,               x2 = 0,                x3 = 1,
                Solution 4 :        x1 = 1.18,            x2 = −0.69,            x3 = 0.70,
                Solution 5 :    x1 = −0.59 − 1.74i,   x2 = −2.35 + 1.03i,    x3 = −5.04 + 3.09i,
                Solution 6 :    x1 = −0.59 + 1.74i,    x2 = 1.35 − 1.03i,    x3 = −1.35 − 3.09i,
                Solution 7 :         x1 = −1,         x3 = −0.5 + 1.32i,      x3 = 0.5 + 1.32i,
                Solution 8 :          x1 = 0,              x2 = −1,               x3 = −1,
                Solution 9 :          x1 = 1,              x2 = −1,                x3 = 0,
                Solution 10 :       x1 = 1.18,            x2 = −0.31,            x3 = 1.08,
                Solution 11 :   x1 = −0.59 − 1.74i,    x2 = 1.35 − 1.03i,    x3 = −1.35 + 1.03i,
                Solution 12 :   x1 = −0.59 + 1.74i,   x2 = −2.35 + 1.03i,    x3 = −5.04 − 1.03i.



3     Application of the Gröbner basis method in the theory of the method of Lya-
      punov functions
3.1   Estimation of the area of attraction
     The set of all initial conditions of a dynamical system, which converge to the same equilibrium state, is
called the area of attraction of this equilibrium state [Forsman91, Sidorov19]. One way to get an estimate of the
domain of attraction is to use the Lyapunov functions.
   The standard result of Lyapunov’s theory is that if x = 0 is an equilibrium point for a system with continuous
time: ẋ = f (x), x ∈ D ⊂ Rn , is a domain containing x = 0 and V : D → R is a continuously differentiable
Lyapunov function such that V (0) = 0 and V (x) > 0, V̇ = Vx f (x) < 0, ∀x ∈ D − {0} ; then the point x = 0
is asymptotically stable. For such a Lyapunov function, consider the sets Ω = {x ∈ Rn : Vx f (x) < 0} and
Bd = {x ∈ Rn : V (x) ≤ d} . If there is a value d > 0 such that Bd ⊂ Ω, then the set Bd is an estimate of the
domain of attraction.
   For polynomial systems with a polynomial Lyapunov function V , the Gröbner basis can be used to determine
Bd . You can determine the largest Bd by finding a d such that Bd ⊂ Ω. For polynomial systems with polynomial
Lyapunov functions, V (x) − d and Vx f (x) are polynomials and the boundaries of the sets Bd and Ω are varieties
Z(V − d) and Z(Vx f (x)), respectively. At the points of contact Z(V − d) and Z(Vx f (x)), the gradients V and
Vx f (x) are parallel [Luenberger16]. Using this information, we obtain a system of n + 2 polynomial equations in
n + 2 variables (x1 , . . . , xn , d, λ), where λ is the Lagrange multiplier (see Appendix):

                                              V − d = 0,
                                              Vx f = 0,                                                       (4)
                                              ∇(Vx f ) − λ∇V = 0.

In the case of the vector of Lagrange multipliers λ = (λ1 , . . . , λm )T ∈ Rm we obtain a system of n + m + 1
equations from n + m + 1 variables x1 , ..., xn , d, λ1 , . . . , λm .
   Calculating the Gröbner basis for this system, where the variable d has the lowest rank in the lex order, we
obtain a polynomial equation for d. The smallest positive solution of this equation (the value dmin > 0), is the
best estimate of the area of attraction.




                                                       4
Example 2
      Consider a second-order system:                     (              )
                                                               −x1
                                           ẋ = f (x) =
                                                           −x2 + 2x1 x22
                                                  (       )                                            (           )
                                                    4 2                                                  8x1 + 4x2
and choose the Lyapunov function V (x) = xT                 x = 4x21 + 4x1 x2 + 3x22 , therefore: Vx =               ;
                                                    2 3                                                  4x1 + 6x2
V̇ = Vx f = −8x21 + 12x1 x32 − 8x1 x2 + 8x21 x22 − 6x22 .
   The fact that the gradients are parallel (∇(Vx f ) − λ · ∇V = 0) gives additional equations:
                             {
                                g1 = 8x1 + 4x2 − λ(−16x1 + 12x32 − 8x2 + 16x1 x22 ),
                                g2 = 4x1 + 6x2 − λ(36x1 x22 − 8x1 + 16x21 x2 − 12x2 ).

   Let us calculate the Gröbner basis for four polynomials {V − d, Vx f, g1 , g2 } in ordering: d ≺ x1 ≺ λ < x2 .
   This{ reduces the system to a polynomial:
                                   }             4d4 − 147d3 + 768d2 + 2048d, which results in the values of the
roots:    0 29.71 −1.92 8.97 . The smallest nonzero positive value of d for which there is a solution to the
system is d ≈ 8.97.
   In the WOLFRAM MATHEMATICA package: V d = 4x12 + 4x1x2 + 3x22 − d;
   V xf = −8x12 − 8x1x2 − 6x22 + 8x12 x22 + 12x1x23 ;
   g1 = 8x1 + 4x2 − lam(−16x1 − 8x2 + 16x1x22 + 12x23 );
   g2 = 4x1 + 6x2 − lam(−8x1 − 12x2 + 16x12 x2 + 36x1x22 ).
   grbas = GroebnerBasis[{V
            {                               } {d, x1, lam, x2}, {x1, lam, x2}],
                               d, V xf, g1, g2},
   grbas = 2048d + 768d2 − 147d3 + 4d4 ,
                                                         {                             }
   Roots[2048d + 768d2 − 147d3 + 4d4 == 0, d] ⇒ d = 0 29.71 −1.92 8.97 . 

3.2     Decomposition of the dynamical system vector field ẋ = f (x)
                                                                                                      ∂
    For a dynamic system: ẋ = f (x); x ∈ Rn , f (x) ∈ Rn , f (0) = 0, form a vector field X = f (x) . Form
                                                                   ⟨         ⟩                        ∂x
                                                                          ∂
the corresponding differential form ω = f (x)dx in the dual basis dxi ,         = δij . Let us construct a scalar
                                                                         ∂xj
potential from the vector field X by using the homotopy operator centered at the point x0 = 0 for the form
ω = f (x)dx :
                                         (     )
                                      ∫1    ∂                      ∫1
                               H(ω) =     x       | (f (λx)dx) dλ = xT f (λx)dλ.
                                      0     ∂x                     0

We will assume that φ(x) = Hω(x) is a scalar potential.

Example 3
      Consider an example of dynamic equations:

                                                      ẋ1 = −x1 + x22 ,
                                                      ẋ2 = −x2 − x21 .

  Let’s construct a dual differential form:

                                        ω = (−x1 + x22 )dx1 + (−x2 − x21 )dx2 ,

to which we apply the homotopy operator with x0 = 0 :
                                                  (   )                             (                     )
                                        ∫1                     ∫1 (             )       −λ̃x1 + λ̃2 x22
                    φ(x) = H(ω(x)) =          T
                                             x f λ̃x dλ̃ =            x1   x2                                 dλ̃ =
                                        0                      0                        −λ̃x2 − λ̃2 x21
                         1              1
                      = − (x21 + x22 ) + (x1 x22 − x2 x21 ).
                         2              3
Let us choose the function as the scalar Lyapunov function:




                                                               5
                                            (        )   (               )
                     V (x) = −6 · ϕ (x) = 3 x21 + x22 + 2 x1 x22 − x2 x21 ,  (           )
                                (                                           ) −x1 + x22
                     V̇ = Vx f = 6x1 + 2x2 − 4x1 x2 6x2 − 2x1 + 4x1 x2
                                             2                    2
                                                                                           =
                                                                               −x2 − x21
                        = −6x1 − 6x2 − 6x2 x1 + 6x1 x2 .
                              2      2     2        2


    Let’s find the solution of the system V d = V (x) − d = 0 in the WOLFRAM MATHEMATICA package:

                      V d = 3x21 + 3x22 + 2x22 x1 − 2x21 x2 − d,
                      Vx f = −6x21 − 6x22 − 6x22 x1 + 6x21 x2 ,
                      g1 = 6x1 + 2x22 − 4x1 x2 + lam(−12x1 − 6x22 + 12x1 x2 ),
                      g2 = 6x2 + 4x2 x1 − 2x21 + lam(−12x2 − 12x2 x1 + 6x21 ),
                      grb = {
                            GroebnerBasis[{V   } d, Vx f, g1 , g2 } , {d, x1 , lam, x2 } , {x1 , lam, x2 }],
                      grb = 54d − 29d2 + d3 .

    The roots of the polynomial 54d − 29d2 + d3 are : d1 = 0, d2 = 2, d3 = 27.
    The smallest nonzero positive value d, for which there is a solution to the system: dmin = 2.

4     Conversions of input-output signals of a nonlinear system
    Consider a differential ring - a ring on which the differentiation operation is defined. It is assumed that
differentiation is carried out with respect to the implicit variable t. A differential ideal is an ideal that is closed
under differentiation.
   A polynomial system in the state space is a system of differential equations:

                                      ẋ1 = f1 (x, u), . . . , ẋn = fn (x, u), y = h(x, u),

where h, fi ∈ R[x, u], ∀i.
  Thus, every polynomial system in the form of a state space corresponds to a differential ideal in R[x, u, y] :

                                                   I = [φ1 , . . . , φn , y − h],

where φi = ẋi − fi (x, u), i = 1, . . . , n.
  The problem of transformation from the state space to the input-output form: let I be a differential ideal;
find a generator for the differential ideal I ∩ R[u, y].

Example 4
    Suppose that it is necessary to find a differential relationship between u and y from the description in the
state space of the system:

                                        ẋ1 = −2x1 + x22 ; ẋ2 = −x1 x2 + u; y = x2 .

   Differentiating the equations of the system with respect to t and replacing ẋi by fi , we get:
   g1 = y − x2 ; g2 = ẏ − (u − x1 x2 ); g3 = ÿ − (u̇ − (x22 − 2x1 )x2 − x1 (u − x1 x2 )).
   Replace y → y0 , ẏ → y1 , ÿ → y2 , u → u0 , u̇ → u1 in gi , calculate the Gröbner basis G for:
(y0 − x1 , y1 − u0 + x1 x2 , y2 − u1 + (x22 − 2x1 )x2 + x1 (u0 − x1 x2 )) relative to lex order: u0 ≺ u1 ≺ y0 ≺ y1 ≺
y2 ≺ x1 ≺ x2 .
   Therefore, the input signals u, u̇ and the output signals y, ẏ, ÿ are related by:

                                (−2u − u̇ + 2ẏ + ÿ + y ẏ)y 3 + (−3u20 + 3uẏ − ẏ 2 )ẏ + u30 .

    In the WOLFRAM MATHEMATICA package:
    g1 = y0 − x1,
    g2 = y1 − u0 + x1 ∗ x2,
    g3 = y2 − u1 + (x22 − 2 ∗ x1) ∗ x2 + x1 ∗ (u0 − x1 ∗ x2).
    grbas = GroebnerBasis[{g1, g2, g3} , {u0, u1, y0, y1, y2, x1, x2}, {x1, x2}],
    grbas = (−2u0 − u1 + 2y1 + y2 + y0 y1 )y03 + (−3u20 + 3u0 y1 − y12 )y1 + u30 .
    




                                                                 6
Conclusion
   The paper considers methods for estimating stability using Lyapunov functions, applied to nonlinear systems.
The canonical relations of a nonlinear system are approximated by polynomials of the components of
the state and control vectors. To assess the stability, Gröbner bases are used. A method for finding the critical
points of a given nonlinear system is proposed. The coordination of input-output signals of the system based on
the construction of Gröbner bases is considered.

Acknowledgements
This work was supported by the Basic Research Program of the Siberian Branch of the Russian Academy of
Sciences No. I.5.1., Project No. 0314-2019-0020.

References
[Krasovsky59] N. N. Krasovsky Problems of the Theory of Stability of Motion Mir, Moscow, 1959

[Papach02] A. Papachristodoulou, S. Prajna On the construction of Lyapunov functions using the sum of squares
        decomposition Proceedings of the 41st IEEE Conference on Decision and Control , 3:3482-3487, 2002

[Forsman91] K. Forsman Construction of Lyapunov functions using Gröbner bases Proceedings of the 30th IEEE
        Conference on Decision and Control, 798-799, 1991

[Chukanov12] S. N. Chukanov , D. V. Ulyanov Decomposition of the vector field of control system by constructing
        a homotopy operator Probl. Upr., 6: 26, 2012

[Edelen85] D. G. B. Edelen Applied Exterior Calculus John Wiley & Sons, Inc., N.-Y., 1985
[Khalil02] H. K. Khalil Nonlinear systems Prentice Hall, 2002

[Luenberger16] D. G. Luenberger, Y. Ye Linear and nonlinear programming Springer, 2016
[Nesic02] D. Nesic, I. M. Y. Mareels, T. Glad, M. Jirstrand The Gröbner basis method in control design: an
         overview Proceedings of the 15th Triennial World Congress of IFAC, 1318, 2002
[Buchberger70] B. Buchberger Ein algorithmisches Kriterium fur die Losbarkeit eines algebraischen Gle-
        ichungssystems Aequationes Muthematicae, 4(3): 374–383, 1970
[Wolfram03] S. Wolfram The Mathematica Book Wolfram Media, Champaign, IL, 2003

[Demenkov15] M. Demenkov A Matlab tool for regions of attraction estimation via numerical algebraic geometry
        In Mechanics - Seventh Polyakhovs Reading, 2015 International Conference, 25, 2015

[Cox13] D. Cox, J. Little, D. OShea Ideals, Varieties, and Algorithms Springer, 2013
[Sidorov19] N. A. Sidorov, D. N. Sidorov, Y. Li Basins of attraction and stability of nonlinear systems equilibrium
         points Differential Equations and Dynamical Systems Springer, 110, 2019




                                                        7