Contents

function [Er, Ar, Br, Cr] = IRKA_rail(k, r, istest)
% Computes a locally H2-optimal reduced order model of order r for
% the selective cooling of Steel profiles application described in
% [1,2,3] via the tangential IRKA method.
%
% Usage: IRKA_rail(k,r,istest)
%
% Inputs:
%
% k           refinement level of the model to use
%             (0 - 5, i.e. 109 - 79841 Dofs)
%             (optional, defaults to 2, i.e. 1357 Dofs)
%
% r           desired dimension of the reduced order model
%             (optional, defaults to 10)
%
% istest      flag to determine whether this demo runs as a CI test or
%             interactive demo
%             (optional, defaults to 0, i.e. interactive demo)
%
% References:
% [1] J. Saak, Effiziente numerische Lösung eines
%     Optimalsteuerungsproblems für die Abkühlung von Stahlprofilen,
%     Diplomarbeit, Fachbereich 3/Mathematik und Informatik, Universität
%     Bremen, D-28334 Bremen (Sep. 2003).
%     https://doi.org/10.5281/zenodo.1187040
%
% [2] P. Benner, J. Saak, A semi-discretized heat transfer model for
%     optimal cooling of steel profiles, in: P. Benner, V. Mehrmann, D.
%     Sorensen (Eds.), Dimension Reduction of Large-Scale Systems, Vol. 45
%     of Lecture Notes in Computational Science and Engineering, Springer-Verlag, Berlin/Heidelberg,
%     Germany, 2005, pp. 353–356. https://doi.org/10.1007/3-540-27909-1_19
%
% [3] J. Saak, Efficient numerical solution of large scale algebraic matrix
%     equations in PDE control and model order reduction, Dissertation,
%     Technische Universität Chemnitz, Chemnitz, Germany (Jul. 2009).
%     URL http://nbn-resolving.de/urn:nbn:de:bsz:ch1-200901642
%
% [4] S. Gugercin, A. C. Antoulas, C. Beattie, H2 model reduction
%     for large-scale linear dynamical systems, SIAM J. Matrix
%     Anal. Appl. 30 (2) (2008) pp. 609–638.
%     https://doi.org/10.1137/060666123

