Anonymous
×
Create a new article
Write your page title here:
We currently have 105 articles on MOR Wiki. Type your article name above or click on one of the titles below and start writing!



MOR Wiki

Difference between revisions of "Inverse Lyapunov Procedure"

(corrections 1)
 
(13 intermediate revisions by 3 users not shown)
Line 1: Line 1:
 
[[Category:benchmark]]
 
[[Category:benchmark]]
  +
[[Category:procedural]]
 
[[Category:linear]]
 
[[Category:linear]]
 
[[Category:first differential order]]
 
[[Category:first differential order]]
 
[[Category:time invariant]]
 
[[Category:time invariant]]
  +
[[Category:MIMO]]
   
 
==Description==
 
==Description==
The '''Inverse Lyapunov Procedure''' (ilp) is a synthetic random linear system generator.
+
The '''Inverse Lyapunov Procedure''' (ILP) is a synthetic random linear system generator.
It is based on reversing the [[Balanced_Truncation|Balanced Truncation]] procedure and was developed in <ref name="smith03">S.C. Smith, J. Fisher, "<span class="plainlinks">[http://ieeexplore.ieee.org/xpl/freeabs_all.jsp?arnumber=1243494 On generating random systems: a gramian approach]</span>", Proceedings of the American Control Conference, 2003.</ref>, where a description of the algorithm is given.
+
It is based on reversing the [[Balanced_Truncation|Balanced Truncation]] procedure and was developed in <ref name="smith03"/>,
  +
where a description of the algorithm is given.
 
In aggregate form, for randomly generated controllability and observability gramians, a balancing transformation is computed.
 
In aggregate form, for randomly generated controllability and observability gramians, a balancing transformation is computed.
The balanced gramian is the basis for an associated state-space system, which is determined by solving a [[wikipedia:Lyapunov_equation|Lyapunov equation]] and then unbalanced.
+
The balanced gramian is the basis for an associated state-space system,
  +
which is determined by solving a [[wikipedia:Lyapunov_equation|Lyapunov equation]] and then unbalanced.
 
A central point is the solution of the Lyapunov equations for the system matrix instead of the gramian matrix.
 
A central point is the solution of the Lyapunov equations for the system matrix instead of the gramian matrix.
 
This is feasable due to the symmetric (semi-)positive definiteness of the gramians and the requirement for a stable system, yet with a non-unique solution.
 
This is feasable due to the symmetric (semi-)positive definiteness of the gramians and the requirement for a stable system, yet with a non-unique solution.
  +
Following, the steps of the ILP are listed:
   
  +
# Sample eigenvalues of controllability and observability Gramians.
==Implementation==
 
  +
# Generate random orthogonal matrices (ie.: SVD of random matrix).
  +
# Compute balancing transformation for these Gramians.
  +
# Sample random input and output matrices.
  +
# Scale output matrix to input matrix.
  +
# Solve Lyapunov equation for system matrix.
  +
# Unbalance system.
   
  +
===Inverse Sylvester Procedure===
An efficient approach to solving the Lyapunov equation is provided by empirical gramians.
 
  +
A variant of the '''Inverse Lyapunov Procedure''' is the inverse Sylvester procedure (ILS) <ref name="himpe17"/>,
  +
which generates only state-space symmetric systems.
  +
Instead of balanced truncation, the [[wikipedia:Cross_Gramian|cross Gramian]] is utilized for the random system generation, and hence a [[wikipedia:Sylvester_equation|Sylvester equation]] is needs to be solved.
  +
The steps for the ILS are listed below:
   
  +
# Sample cross Gramian eigenvalues.
==Usage==
 
  +
# Sample random input matrix, and set output matrix as its transpose.
  +
# Solve Sylvester equation for system matrix.
  +
# Sample orthogonal unbalancing transformation (QR of random matrix).
  +
# Unbalance system.
   
  +
Even though the ILS is more limited than the ILP, for large systems it can be more efficient.
Use the following matlab code to generate a random system as described above:
 
   
  +
==Data==
  +
This benchmark is procedural and the input, state and output dimensions can be chosen.
  +
