ml_ct_dss_bst
Balanced stochastic truncation for descriptor systems.
Contents
Syntax
[Ar, Br, Cr, Dr, Er, info] = ml_ct_dss_bst(A, B, C, D, E) [Ar, Br, Cr, Dr, Er, info] = ml_ct_dss_bst(A, B, C, D, E, opts)
[rom, info] = ml_ct_dss_bst(sys) [rom, info] = ml_ct_dss_bst(sys, opts)
Description
This function computes the generalized balanced stochastic truncation for a descriptor system of the form
E*x'(t) = A*x(t) + B*u(t), (1) y(t) = C*x(t) + D*u(t), (2)
Therefore, first an additive decomposition of the system is performed using the matrix disk function, such that
[ Ei 0 ] [ Ai 0 ] [ Ci ] E2 = [ ], A2 = [ ], B2 = [ Bi, Bp ], C2 = [ ], [ 0 Ep ] [ 0 Ap ] [ Cp ]
with (Ei, Ai, Bi, Ci, D) belonging to the polynomial part and (Ep, Ap, Bp, Cp, 0) belonging to the strictly proper part. Now, for the reduction of the strictly proper part the generalized continuous-time Lyapunov equation
Ap*Pp*Ep' + Ep*Pp*Ap' + Bp*Bp' = 0
is solved for Pp and then, the corresponding generalized Riccati equation
Ap'*Qp*Ep + Ep'*Qp*Ap + (Cp - Bw' * Qp * Ep)' * inv(M*M') * (Cp - Bw' * Qp * Ep) = 0
is solved for Qp, with M = D - Ci * inv(Ai) * Bi. Also, the the two generalized discrete-time Lyapunov equations
Ai*Pi*Ai' - Ei*Pi*Ei' - Bi*Bi' = 0, Ai'*Qi*Ai - Ei'*Qi*Ei - Ci'*Ci = 0
are solved for the reduction of the polynomial part. As result, a reduced-order system of the form
Er*x'(t) = Ar*x(t) + Br*u(t), (3) y(t) = Cr*x(t) + Dr*u(t) (4)
is computed, such that for the original proper transfer function G and the reduced-order transfer function Gr with an r-th order strictly proper part it holds
||G - Gr||_{\infty} / ||G||_{\infty} <= ((1 + Hsvp(r+1))/(1 - Hsvp(r+1)) * ... * (1 + Hsvp(n))/(1 - Hsvp(n))) + 1,
with Hsvp, a vector containing the proper characteristic stochastic singular values of the system. If the transfer function is invertible it holds
||inv(G)*(G - Gr)||_{\infty} <= ((1 + Hsvp(r+1))/(1 - Hsvp(r+1)) * ... * (1 + Hsvp(n))/(1 - Hsvp(n))) + 1.
Notes: 1) The equations above are defined for the case of p < m. Otherwise the system is transposed, then reduced and back transposed. 2) In case of a rank-deficient M term, an epsilon regularization is performed, which replaces the M during the computations with an identity matrix scaled by a given epsilon. 3) For unstable systems, first an additive decomposition into the stable and anti-stable parts is performed and then only the stable part will be reduced.
Input
A - matrix from (1) with dimensions n x n B - matrix from (1) with dimensions n x m C - matrix from (2) with dimensions p x n D - matrix from (2) with dimensions p x m E - matrix from (1) with dimensions n x n sys - structure or state-space object, containing the descriptor system's matrices:
Entry | Meaning |
A | matrix from (1) with dimensions n x n |
B | matrix from (1) with dimensions n x m |
C | matrix from (2) with dimensions p x n |
D | matrix from (2) with dimensions p x m |
E | matrix from (1) with dimensions n x n |
opts - structure, containing the following optional entries:
Parameter | Meaning |
DecompEig | positive scalar, overestimation of the absolute value of the largest finite eigenvalue of s*E - A, if set, replaces the computation with DecompTol default: [] |
DecompTol | nonnegative scalar, tolerance multiplied with the largest singular value of E to determine the smallest non-quasi-zero singular value of E default: log(n)*eps |
Epsilon | positive scalar, used in the case of a non-full-rank M + M' term for epsilon regularization by multiplying with an identity matrix of appropriate size default: 1.0e-03 |
gdlyapopts | structure, containing the optional parameters for the computation of the generalized discrete-time Lyapunov equations, see ml_gdlyap_smith_fac default: struct() |
ImproperTrunc | nonnegative scalar, tolerance multiplied with the largest proper Hankel singular value of the system to truncate the improper part default: log(n)*eps |
Index | nonnegative integer, index of the descriptor system used to set an upper bound on the size of the reduced improper part, Inf if unknown default: Inf |
infdecopts | structure, containing the optional parameters for the decomposition of the finite and infinite parts of the system using the disk function and subspace extraction method, see ml_disk and ml_getqz default: struct() |
lyapopts | structure, containing the optional parameters for the computation of the continuous-time algebraic Lyapunov equation, see ml_lyap_sgn_fac default: struct() |
Method | character array, determining algorithm for the computation of the reduced-order model
|
Order | positive integer, order of the resulting reduced-order model chosen by the user if 'order' is set for OrderComputation default: min(10,length(Hsvp)) + Nu + Ni |
OrderComputation | character array, determining the method for the computation of the size of the reduced-order model
|
pcareopts | structure, containing the optional parameters for the computation of the continuous-time algebraic positive Riccati equation, see ml_pcare_nwt_fac default: struct() |
stabdecopts | structure, containing the optional parameters for the decomposition of the stable and unstable parts of the system using the disk function and subspace extraction method, see ml_disk and ml_getqz default: struct() |
Tolerance | nonnegative scalar, tolerance used for the computation of the size of the reduced-order model by an absolute error bound if 'tolerance' is set for OrderComputation default: 1.0e-02 |
Output
Ar - matrix of (3) with dimensions r x r Br - matrix of (3) with dimensions r x m Cr - matrix of (4) with dimensions p x r Dr - matrix of (4) with dimensions p x m Er - matrix of (3) with dimensions r x r rom - structure or state-space object, containing the reduced-order descriptor system:
Entry | Meaning |
A | matrix from (3) with dimensions r x r |
B | matrix from (3) with dimensions r x m |
C | matrix from (4) with dimensions p x r |
D | matrix from (4) with dimensions p x m |
E | matrix from (3) with dimensions r x r |
info - structure, containing the following information:
Entry | Meaning |
AbsErrBound | computed error bound for the absolute error of the reduced-order model in H-infinity norm |
Hsvi | a vector, containing the computed Hankel singular values of the improper part of the system |
Hsvp | a vector, containing the computed characteristic stochastic singular values of the proper part of the system |
infoGADTF | structure, containing information about the additive decomposition of the system into its infinite, finite stable and finite anti-stable parts, see ml_ct_dss_adtf |
infoGDLYAP_C | structure, containing information about the generalized discrete-time Lyapunov equation solver for the improper controllability Gramian, see ml_gdlyap_smith_fac |
infoGDLYAP_O | structure, containing information about the generalized discrete-time Lyapunov equation solver for the improper observability Gramian, see ml_gdlyap_smith_fac |
infoLYAP | structure, containing information about the continuous-time Lyapunov equation sovler for the controllability Gramian, ml_lyap_sgn_fac |
infoPCARE | structure, containing information about the continuous-time algebraic positive Riccati equation for the observability Gramian, see ml_pcare_nwt_fac |
Ni | Dimension of the improper part in the reduced- order model |
Np | Dimension of the proper part in the reduced-order model |
Nu | Dimension of the unstable part in the reduced- order model |
RelErrBound | computed error bound for the relative error of the of the reduced-order model in H-infinity norm |
Reference
P. Benner, T. Stykel, Model order reduction for differential-algebraic equations: A survey, in: A. Ilchmann, T. Reis (Eds.), Surveys in Differential-Algebraic Equations IV, Differential-Algebraic Equations Forum, Springer International Publishing, Cham, 2017, pp. 107--160. https://doi.org/10.1007/978-3-319-46618-7_3
See Also