Inverse Lyapunov Procedure

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 , 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.

Inverse Sylvester Procedure

A variant of the Inverse Lyapunov Procedure is the inverse Sylvester procedure .

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 = {http://modelreduction.org/index.php/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}
    }

Contact

Christian Himpe