%
% This file is part of the M-M.E.S.S. project
% (http://www.mpi-magdeburg.mpg.de/projects/mess).
% Copyright (c) 2009-2025 Jens Saak, Martin Koehler, Peter Benner and others.
% All rights reserved.
% License: BSD 2-Clause License (see COPYING)
%

Load matrix data

if nargin < 1
    k = 2;
end
eqn = mess_get_linear_rail(k);

register corresponding usfs

opts = struct();
[oper, opts] = operatormanager(opts, 'default');

collect IRKA parameters

if nargin < 2

    opts.irka.r = 20;

else

    opts.irka.r = r;

end

if nargin < 3
    istest = false;
end

opts.irka.maxiter   = 100;
opts.irka.shift_tol = 1e-3;
opts.irka.h2_tol    = 1e-6;

if istest

    opts.irka.info = 1;

else

    opts.irka.info = 3;

end

opts.irka.init = 'logspace';

Run IRKA

[Er, Ar, Br, Cr, ~, outinfo] = mess_tangential_irka(eqn, opts, oper);

if istest && isequal(outinfo.term_flag, 'maxiter')
    mess_err(opts, 'convergence', 'IRKA converged unexpectedly slow');
end
Warning: IRKA step 1 : 1 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step   1, rel. chg. shifts = 1.466613e+02 , rel. H2-norm chg. ROM = 1.000000e+00  
Warning: IRKA step 2 : 2 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step   2, rel. chg. shifts = 1.442353e+03 , rel. H2-norm chg. ROM = Inf  
IRKA step   3, rel. chg. shifts = 4.829400e+01 , rel. H2-norm chg. ROM = Inf  
Warning: IRKA step 4 : 4 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step   4, rel. chg. shifts = 5.293215e+01 , rel. H2-norm chg. ROM = Inf  
IRKA step   5, rel. chg. shifts = 4.453988e+01 , rel. H2-norm chg. ROM = 4.162903e-01  
Warning: IRKA step 6 : 2 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step   6, rel. chg. shifts = 2.502457e+01 , rel. H2-norm chg. ROM = Inf  
Warning: IRKA step 7 : 1 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step   7, rel. chg. shifts = 5.545525e+00 , rel. H2-norm chg. ROM = Inf  
Warning: IRKA step 8 : 1 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step   8, rel. chg. shifts = 6.424616e+01 , rel. H2-norm chg. ROM = Inf  
Warning: IRKA step 9 : 2 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step   9, rel. chg. shifts = 1.958368e+03 , rel. H2-norm chg. ROM = Inf  
Warning: IRKA step 10 : 1 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step  10, rel. chg. shifts = 2.945539e+02 , rel. H2-norm chg. ROM = Inf  
IRKA step  11, rel. chg. shifts = 4.064371e+00 , rel. H2-norm chg. ROM = 4.277725e-01  
IRKA step  12, rel. chg. shifts = 3.518750e+00 , rel. H2-norm chg. ROM = 8.048146e-02  
IRKA step  13, rel. chg. shifts = 4.308130e+03 , rel. H2-norm chg. ROM = 2.093892e-02  
IRKA step  14, rel. chg. shifts = 4.747452e+00 , rel. H2-norm chg. ROM = 1.509600e-02  
IRKA step  15, rel. chg. shifts = 5.807931e-01 , rel. H2-norm chg. ROM = 1.256206e-02  
IRKA step  16, rel. chg. shifts = 4.244054e-01 , rel. H2-norm chg. ROM = 1.241050e-02  
IRKA step  17, rel. chg. shifts = 3.801743e-01 , rel. H2-norm chg. ROM = 1.231327e-02  
IRKA step  18, rel. chg. shifts = 1.434028e+02 , rel. H2-norm chg. ROM = 1.259667e-02  
IRKA step  19, rel. chg. shifts = 5.883020e+01 , rel. H2-norm chg. ROM = 1.938550e-02  
Warning: IRKA step 20 : 2 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step  20, rel. chg. shifts = 1.106745e+01 , rel. H2-norm chg. ROM = Inf  
Warning: IRKA step 21 : 1 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step  21, rel. chg. shifts = 1.789179e+01 , rel. H2-norm chg. ROM = Inf  
IRKA step  22, rel. chg. shifts = 6.924738e+00 , rel. H2-norm chg. ROM = 3.673466e-02  
Warning: IRKA step 23 : 1 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step  23, rel. chg. shifts = 1.332098e+02 , rel. H2-norm chg. ROM = Inf  
IRKA step  24, rel. chg. shifts = 1.033895e+02 , rel. H2-norm chg. ROM = 3.249677e-02  
IRKA step  25, rel. chg. shifts = 4.298765e+02 , rel. H2-norm chg. ROM = 5.132754e-02  
IRKA step  26, rel. chg. shifts = 2.478160e+01 , rel. H2-norm chg. ROM = 3.670859e-02  
IRKA step  27, rel. chg. shifts = 1.309895e+00 , rel. H2-norm chg. ROM = 3.685648e-02  
Warning: IRKA step 28 : 1 non-stable reduced eigenvalues detected.
 
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/package/package.m',13)">package (line 13)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/_release/publish_demos.m',18)">publish_demos (line 18)</a>
↳ In <a href="matlab:opentoline('/matlab/R2020b/toolbox/matlab/codetools/publish.p',0)">publish (line 0)</a>
↳ In <a href="matlab:opentoline('',21)">evalmxdom (line 21)</a>
↳ In <a href="matlab:opentoline('',109)">instrumentAndRun (line 109)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/DEMOS/Rail/IRKA_rail.m',95)">IRKA_rail (line 95)</a>
↳ In <a href="matlab:opentoline('/builds/mess/mmess/mor/mess_tangential_irka.m',369)">mess_tangential_irka (line 369)</a>
IRKA step  28, rel. chg. shifts = 2.520936e+01 , rel. H2-norm chg. ROM = Inf  
IRKA step  29, rel. chg. shifts = 1.137564e+01 , rel. H2-norm chg. ROM = 4.627848e-01  
IRKA step  30, rel. chg. shifts = 1.388013e+01 , rel. H2-norm chg. ROM = 1.052623e-01  
IRKA step  31, rel. chg. shifts = 3.837583e+00 , rel. H2-norm chg. ROM = 4.553695e-02  
IRKA step  32, rel. chg. shifts = 3.282754e-01 , rel. H2-norm chg. ROM = 4.687208e-02  
IRKA step  33, rel. chg. shifts = 5.495177e+00 , rel. H2-norm chg. ROM = 3.180306e-02  
IRKA step  34, rel. chg. shifts = 5.878087e+01 , rel. H2-norm chg. ROM = 2.013201e-02  
IRKA step  35, rel. chg. shifts = 2.077123e-01 , rel. H2-norm chg. ROM = 1.134691e-02  
IRKA step  36, rel. chg. shifts = 5.032793e+00 , rel. H2-norm chg. ROM = 6.892876e-03  
IRKA step  37, rel. chg. shifts = 2.334564e+01 , rel. H2-norm chg. ROM = 3.767834e-03  
IRKA step  38, rel. chg. shifts = 3.247672e-02 , rel. H2-norm chg. ROM = 1.902138e-03  
IRKA step  39, rel. chg. shifts = 2.395999e-02 , rel. H2-norm chg. ROM = 1.321210e-03  
IRKA step  40, rel. chg. shifts = 1.021422e-02 , rel. H2-norm chg. ROM = 1.032384e-03  
IRKA step  41, rel. chg. shifts = 1.093718e-02 , rel. H2-norm chg. ROM = 8.175446e-04  
IRKA step  42, rel. chg. shifts = 7.017752e-03 , rel. H2-norm chg. ROM = 6.700431e-04  
IRKA step  43, rel. chg. shifts = 5.507644e-03 , rel. H2-norm chg. ROM = 5.561991e-04  
IRKA step  44, rel. chg. shifts = 4.064457e-03 , rel. H2-norm chg. ROM = 4.636019e-04  
IRKA step  45, rel. chg. shifts = 2.858886e-03 , rel. H2-norm chg. ROM = 3.878166e-04  
IRKA step  46, rel. chg. shifts = 2.139374e-03 , rel. H2-norm chg. ROM = 3.248186e-04  
IRKA step  47, rel. chg. shifts = 1.540009e-03 , rel. H2-norm chg. ROM = 2.725457e-04  
IRKA step  48, rel. chg. shifts = 1.346657e-03 , rel. H2-norm chg. ROM = 2.287556e-04  
IRKA step  49, rel. chg. shifts = 8.596453e-04 , rel. H2-norm chg. ROM = 1.921911e-04  
IRKA terminated due to relative change of shifts criterion.

