Simplification of descriptor system models

sysr = gminreal(sys; contr = true, obs = true, noseig = true, prescale, fast = true, 
                atol = 0, atol1 = atol, atol2 = atol, rtol = nϵ)

Compute for a descriptor system sys = (A-λE,B,C,D) of order n a reduced order descriptor system sysr = (Ar-λEr,Br,Cr,Dr) of order nr ≤ n such that sys and sysr have the same transfer function matrix, i.e.,

         -1                    -1
 C*(λE-A)  *B + D = Cr*(λEr-Ar)  *Br + Dr .

If prescale = true, a preliminary balancing of the descriptor system matrices is performed. The default setting is prescale = gbalqual(sys) > 10000, where gbalqual(sys) is the scaling quality of the descriptor system model sys (see gbalqual).

The least possible order nr is achieved if contr = true, obs = true and nseig = true. Such a realization is called minimal and satisfies:

 (1) rank[Br Ar-λEr] = nr for all finite λ (finite controllability)

 (2) rank[Br Er] = nr (infinite controllability)

 (3) rank[Ar-λEr; Cr] = nr for all finite λ (finite observability)

 (4) rank[Er; Cr] = nr (infinite observability)

 (5) Ar-λEr has no simple infinite eigenvalues

A realization satisfying only conditions (1)-(4) is called irreducible.

Some reduction steps can be skipped by appropriately selecting the keyword arguments contr, obs and nseig.

If contr = false, then the controllability conditions (1) and (2) are not enforced.

If obs = false, then observability condition (3) and (4) are not enforced.

If nseig = false, then condition (5) on the lack of simple infinite eigenvalues is not enforced.

To enforce conditions (1)-(4), orthogonal similarity transformations are performed on the matrices of the original realization (A-λE,B,C,D) to obtain an irreducible realization using structured pencil reduction algorithms, as the fast versions of the reduction techniques of the full row rank pencil [B A-λE] and full column rank pencil [A-λE;C] proposed in [1]. To enforce condition (5), residualization formulas (see, e.g., [2, page 329]) are employed which involves matrix inversions.

The underlying pencil manipulation algorithms employ rank determinations based on either the use of rank revealing QR-decomposition with column pivoting, if fast = true, or the SVD-decomposition. The rank decision based on the SVD-decomposition is generally more reliable, but the involved computational effort is higher.

The keyword arguments atol1, atol2, and rtol, specify, respectively, the absolute tolerance for the nonzero elements of matrices A, B, C, D, the absolute tolerance for the nonzero elements of E, and the relative tolerance for the nonzero elements of A, B, C, D and E. The default relative tolerance is nϵ, where ϵ is the working machine epsilon and n is the order of the system sys. The keyword argument atol can be used to simultaneously set atol1 = atol and atol2 = atol.

[1] P. Van Dooreen, The generalized eigenstructure problem in linear system theory, IEEE Transactions on Automatic Control, vol. AC-26, pp. 111-129, 1981.

[2] A. Varga, Solving Fault Diagnosis Problems - Linear Synthesis Techniques, Springer Verlag, 2017.

sysr = gir(sys; finite = true, infinite = true, contr = true, obs = true, noseig = false,
           fast = true, atol = 0, atol1 = atol, atol2 = atol, rtol = nϵ)

Compute for a descriptor system sys = (A-λE,B,C,D) of order n a reduced order descriptor system sysr = (Ar-λEr,Br,Cr,Dr) of order nr ≤ n such that sys and sysr have the same transfer function matrix, i.e.,

         -1                    -1
 C*(λE-A)  *B + D = Cr*(λEr-Ar)  *Br + Dr .

The least possible order nr is achieved if finite = true, infinite = true, contr = true, obs = true and nseig = true. Such a realization is called minimal and satisfies:

 (1) rank[Br Ar-λEr] = nr for all finite λ (finite controllability)

 (2) rank[Br Er] = nr (infinite controllability)

 (3) rank[Ar-λEr; Cr] = nr for all finite λ (finite observability)

 (4) rank[Er; Cr] = nr (infinite observability)

 (5) Ar-λEr has no simple infinite eigenvalues

A realization satisfying only conditions (1)-(4) is called irreducible and is computed by default.

Some reduction steps can be skipped by appropriately selecting the keyword arguments contr, obs, finite, infinite and nseig.

If contr = false, then the controllability conditions (1) and (2) are not enforced. If contr = true and finite = true, then the finite controllability condition (1) is enforced. If contr = true and finite = false, then the finite controllability condition (1) is not enforced. If contr = true and infinite = true, then the infinite controllability condition (2) is enforced. If contr = true and infinite = false, then the infinite controllability condition (2) is not enforced.

