Control theory

The following algorithms draw heavily from the second chapter of Nicholas J. Higham’s “Functions of Matrices: Theory and Computation” and depend heavily on the matrix sign function. They have only undergone cursory testing.

Sylvester

As long as both \(A\) and \(B\) only have eigenvalues in the open right-half plane, the following routines solve for \(X\) in the Sylvester equation

\[A X + X B = C\]

via computing \(\text{sgn}(W)\), where

\[\begin{split}W = \begin{pmatrix} A & -C \\ 0 & -B \end{pmatrix}.\end{split}\]

int Sylvester(const Matrix<F> &A, const Matrix<F> &B, const Matrix<F> &C, Matrix<F> &X)

int Sylvester(const DistMatrix<F> &A, const DistMatrix<F> &B, const DistMatrix<F> &C, DistMatrix<F> &X)

One may also directly pass in \(W\) in order to save space.

int Sylvester(int m, Matrix<F> &W, Matrix<F> &X)

int Sylvester(int m, DistMatrix<F> &W, DistMatrix<F> &X)

Lyapunov

A special case of the Sylvester solver, where \(B = A^H\).

int Lyapunov(const Matrix<F> &A, const Matrix<F> &C, Matrix<F> &X)

int Lyapunov(const DistMatrix<F> &A, const DistMatrix<F> &C, DistMatrix<F> &X)

Algebraic Ricatti

Under suitable conditions, the following routines solve for \(X\) in the algebraic Ricatti equation

\[X K X - A^H X - X A = L,\]

where both \(K\) and \(L\) are Hermitian. In each case, the number of Newton iterations required for convergence of \(\text{sgn}(W)\), where

\[\begin{split}W = \begin{pmatrix} A^H & L \\ K & -A \end{pmatrix},\end{split}\]

is returned.

int Ricatti(UpperOrLower uplo, const Matrix<F> &A, const Matrix<F> &K, const Matrix<F> &L, Matrix<F> &X)

int Ricatti(UpperOrLower uplo, const DistMatrix<F> &A, const DistMatrix<F> &K, const DistMatrix<F> &L, DistMatrix<F> &X)

Alternatively, one may directly fill the matrix \(W\).

int Ricatti(Matrix<F> &W, Matrix<F> &X)

int Ricatti(DistMatrix<F> &W, DistMatrix<F> &X)