Computing TFMs of original and reduced order systems and MOR errors
 Step  10 / 100 Step  20 / 100 Step  30 / 100 Step  40 / 100 Step  50 / 100 Step  60 / 100 Step  70 / 100 Step  80 / 100 Step  90 / 100 Step 100 / 100

In case this is a CI test and the sparss class is available (recent MATLAB)

try the same again with that.
if istest && exist('sparss', 'class')

    sys = sparss(eqn.A_, eqn.B, eqn.C, [], eqn.E_);

    [Er, Ar, Br, Cr, ~, outinfo] = mess_tangential_irka(sys, opts);

    if istest && isequal(outinfo.term_flag, 'maxiter')

        mess_err(opts, 'convergence', 'IRKA converged unexpectedly slow');

    end
end
ans =

  Columns 1 through 13

     1     0     0     0     0     0     0     0     0     0     0     0     0
     0     1     0     0     0     0     0     0     0     0     0     0     0
     0     0     1     0     0     0     0     0     0     0     0     0     0
     0     0     0     1     0     0     0     0     0     0     0     0     0
     0     0     0     0     1     0     0     0     0     0     0     0     0
     0     0     0     0     0     1     0     0     0     0     0     0     0
     0     0     0     0     0     0     1     0     0     0     0     0     0
     0     0     0     0     0     0     0     1     0     0     0     0     0
     0     0     0     0     0     0     0     0     1     0     0     0     0
     0     0     0     0     0     0     0     0     0     1     0     0     0
     0     0     0     0     0     0     0     0     0     0     1     0     0
     0     0     0     0     0     0     0     0     0     0     0     1     0
     0     0     0     0     0     0     0     0     0     0     0     0     1
     0     0     0     0     0     0     0     0     0     0     0     0     0
     0     0     0     0     0     0     0     0     0     0     0     0     0
     0     0     0     0     0     0     0     0     0     0     0     0     0
     0     0     0     0     0     0     0     0     0     0     0     0     0
     0     0     0     0     0     0     0     0     0     0     0     0     0
     0     0     0     0     0     0     0     0     0     0     0     0     0
     0     0     0     0     0     0     0     0     0     0     0     0     0

  Columns 14 through 20

     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     0     0     0     0     0     0     0
     1     0     0     0     0     0     0
     0     1     0     0     0     0     0
     0     0     1     0     0     0     0
     0     0     0     1     0     0     0
     0     0     0     0     1     0     0
     0     0     0     0     0     1     0
     0     0     0     0     0     0     1