If obs = false, then observability condition (3) and (4) are not enforced. If obs = true and finite = true, then the finite observability condition (3) is enforced. If obs = true and finite = false, then the finite observability condition (3) is not enforced. If obs = true and infinite = true, then the infinite observability condition (4) is enforced. If obs = true and infinite = false, then the infinite observability condition (4) is not enforced.

If nseig = true, then condition (5) on the lack of simple infinite eigenvalues is also enforced.

To enforce conditions (1)-(4), the Procedure GIR in [1, page 328] is employed, which performs orthogonal similarity transformations on the matrices of the original realization (A-λE,B,C,D) to obtain an irreducible realization using structured pencil reduction algorithms. To enforce condition (5), residualization formulas (see, e.g., [1, page 329]) are employed which involves matrix inversions.

The underlying pencil manipulation algorithms employ rank determinations based on either the use of rank revealing QR-decomposition with column pivoting, if fast = true, or the SVD-decomposition. The rank decision based on the SVD-decomposition is generally more reliable, but the involved computational effort is higher.

The keyword arguments atol1, atol2, and rtol, specify, respectively, the absolute tolerance for the nonzero elements of matrices A, B, C, D, the absolute tolerance for the nonzero elements of E, and the relative tolerance for the nonzero elements of A, B, C, D and E. The default relative tolerance is nϵ, where ϵ is the working machine epsilon and n is the order of the system sys. The keyword argument atol can be used to simultaneously set atol1 = atol and atol2 = atol.

[1] A. Varga, Solving Fault Diagnosis Problems - Linear Synthesis Techniques, Springer Verlag, 2017.

gir_lrtran(sys; ltran = false, rtran = false, finite = true, infinite = true, contr = true, obs = true, 
         noseig = false, fast = true, atol = 0, atol1 = atol, atol2 = atol, rtol = nϵ) -> (sysr, Q, Z)

This is a special version of the function gir to additionally determine the left and right transformation matrices Q = [Q1 Q2] and Z = [Z1 Z2], respectively, such that the matrices Ar, Er, Br and Cr of the resulting descriptor system sysr = (Ar-λEr,Br,Cr,Dr) are given by Ar = Q1'*A*Z1, Er = Q1'*E*Z1, Br = Q1'*B, Cr = C*Z1, where the number of columns of Q1 and Z1 is equal to the order of matrix Ar. Q and Z result orthogonal if noseig = false. Q = nothing if ltran = false and Z = nothing if rtran = false. See gir for details on the rest of keyword parameters.

gbalmr(sys, balance = false, matchdc = false, ord = missing, offset = √ϵ,
       atolhsv = 0, rtolhsv = nϵ, atolmin = atolhsv, rtolmin = rtolhsv, 
       atol = 0, atol1 = atol, atol2 = atol, rtol, fast = true) -> (sysr, hs)

Compute for a proper and stable descriptor system sys = (A-λE,B,C,D) with the transfer function matrix G(λ), a reduced order realization sysr = (Ar-λEr,Br,Cr,Dr) and the vector hs of decreasingly ordered Hankel singular values of the system sys. If balance = true, a balancing-based approach is used to determine a reduced order minimal realization of the form sysr = (Ar-λI,Br,Cr,Dr). For a continuous-time system sys, the resulting realization sysr is balanced, i.e., the controllability and observability grammians are equal and diagonal. If additonally matchdc = true, the resulting sysr is computed using state rezidualization formulas (also known as singular perturbation approximation) which additionally preserves the DC-gain of sys. In this case, the resulting realization sysr is balanced (for both continuous- and discrete-time systems). If balance = false, an enhanced accuracy balancing-free approach is used to determine the reduced order system sysr.

If ord = nr, the resulting order of sysr is min(nr,nrmin), where nrmin is the order of a minimal realization of sys determined as the number of Hankel singular values exceeding max(atolmin,rtolmin*HN), with HN, the Hankel norm of G(λ). If ord = missing, the resulting order is chosen as the number of Hankel singular values exceeding max(atolhsv,rtolhsv*HN).

To check the stability of the eigenvalues of the pencil A-λE, the real parts of eigenvalues must be less than -β for a continuous-time system or the moduli of eigenvalues must be less than 1-β for a discrete-time system, where β is the stability domain boundary offset. The offset β to be used can be specified via the keyword parameter offset = β. The default value used for β is sqrt(ϵ), where ϵ is the working machine precision.

The keyword arguments atol1, atol2, and rtol, specify, respectively, the absolute tolerance for the nonzero elements of A, B, C, D, the absolute tolerance for the nonzero elements of E, and the relative tolerance for the nonzero elements of A, B, C, D and E. The default relative tolerance is nϵ, where ϵ is the working machine epsilon and n is the order of the system sys. The keyword argument atol can be used to simultaneously set atol1 = atol and atol2 = atol.

