Contents
MORLAB Demo: Model Reduction for Continuous-Time Second-Order Systems
This demo script contains the application of the model reduction methods for continuous-time second-order systems and some hints for proper usage.
% % This program is free software: you can redistribute it and/or modify % it under the terms of the GNU Affero General Public License as published % by the Free Software Foundation, either version 3 of the License, or % (at your option) any later version. % % This program is distributed in the hope that it will be useful, % but WITHOUT ANY WARRANTY; without even the implied warranty of % MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the % GNU Affero General Public License for more details. % % You should have received a copy of the GNU Affero General Public License % along with this program. If not, see <http://www.gnu.org/licenses/>. % % Copyright (C) 2006-2019 Peter Benner, Steffen W. R. Werner %
General Handling
The MORLAB toolbox implements structure-preserving model reduction methods for continuous-time second-order systems, i.e., for systems of the form
M*x''(t) = -K*x(t) - E*x'(t) + Bu*u(t), y(t) = Cp*x(t) + Cv*x'(t) + D*u(t).
In practice, those systems are described by a large number of differential equations, which makes the usage in optimization and controller design difficult. The goal is to construct an approximation given by
Mr*xr''(t) = -Kr*xr(t) - Er*xr'(t) + Bur*u(t), yr(t) = Cpr*xr(t) + Cvr*xr'(t) + Dr*u(t),
which is described by less differential equations than the original system. The new system has to approximate the input-output behavior of the original system, i.e., the same inputs u(t) the outputs are close in a certain sense |y(t) - yr(t)| < tol.
To get started with the model reduction in MORLAB, we load a prepared data file containing a general stable second-order system.
if exist('OCTAVE_VERSION', 'builtin') orig_warn = warning('off', 'Octave:data-file-in-path'); load morlab_data_so_stab.mat; warning(orig_warn); else load morlab_data_so_stab.mat; end
For the application of the model reduction methods, we just need to give the system to the model reduction function we want to use. This can be done in different ways. First, MORLAB supports the single matrices given to the function. As example we want to apply the second-order balanced truncation method to the system matrices.
[Mr, Er, Kr, Bur, Cpr, Cvr, Dr, info] = ml_ct_soss_bt(M, E, K, Bu, Cp, Cv, D);
The output matrices Ar, Br, Cr, Dr are now the system matrices of the approximation. The info struct contains general information about the underlying used methods.
disp(info);
Hsv: [1x1 struct] infoLYAPDL: [1x1 struct] N: 6 T: [] W: []
A more advanced way of using the model reduction function is to construct first the system as a struct. This is recommended, since most of the evaluation and analysis tools in MORLAB assume the system to be in a struct.
sys = struct( ... 'M' , M, ... 'E' , E, ... 'K' , K, ... 'Bu', Bu, ... 'Cp', Cp, ... 'Cv', Cv, ... 'D' , D);
Now we want to apply the frequency-limited second-order balanced truncation method to this system.
rom = ml_ct_soss_flbt(sys);
The reduced-order model is now given in the rom struct with the same naming of the matrices as in sys.
disp(rom);
M: [6x6 double] E: [6x6 double] K: [6x6 double] Bu: [6x3 double] Cp: [2x6 double] Cv: [2x6 double] D: [2x3 double]
To take a closer look on the result, we can use the analysis tools of MORLAB, e.g., let's compute the transfer functions of the original system and the reduction in the frequency domain by the sigma plot function. Note that it is not necessary to compute the first-order realization of the system, since the methods in MORLAB can handle the second-order structure.
ml_sigmaplot(sys, rom);

Or we can compare the two systems in a time simulation.
ml_ct_soss_simulate_ss22(sys, rom);

