ml_ct_ss_brbt

Bounded-real balanced truncation for standard systems.

Contents

Syntax

[Ar, Br, Cr, Dr, info] = ml_ct_ss_brbt(A, B, C, D)
[Ar, Br, Cr, Dr, info] = ml_ct_ss_brbt(A, B, C, D, opts)
[rom, info] = ml_ct_ss_brbt(sys)
[rom, info] = ml_ct_ss_brbt(sys, opts)

Description

This function computes the bounded-real balanced truncation for a bounded-real standard system of the form

   x'(t) = A*x(t) + B*u(t),                                        (1)
    y(t) = C*x(t) + D*u(t).                                        (2)

Therefore, the two bounded-real Riccati equations

   A*P  + P*A' + B*B' + (P*C' + B*D') * inv(Rb) * (P*C' + B*D')' = 0,
   A'*Q + Q*A  + C'*C + (B'*Q + D'*C)' * inv(Rc) * (B'*Q + D'*C) = 0,

are solved for the Gramians P and Q, with

   Rb = I - D*D'   and
   Rc = I - D'*D.

As result, a reduced-order bounded-real system of the form

   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 transfer function G with the spectral factors W'*W = I - G'*G and V*V' = I - G*G', and the r-th order transfer function Gr with the spectral factors Wr'*Wr = I - Gr'*Gr and Vr*Vr' = I - Gr*Gr' it holds

   max(||[G - Gr; W - Wr]||_{\infty}, ||[G - Gr; V - Vr]||_{\infty})
   <= 2 * (Hsv(r+1) + ... + Hsv(n)),

with Hsv, a vector containing the characteristic bounded-real singular values of the system.

Input

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

Parameter
Meaning
caredlopts
structure, containing the optional parameters for the Riccati equation sign function solver, only used if RiccatiSolver = 'sign', see ml_caredl_sgn_fac
default: struct()
Method
{!}
character array, determining algorithm for the computation of the reduced-order model
  • 'sr' - square-root method
  • 'bfsr' - balancing-free square-root method
default: 'sr'
Order
{!}
positive integer, order of the resulting reduced-order model chosen by the user if 'order' is set for OrderComputation
default: min(10,length(Hsv)) + Nu
OrderComputation
{!}
character array, determining the method for the computation of the size of the reduced-order model
  • 'order' - take explicit order
  • 'tolerance' - using absolute error bound
default: 'tolerance'
pcareopts
structure, containing the optional parameters for the computation of the continuous-time algebraic positive Riccati equation, only used if RiccatiSolver = 'newton', see ml_pcare_nwt_fac
default: struct()
RiccatiSolver
character array, determining the solver for the dual Riccati equations
  • 'newton' - Newton iteration
  • 'sign' - dual sign function method
default: 'sign'
StoreProjection
{0, 1}, used to disable/enable storing of the computed projection matrices W and T
default: 0
Tolerance
{!}
nonnegative scalar, tolerance used for the computation of the size of the reduced-order model by an relative error bound if 'tolerance' is set for OrderComputation
default: 1.0e-02

Note: Parameters marked with {!} may also be a cell array containing multiple arguments. In this case an cell array of the same size is returned with one entry computed for each input argument and the marked fields of the info struct are cells as well. When multiple arguments are given as cells, they are expected to have the same length.

Output

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

Entry
Meaning
AbsErrBound
{!}
computed error bound for the absolute error of the reduced-order model (and the spectral factors of I - G*G' and I - G'*G) in H-infinity norm
Hsv
a vector, containing the computed characteristic bounded-real singular values
infoCAREDL
structure, containing information about the sign function solver for the dual Riccati equations, see ml_caredl_sgn_fac
infoPCARE_C
structure, containing information about the continuous-time algebraic positive Riccati equation solver for the controllability Gramian, see ml_pcare_nwt_fac
infoPCARE_O
structure, containing information about the continuous-time algebraic positive Riccati equation solver for the observability Gramian, see ml_pcare_nwt_fac
N
{!}
Dimension of the reduced-order model
T
{!}
projection matrix used as right state-space transformation to obtain the resulting block system, if opts.StoreProjection == 1
W
{!}
projection matrix used as left state-space transformation to obtain the resulting block system, if opts.StoreProjection == 1

Reference

S. Gugercin, A. C. Antoulas, A survey of model reduction by balanced truncation and some new results, Internat. J. Control 77 (8) (2004) 748--766. https://doi.org/10.1080/00207170410001713448

See Also

ml_ct_dss_brbt | ml_morlabopts