If E is singular, the uncontrollable infinite eigenvalues of the pair (A,E) and the non-dynamic modes are elliminated using minimal realization techniques. The rank determinations in the performed reductions are based on rank revealing QR-decompositions with column pivoting if fast = true or the more reliable SVD-decompositions if fast = false.

Method: For the order reduction of a standard system, the balancing-free method of [1] or the balancing-based method of [2] are used. For a descriptor system the balancing related order reduction methods of [3] are used. To preserve the DC-gain of the original system, the singular perturbation approximation method of [4] is used in conjunction with the balancing-based or balancing-free approach of [5].

References

[1] A. Varga. Efficient minimal realization procedure based on balancing. In A. El Moudni, P. Borne, and S.G. Tzafestas (Eds.), Prepr. of the IMACS Symp. on Modelling and Control of Technological Systems, Lille, France, vol. 2, pp.42-47, 1991.

[2] M. S. Tombs and I. Postlethwaite. Truncated balanced realization of a stable non-minimal state-space system. Int. J. Control, vol. 46, pp. 1319–1330, 1987.

[3] T. Stykel. Gramian based model reduction for descriptor systems. Mathematics of Control, Signals, and Systems, 16:297–319, 2004.

[4] Y. Liu Y. and B.D.O. Anderson Singular Perturbation Approximation of Balanced Systems, Int. J. Control, Vol. 50, pp. 1379-1405, 1989.

[5] Varga A. Balancing-free square-root algorithm for computing singular perturbation approximations. Proc. 30-th IEEE CDC, Brighton, Dec. 11-13, 1991, Vol. 2, pp. 1062-1065.

gss2ss(sys; Eshape = "ident", atol = 0, atol1 = atol, atol2 = atol, rtol = nϵ) -> (sysr, r)

Convert the descriptor system sys = (A-λE,B,C,D) to an input-output equivalent descriptor system realization sysr = (Ar-λEr,Br,Cr,Dr) without non-dynamic modes and having the same transfer function matrix. The resulting Er is in the SVD-like form Er = blockdiag(E1,0), with E1 an r × r nonsingular matrix, where r is the rank of E.

The keyword argument Eshape specifies the shape of E1 as follows:

if Eshape = "ident" (the default option), E1 is an identity matrix of order r and if E is nonsingular, then the resulting system sysr is a standard state-space system;

if Eshape = "diag", E1 is a diagonal matrix of order r, where the diagonal elements are the decreasingly ordered nonzero singular values of E;

if Eshape = "triu", E1 is an upper triangular nonsingular matrix of order r.

The keyword arguments atol1, atol2, and rtol, specify, respectively, the absolute tolerance for the nonzero elements of A, the absolute tolerance for the nonzero elements of E, and the relative tolerance for the nonzero elements of A and E. The keyword argument atol can be used to simultaneously set atol1 = atol and atol2 = atol.

If Eshape = "triu", the reductions of E and A are performed using rank decisions based on rank revealing QR-decompositions with column pivoting. If Eshape = "ident" or Eshape = "diag" the reductions are performed using the more reliable SVD-decompositions.

 dss2ss(sys[, x0 = 0]; state-mapping = false, simple_infeigs = true, fast = true, atol1, atol2, rtol) 
           -> (sysr, xr0, Mx, Mu)

Return for a proper descriptor system sys = (A-λE,B,C,D) and initial state x0, the equivalent reduced order standard system sysr = (Ar-λI,Br,Cr,Dr) and the corresponding reduced consistent initial state xr0.

If state_mapping = true, the state mapping matrices Mx and Mu are also determined such that the values x(t) and xr(t) of the state vectors of the systems sys and sysr, respectively, and the input vector u(t) are related as x(t) = Mx*xr(t)+Mu*u(t). In this case, higher order uncontrollable infinite eigenvalues can be eliminated if simple_infeigs = false.

By default, state_mapping = false and Mx = nothing and Mu = nothing. In this case, higher order uncontrollable or unobservable infinite eigenvalues can be eliminated if simple_infeigs = false.

By default, simple_infeigs = true, and simple infinite eigenvalues for the pair (A,E) are assumed and eliminated.

The underlying pencil manipulation algorithms employ rank determinations based on either the use of rank revealing QR-decomposition with column pivoting, if fast = true (default), or the SVD-decomposition, if fast = false. The rank decision based on the SVD-decomposition is generally more reliable, but the involved computational effort is higher.

The keyword arguments atol1, atol2 and rtol specify, respectively, the absolute tolerance for the nonzero elements of A, B, C, D, the absolute tolerance for the nonzero elements of E, and the relative tolerance for the nonzero elements of A, B, C, D and E. The default relative tolerance is n*ϵ, where n is the order of the square matrices A and E, and ϵ is the working machine epsilon.