For creating other plots later, we are closing the subplot above again.
close;
All model reduction methods allow customization by optional parameters. option contructor. Here we will just call the function for simplicity and adjust the parameters as we want to.
opts = ml_morlabopts('ml_ct_soss_bt');
Now opts is an empty option struct with all the fields that the balanced truncation routine accepts. Note that MORLAB functions doesn't need the full option struct in general, non-existent or empty fields are replaced by default values.
disp(opts);
BalanceType: [] lyapdlopts: [1x1 struct] Method: [] Order: [] OrderComputation: [] StoreProjection: [] Tolerance: []
As a simple start, we want to compute a reduced-order model of order 5. Therefor, we are arranging the opts struct in the following way.
opts.Order = 5;
opts.OrderComputation = 'Order';
Now we can start the model reduction method again but with opts as additional input.
rom = ml_ct_soss_bt(sys, opts); disp(rom);
M: [5x5 double] E: [5x5 double] K: [5x5 double] Bu: [5x3 double] Cp: [2x5 double] Cv: [2x5 double] D: [2x3 double]
And we see that rom now contains the matrices of an order 5 system.
Special Features: Vectorization and Projection Matrices
An important feature of many model reduction methods in MORLAB is the vectorization, i.e., one can compute a bunch of reduced-order models efficiently by one call of the model reduction routine.
As example, we want to compute 3 different reduced-order models using the second-order balanced truncation method. Again we have to adjust the option struct to tell the method how the different reduced-order models should be constructed. For simplicity, we just set 3 different orders.
opts = ml_morlabopts('ml_ct_soss_bt'); opts.OrderComputation = 'Order'; opts.Order = {1, 3, 5}; rom = ml_ct_soss_bt(sys, opts); disp(rom);
[1x1 struct] [1x1 struct] [1x1 struct]
rom is now a cell array with the 3 entries corresponding to our chosen order.
disp(rom{1}); disp(rom{2}); disp(rom{3});
M: 1 E: 5.2158 K: 102.7165 Bu: [-16.1954 10.0263 -15.4432] Cp: [2x1 double] Cv: [2x1 double] D: [2x3 double] M: [3x3 double] E: [3x3 double] K: [3x3 double] Bu: [3x3 double] Cp: [2x3 double] Cv: [2x3 double] D: [2x3 double] M: [5x5 double] E: [5x5 double] K: [5x5 double] Bu: [5x3 double] Cp: [2x5 double] Cv: [2x5 double] D: [2x3 double]
Note that the evaluation tools also allow cell arrays of structs and ss objects to be evaluated. So we could just call
ml_sigmaplot(rom);

To see the sigma plots of all 3 computed roms.
The second feature we want to emphasize here is the possibility to get the projection matrices. This feature can be activated for projection-based methods by the optional parameter StoreProjection.
As example, we want to get the projection matrices for using the second-order frequency-limited balanced truncation method. Since we only need this one field, we construct the option struct this time by hand. Remember that the non-existent fields in the struct are replaced by default values. Also we change the default balancing method to position-velocity balancing, since the default second-order balancing is not a second-order projection method.
opts = struct('StoreProjection', 1, 'BalanceType', 'pv'); [rom, info] = ml_ct_soss_flbt(sys, opts); disp(rom); disp(info);
M: [7x7 double] E: [7x7 double] K: [7x7 double] Bu: [7x3 double] Cp: [2x7 double] Cv: [2x7 double] D: [2x3 double] Hsv: [36x1 double] infoLYAPDL: [1x1 struct] N: 7 T: [100x7 double] W: [7x100 double]
Info now contains the matrices T and W, that can be used to contruct the reduced-order model from the original data:
rom.M = info.W * sys.M * info.T rom.E = info.W * sys.E * info.T rom.K = info.W * sys.K * info.T rom.Bu = info.W * sys.Bu rom.Cp = sys.Cp * info.T rom.Cv = sys.Cv * info.T rom.D = sys.D
Remarks
- All model reduction routines are more or less equivalent to each other, which means that all usage examples here also apply to the other routines.
See Also