Main Content

stabsep

Stable-unstable decomposition

    Description

    [G1,G2] = stabsep(G) decomposes the linear time-invariant (LTI) model G into its stable and unstable parts.

    G=G1+G2

    G1 contains all stable modes (poles) that can be separated from the unstable modes in a numerically stable way, and G2 contains the remaining modes. G1 is always proper and contains the feedthrough term (if any).

    example

    [G1,G2] = stabsep(G,Offset) shifts the stability boundary so that G1 includes all poles satisfying: (since R2023b)

    • Re(s) < -Offset × max(1,|Im(s)|) (continuous time)

    • |z| < 1 - Offset (discrete time)

    example

    [G1,G2] = stabsep(G,___,Name=Value) specifies additional options for controlling accuracy or focusing on unstable modes. (since R2023b)

    [G1,G2,info] = stabsep(___) returns a structure info containing the block-diagonalizing transformation matrices. (since R2023b)

    Examples

    collapse all

    This example shows how to separate a zero-pole-gain model into its stable and unstable parts.

    Consider the following system.

    H(s)=0.1(s+2)(s+1)(s-1)(s+0.001)

    The system contains three stable and one unstable poles. One stable pole is really close to the stability boundary.

    Create a zpk model of this system and obtain the stable-unstable decomposition.

    h = zpk(1,[-2 -1 1 -0.001],0.1);
    [h1,h2] = stabsep(h)
    h1 =
     
               0.1
      ---------------------
      (s+0.001) (s+1) (s+2)
     
    Continuous-time zero/pole/gain model.
    
    h2 =
     
      -2.8326e-18
      -----------
         (s-1)
     
    Continuous-time zero/pole/gain model.
    

    The stable part h1 contains all the stable poles. You can use the Offset argument to shift the stability boundary and exclude the nominally stable pole (s = –0.001) from the stable part of the decomposition.

    Compute the decomposition with the shifted boundary.

    [h1,h2] = stabsep(h,0.1)
    h1 =
     
      -0.050075 (s+2.999)
      -------------------
          (s+1) (s+2)
     
    Continuous-time zero/pole/gain model.
    
    h2 =
     
      0.050075 (s-1)
      ---------------
      (s+0.001) (s-1)
     
    Continuous-time zero/pole/gain model.
    

    For this decomposition, the nominally stable pole is in the unstable part h2.

    By default, the function provides the stable part in the first output argument and the unstable part in the second output argument. You can use the Focus option to reverse the order.

    [h1,h2] = stabsep(h,Focus="unstable")
    h1 =
     
      -2.1943e-18
      -----------
         (s-1)
     
    Continuous-time zero/pole/gain model.
    
    h2 =
     
               0.1
      ---------------------
      (s+0.001) (s+1) (s+2)
     
    Continuous-time zero/pole/gain model.
    

    For this decomposition, h1 contains the unstable pole and the nominally stable pole.

    Input Arguments

    collapse all

    Dynamic system to decompose, specified as a numeric LTI model, such as a ss or tf model.

    Since R2023b

    Offset for stable/unstable boundary, specified as a scalar. In the stable/unstable decomposition, the stable term includes only poles satisfying:

    • Re(s) < -Offset × max(1,|Im(s)|) (continuous time)

    • |z| < 1 - Offset (discrete time)

    Increase the value of Offset to exclude poles close to the stability boundary from G1.

    Name-Value Arguments

    Specify optional pairs of arguments as Name1=Value1,...,NameN=ValueN, where Name is the argument name and Value is the corresponding value. Name-value arguments must appear after other arguments, but the order of the pairs does not matter.

    Example: [G1,G2] = stabsep(G,0.1,SepTol=100)

    Since R2023b

    Focus of decomposition, specified as either "stable" or "unstable".

    • "stable" — The first argument G1 includes all poles satisfying

      • Re(s) < -Offset × max(1,|Im(s)|) (continuous time)

      • |z| < 1 - Offset (discrete time)

    • "unstable" — The first argument G1 includes all poles satisfying

      • Re(s) > Offset × max(1,|Im(s)|) (continuous time)

      • |z| > 1 + Offset (discrete time)

    When an exact stable/unstable split is numerically impossible, G1 includes a maximum subset of the desired modes and G2 contains the remaining mix of stable and unstable modes.

    Since R2023b

    Accuracy loss factor for stable/unstable decomposition, specified as a positive scalar value. When decomposing a model G(s), stabsep ensures that the frequency responses of G and G1 + G2 differ by no more than SepTol times the absolute accuracy of the computed value of G(s). Increasing SepTol helps separate nearby stable and unstable modes at the expense of accuracy.

    Output Arguments

    collapse all

    Stable part of the decomposition, returned as a model of same type as G.

    When Focus is set to "unstable", the function returns the unstable part of the decomposition in the G1 argument.

    Unstable part of the decomposition, returned as a model of same type as G.

    When Focus is set to "unstable", the function returns the stable part of the decomposition in the G2 argument.

    Since R2023b

    Additional information about the decomposition, returned as structure with these fields.

    FieldDescription
    TLLeft-side matrix of the block-diagonalizing transformation, returned as a matrix with dimensions Nx-by-Nx, where Nx is the number of states in the model G.
    TRRight-side matrix of the block-diagonalizing transformation, returned as a matrix with dimensions Nx-by-Nx, where Nx is the number of states in the model G.

    The algorithm transforms the state-space realization (A, B, C, D) of a model to

    TLATR=(A100A2),TLB=(B1B2),CTR=(C1C2)

    The transformation scales the system and separates the stable and unstable parts. All modes of A1 are unstable and all modes of A2 are stable.

    The function returns an empty value [] for this argument when the input model G is not a state-space model.

    References

    [1] Bartels, R. H., and G. W. Stewart. “Algorithm 432 [C2]: Solution of the Matrix Equation AX + XB = C [F4].” Communications of the ACM 15, no. 9 (September 1972): 820–26. https://doi.org/10.1145/361573.361582.

    [2] Bavely, Connice A., and G. W. Stewart. “An Algorithm for Computing Reducing Subspaces by Block Diagonalization.” SIAM Journal on Numerical Analysis 16, no. 2 (April 1979): 359–67. https://doi.org/10.1137/0716028.

    [3] Demmel, James. “The Condition Number of Equivalence Transformations That Block Diagonalize Matrix Pencils.” SIAM Journal on Numerical Analysis 20, no. 3 (June 1983): 599–610. https://doi.org/10.1137/0720040.

    Version History

    Introduced before R2006a

    expand all

    See Also

    |