Use the following [http://matlab.com MATLAB] code to generate a random system as described above:
  +
  +
<div class="thumbinner" style="width:20%;text-align:left;"><!--[[Media:ilp.m|-->
 
<source lang="matlab">
 
<source lang="matlab">
   
function [A B C] = ilp(J,N,O,s,r)
+
function [A,B,C] = ilp(M,N,Q,s,r)
 
% ilp (inverse lyapunov procedure)
 
% ilp (inverse lyapunov procedure)
% by Christian Himpe, 2013 ( http://gramian.de )
+
% by Christian Himpe, 2013--2018
% released under BSD 2-Clause License ( http://gramian.de/#license )
+
% released under BSD 2-Clause License
 
%*
 
%*
  +
if(nargin==5)
  +
rand('seed',r);
  +
randn('seed',r);
  +
end;
   
  +
% Gramian Eigenvalues
if(exist('emgr')~=2) disp('emgr framework is required. Download at http://gramian.de/emgr.m'); return; end
 
  +
WC = exp(0.5*rand(N,1));
  +
WO = exp(0.5*rand(N,1));
   
  +
% Gramian Eigenvectors
if(nargin==5) rand('seed',r); randn('seed',r); end;
 
  +
[P,S,R] = svd(randn(N));
   
  +
% Balancing Transformation
%% Gramian Eigenvalues
 
WC = exp(-N + N*rand(N,1));
+
WC = P*diag(sqrt(WC))*P';
WO = exp(-N + N*rand(N,1));
+
WO = R*diag(sqrt(WO))*R';
  +
[U,D,V] = svd(WC*WO);
   
  +
% Input and Output
%% Gramian Eigenvectors
 
X = randn(N,N);
+
B = randn(N,M);
[U E V] = svd(X);
 
   
  +
if(nargin>=4 && s~=0)
%% Balancing Trafo
 
  +
C = B';
[P D Q] = svd(diag(WC.*WO));
 
W = -D;
+
else
  +
C = randn(Q,N);
  +
end
   
%% Input and Output
+
% Scale Output Matrix
  +
BB = sum(B.*B,2); % = diag(B*B')
B = randn(N,J);
 
  +
CC = sum(C.*C,1)'; % = diag(C'*C)
  +
C = bsxfun(@times,C,sqrt(BB./CC)');
   
  +
% Solve System Matrix
if(nargin<4 || s==0)
 
C = randn(O,N);
+
A = -sylvester(D,D,B*B');
else
 
C = B';
 
end
 
   
  +
% Unbalance System
%% Scale Output Matrix
 
  +
A = V*A*U';
BB = sum(B.*B,2); % = diag(B*B')
 
  +
B = V*B;
CC = sum(C.*C,1)'; % = diag(C'*C)
 
  +
C = C*U';
C = bsxfun(@times,C,sqrt(BB./CC)');
 
 
%% Solve System Matrix
 
f = @(x,u,p) W*x+B*u;
 
g = @(x,u,p) C*x;
 
A = -emgr(f,g,[J N O],0,[0 0.01 1],'c');
 
 
%% Unbalance System
 
T = U'*P';
 
A = T*A*T';
 
B = T*B;
 
C = C*T';
 
   
  +
end
 
</source>
 
</source>
  +
<!--]]--></div>
   
The function call requires three parameters; the number of inputs <math>J</math>, of states <math>N</math> and outputs <math>O</math>.
+
The function call requires three parameters; the number of inputs <math>M</math>, of states <math>N</math> and outputs <math>Q</math>.
Optionally, a symmetric system can be enforced with the parameter <math>s=1</math>.
+
Optionally, a symmetric system can be enforced with the parameter <math>s \neq 0</math>.
 
For reproducibility, the random number generator seed can be controlled by the parameter <math>r \in \mathbb{N}</math>.
 
For reproducibility, the random number generator seed can be controlled by the parameter <math>r \in \mathbb{N}</math>.
 
The return value consists of three matrices; the system matrix <math>A</math>, the input matrix <math>B</math> and the output matrix <math>C</math>.
 
The return value consists of three matrices; the system matrix <math>A</math>, the input matrix <math>B</math> and the output matrix <math>C</math>.
   
 
:<source lang="matlab">
 
:<source lang="matlab">
[A,B,C] = ilp(J,N,O,s,r);
+
[A,B,C] = ilp(M,N,Q,s,r);
 
</source>
 
</source>
   
  +
A variant of the above code using empirical Gramians instead of a matrix equation solution can be found at http://gramian.de/utils/ilp.m , which may yield preferable results.
'''ilp''' is compatible with [[wikipedia:MATLAB|MATLAB]] and [[wikipedia:GNU_Octave|OCTAVE]].
 
  +
The matlab code can be downloaded: [[Media:ilp.m.tar.gz|ilp.m]].
 
  +
==Dimensions==
The required [[Emgr|Empirical Gramian Framework]] can be obtained from [http://gramian.de/emgr.m http://gramian.de].
 
  +
  +
System structure:
  +
:<math>
  +
\begin{align}
  +
\dot{x}(t) &= Ax(t) + Bu(t) \\
  +
y(t) &= Cx(t)
  +
\end{align}
  +
</math>
  +
  +
System dimensions:
  +
  +
<math>A \in \mathbb{R}^{N \times N}</math>,
  +
<math>B \in \mathbb{R}^{N \times M}</math>,
  +
<math>C \in \mathbb{R}^{Q \times N}</math>.
  +
  +
==Citation==
  +
To cite this benchmark, use the following references:
  +
  +
* For the benchmark itself and its data:
  +
::The MORwiki Community, '''Inverse Lyapunov Procedure'''. MORwiki - Model Order Reduction Wiki, 2018. http://modelreduction.org/index.php/Inverse_Lyapunov_Procedure
  +
  +
@MISC{morwiki-invlyapproc,
  +
author = <nowiki>{{The MORwiki Community}}</nowiki>,
  +
title = {Inverse Lyapunov Procedure},
  +
howpublished = {{MORwiki} -- Model Order Reduction Wiki},
  +
url = <nowiki>{https://modelreduction.org/morwiki/Inverse_Lyapunov_Procedure}</nowiki>,
  +
year = 2018
  +
}
  +
  +
* For the background on the benchmark:
  +
  +
@INPROCEEDINGS{SmiF03,
  +
author = {Smith, S.~C. and Fisher, J.},
  +
title = {On generating random systems: a gramian approach},
  +
booktitle = {Proc. Am. Control. Conf.},
  +
volume = 3,
  +
pages = {2743--2748},
  +
year = 2003,
  +
doi = {10.1109/ACC.2003.1243494}
  +
}
   
 
==References==
 
==References==
<references />
+
<references>
  +
  +
<ref name="smith03">S.C. Smith, J. Fisher, "<span class="plainlinks">[https://doi.org/10.1109/ACC.2003.1243494 On generating random systems: a gramian approach]</span>", Proceedings of the American Control Conference, 3: 2743--2748, 2003.</ref>
   
  +
<ref name="himpe17">C. Himpe, M. Ohlberger, "<span class="plainlinks">[https://doi.org/10.1007/978-3-319-58786-8_17 Cross-Gramian-Based Model Reduction: A Comparison]</span>", In: Model Reduction of Parametrized Systems, Modeling, Simulation and Applications, vol. 17: 271--283, 2017.</ref>
   
  +
</references>
   
 
==Contact==
 
==Contact==

Latest revision as of 07:32, 17 June 2025


Description

The Inverse Lyapunov Procedure (ILP) is a synthetic random linear system generator. It is based on reversing the Balanced Truncation procedure and was developed in [1], where a description of the algorithm is given. In aggregate form, for randomly generated controllability and observability gramians, a balancing transformation is computed. The balanced gramian is the basis for an associated state-space system, which is determined by solving a Lyapunov equation and then unbalanced. A central point is the solution of the Lyapunov equations for the system matrix instead of the gramian matrix. This is feasable due to the symmetric (semi-)positive definiteness of the gramians and the requirement for a stable system, yet with a non-unique solution. Following, the steps of the ILP are listed:

  1. Sample eigenvalues of controllability and observability Gramians.
  2. Generate random orthogonal matrices (ie.: SVD of random matrix).
  3. Compute balancing transformation for these Gramians.
  4. Sample random input and output matrices.
  5. Scale output matrix to input matrix.
  6. Solve Lyapunov equation for system matrix.
  7. Unbalance system.

Inverse Sylvester Procedure

A variant of the Inverse Lyapunov Procedure is the inverse Sylvester procedure (ILS) [2], which generates only state-space symmetric systems. Instead of balanced truncation, the cross Gramian is utilized for the random system generation, and hence a Sylvester equation is needs to be solved. The steps for the ILS are listed below:

  1. Sample cross Gramian eigenvalues.
  2. Sample random input matrix, and set output matrix as its transpose.
  3. Solve Sylvester equation for system matrix.
  4. Sample orthogonal unbalancing transformation (QR of random matrix).
  5. Unbalance system.

Even though the ILS is more limited than the ILP, for large systems it can be more efficient.

Data

This benchmark is procedural and the input, state and output dimensions can be chosen. Use the following MATLAB code to generate a random system as described above:

function [A,B,C] = ilp(M,N,Q,s,r)
% ilp (inverse lyapunov procedure)
% by Christian Himpe, 2013--2018
% released under BSD 2-Clause License
%*
  if(nargin==5)
    rand('seed',r);
    randn('seed',r);
  end;

% Gramian Eigenvalues
  WC = exp(0.5*rand(N,1));
  WO = exp(0.5*rand(N,1));

% Gramian Eigenvectors
  [P,S,R] = svd(randn(N));

% Balancing Transformation
  WC = P*diag(sqrt(WC))*P';
  WO = R*diag(sqrt(WO))*R';
  [U,D,V] = svd(WC*WO);

% Input and Output
  B = randn(N,M);

  if(nargin>=4 && s~=0)
    C = B';
  else
    C = randn(Q,N);
  end

% Scale Output Matrix
  BB = sum(B.*B,2);  % = diag(B*B')
  CC = sum(C.*C,1)'; % = diag(C'*C)
  C = bsxfun(@times,C,sqrt(BB./CC)');

% Solve System Matrix
  A = -sylvester(D,D,B*B');

% Unbalance System
  A = V*A*U';
  B = V*B;
  C = C*U';

end

The function call requires three parameters; the number of inputs M, of states N and outputs Q. Optionally, a symmetric system can be enforced with the parameter s \neq 0. For reproducibility, the random number generator seed can be controlled by the parameter r \in \mathbb{N}. The return value consists of three matrices; the system matrix A, the input matrix B and the output matrix C.

[A,B,C] = ilp(M,N,Q,s,r);

A variant of the above code using empirical Gramians instead of a matrix equation solution can be found at http://gramian.de/utils/ilp.m , which may yield preferable results.

Dimensions

System structure:


\begin{align}
\dot{x}(t) &= Ax(t) + Bu(t) \\
y(t) &= Cx(t)
\end{align}

System dimensions:

A \in \mathbb{R}^{N \times N}, B \in \mathbb{R}^{N \times M}, C \in \mathbb{R}^{Q \times N}.

Citation

To cite this benchmark, use the following references:

  • For the benchmark itself and its data:
The MORwiki Community, Inverse Lyapunov Procedure. MORwiki - Model Order Reduction Wiki, 2018. http://modelreduction.org/index.php/Inverse_Lyapunov_Procedure
@MISC{morwiki-invlyapproc,
  author =       {{The MORwiki Community}},
  title =        {Inverse Lyapunov Procedure},
  howpublished = {{MORwiki} -- Model Order Reduction Wiki},
  url =          {https://modelreduction.org/morwiki/Inverse_Lyapunov_Procedure},
  year =         2018
}
  • For the background on the benchmark:
@INPROCEEDINGS{SmiF03,
  author =       {Smith, S.~C. and Fisher, J.},
  title =        {On generating random systems: a gramian approach},
  booktitle =    {Proc. Am. Control. Conf.},
  volume =       3,
  pages =        {2743--2748},
  year =         2003,
  doi =          {10.1109/ACC.2003.1243494}
}

References

  1. S.C. Smith, J. Fisher, "On generating random systems: a gramian approach", Proceedings of the American Control Conference, 3: 2743--2748, 2003.
  2. C. Himpe, M. Ohlberger, "Cross-Gramian-Based Model Reduction: A Comparison", In: Model Reduction of Parametrized Systems, Modeling, Simulation and Applications, vol. 17: 271--283, 2017.

Contact

Christian Himpe