Model Reduction Using Schur and Square-Root Methods

Make sure the application is loaded.

<< NumericalControlMethods`

Load the collection of test examples.

<< NumericalControlMethods`Examples`

This experiment uses the state-space model of a servomechanism.

ss = ServoJ100[];

Here are the numbers of states, inputs, and outputs of the system.

{CountStates[ss], CountInputs[ss], CountOutputs[ss]}

Here are the ratios of each Hankel singular value to the largest one. The Hankel singular values are computed using the Cholesky factors of the controllability and observability Gramians.

hsv = SingularValues[
            CholeskyFactorControllabilityGramian[ss] . 
              Transpose[
                CholeskyFactorObservabilityGramian[
                  ss]]][[2]];

(hsv/Max)[hsv]

Selecting the threshold (10^(-3)) in the ratios retains only the first five strongest modes. This is the corresponding reduced-order model obtained by the Schur method.

NumberForm[
      rSD = DominantSubsystem[ss, ReductionMethod -> SchurDecomposition, 
          RejectionLevel -> (10^(-3))], 3]

The reduced-order model obtained by the Schur method is not internally balanced. The controllability Gramian is not diagonal.

ControllabilityGramian[rSD]

Here is the reduced-order model of the servomechanism system obtained by the square-root method.

NumberForm[
      rSR = DominantSubsystem[ss, ReductionMethod -> SquareRoot, 
          RejectionLevel -> (10^(-3))], 3]

The reduced-order model obtained by the square-root method is internally balanced.

Chop[ControllabilityGramian[rSR], (10^(-7))]

Chop[ObservabilityGramian[rSR], (10^(-7))]

The frequency response from the first input of the original system (solid lines) is compared with that of each of the reduced-order models obtained by the Schur method (dotted lines) and the square-root method (thick dotted lines). The graphs show that the transfer functions of the reduced-order models in both cases match very closely at the low frequencies, while they differ considerably at the higher frequencies.

ins = {2}; outs = {1, 3}; freq = {(10^(-2)), (10^3)};
DisplayTogetherGraphicsArray[
        BodePlot[Subsystem[ss, ins, outs], freq], 
        BodePlot[Subsystem[rSR, ins, outs], freq, 
          PlotStyle -> Dashing[{ .003, 0.007}]], 
        BodePlot[Subsystem[rSD, ins, outs], freq, 
          PlotStyle -> {{Dashing[{0.002,  .025}], 
                Thickness[ .01]}}]];