Lyapunov equation

In the following theorems $A,\; P,\; Q\; \backslash in\; \backslash mathbb\{R\}^\{n\; \backslash times\; n\}$, and $P$ and $Q$ are symmetric. The notation $P>0$ means that the matrix $P$ is positive definite.

Theorem (continuous time version). Given any $Q>0$, there exists a unique $P>0$ satisfying $A^T\; P\; +\; P\; A\; +\; Q\; =\; 0$ if and only if the linear system $\backslash dot\{x\}=A\; x$ is globally asymptotically stable. The quadratic function $V(x)=x^T\; P\; x$ is a Lyapunov function that can be used to verify stability.

Theorem (discrete time version). Given any $Q>0$, there exists a unique $P>0$ satisfying $A^T\; P\; A\; -P\; +\; Q\; =\; 0$ if and only if the linear system $x\_\{t+1\}=A\; x\_\{t\}$ is globally asymptotically stable. As before, $x^T\; P\; x$ is a Lyapunov function.

The Lyapunov equation is linear; therefore, if $X$ contains $n$ entries, the equation can be solved in $\backslash mathcal\; O(n^3)$ time using standard matrix factorization methods.

However, specialized algorithms are available which can yield solutions much quicker owing to the specific structure of the Lyapunov equation. For the discrete case, the Schur method of Kitagawa is often used.{{cite journal |last=Kitagawa |first=G. |title=An Algorithm for Solving the Matrix Equation X = F X F' + S |journal=International Journal of Control |volume=25 |issue=5 |pages=745–753 |year=1977 |doi=10.1080/00207177708922266 }} For the continuous Lyapunov equation the Bartels–Stewart algorithm can be used.{{cite journal |first=R. H. |last=Bartels |first2=G. W. |last2=Stewart |title=Algorithm 432: Solution of the matrix equation AX + XB = C |journal=Comm. ACM |volume=15 |year=1972 |issue=9 |pages=820–826 |doi=10.1145/361573.361582 |doi-access=free }}

Defining the vectorization operator $\backslash operatorname\{vec\}\; (A)$ as stacking the columns of a matrix $A$ and $A\; \backslash otimes\; B$ as the Kronecker product of $A$ and $B$, the continuous time and discrete time Lyapunov equations can be expressed as solutions of a matrix equation. Furthermore, if the matrix $A$ is "stable", the solution can also be expressed as an integral (continuous time case) or as an infinite sum (discrete time case).

Using the result that $\backslash operatorname\{vec\}(ABC)=(C^\{T\}\; \backslash otimes\; A)\backslash operatorname\{vec\}(B)$, one has

:$(I\_\{n^2\}-\backslash bar\{A\}\; \backslash otimes\; A)\backslash operatorname\{vec\}(X)\; =\; \backslash operatorname\{vec\}(Q)$

where $I\_\{n^2\}$ is a conformable identity matrix and $\backslash bar\{A\}$ is the element-wise complex conjugate of $A$.{{cite book |first=J. |last=Hamilton |year=1994 |title=Time Series Analysis |at=Equations 10.2.13 and 10.2.18 |publisher=Princeton University Press |isbn=0-691-04289-6 }} One may then solve for $\backslash operatorname\{vec\}(X)$ by inverting or solving the linear equations. To get $X$, one must just reshape $\backslash operatorname\{vec\}\; (X)$ appropriately.

Moreover, if $A$ is stable (in the sense of Schur stability, i.e., having eigenvalues with magnitude less than 1), the solution $X$ can also be written as

:$X\; =\; \backslash sum\_\{k=0\}^\{\backslash infty\}\; A^\{k\}\; Q\; (A^\{H\})^k$.

For comparison, consider the one-dimensional case, where this just says that the solution of $(1\; -\; a^2)\; x\; =\; q$ is

:$x\; =\; \backslash frac\{q\}\{1-a^2\}\; =\; \backslash sum\_\{k=0\}^\{\backslash infty\}\; qa^\{2k\}$.

Using again the Kronecker product notation and the vectorization operator, one has the matrix equation

:$(I\_n\; \backslash otimes\; A\; +\; \backslash bar\{A\}\; \backslash otimes\; I\_n)\; \backslash operatorname\{vec\}X\; =\; -\backslash operatorname\{vec\}Q,$

where $\backslash bar\{A\}$ denotes the matrix obtained by complex conjugating the entries of $A$.

Similar to the discrete-time case, if $A$ is stable (in the sense of Hurwitz stability, i.e., having eigenvalues with negative real parts), the solution $X$ can also be written as

:$X\; =\; \backslash int\_0^\backslash infty\; \{e\}^\{A\; \backslash tau\}\; Q\; \backslash mathrm\{e\}^\{A^H\; \backslash tau\}\; d\backslash tau$,

which holds because

:

=&\int_0^\infty \frac{d}{d\tau} {e}^{A \tau} Q \mathrm{e}^{A^H \tau} d\tau \\

=& {e}^{A \tau} Q \mathrm{e}^{A^H \tau} \bigg|_0^\infty\\

=& -Q.

\end{align}

For comparison, consider the one-dimensional case, where this just says that the solution of $2ax\; =\; -\; q$ is

:$x\; =\; \backslash frac\{-q\}\{2a\}\; =\; \backslash int\_0^\{\backslash infty\}\; q\{e\}^\{2\; a\; \backslash tau\}\; d\backslash tau$.

We start with the continuous-time linear dynamics:

:$\backslash dot\{\backslash mathbf\{x\}\}\; =\; \backslash mathbf\{A\}\backslash mathbf\{x\}$.

And then discretize it as follows:

:$\backslash dot\{\backslash mathbf\{x\}\}\; \backslash approx\; \backslash frac\{\backslash mathbf\{x\}\_\{t+1\}-\backslash mathbf\{x\}\_\{t\}\}\{\backslash delta\}$

Where $\backslash delta\; >\; 0$ indicates a small forward displacement in time. Substituting the bottom equation into the top and shuffling terms around, we get a discrete-time equation for $\backslash mathbf\{x\}\_\{t+1\}$.

$\backslash mathbf\{x\}\_\{t+1\}\; =\; \backslash mathbf\{x\}\_t\; +\; \backslash delta\; \backslash mathbf\{A\}\; \backslash mathbf\{x\}\_t\; =\; (\backslash mathbf\{I\}\; +\; \backslash delta\backslash mathbf\{A\})\backslash mathbf\{x\}\_t\; =\; \backslash mathbf\{B\}\backslash mathbf\{x\}\_t$

Where we've defined $\backslash mathbf\{B\}\; \backslash equiv\; \backslash mathbf\{I\}\; +\; \backslash delta\backslash mathbf\{A\}$. Now we can use the discrete time Lyapunov equation for $\backslash mathbf\{B\}$:

$\backslash mathbf\{B\}^T\backslash mathbf\{M\}\backslash mathbf\{B\}\; -\; \backslash mathbf\{M\}\; =\; -\backslash delta\backslash mathbf\{Q\}$

Plugging in our definition for $\backslash mathbf\{B\}$, we get:

$(\backslash mathbf\{I\}\; +\; \backslash delta\; \backslash mathbf\{A\})^T\backslash mathbf\{M\}(\backslash mathbf\{I\}\; +\; \backslash delta\; \backslash mathbf\{A\})\; -\; \backslash mathbf\{M\}\; =\; -\backslash delta\; \backslash mathbf\{Q\}$

Expanding this expression out yields:

$(\backslash mathbf\{M\}\; +\; \backslash delta\; \backslash mathbf\{A\}^T\backslash mathbf\{M\})\; (\backslash mathbf\{I\}\; +\; \backslash delta\; \backslash mathbf\{A\})\; -\; \backslash mathbf\{M\}\; =\; \backslash delta(\backslash mathbf\{A\}^T\backslash mathbf\{M\}\; +\; \backslash mathbf\{M\}\backslash mathbf\{A\})\; +\; \backslash delta^2\; \backslash mathbf\{A\}^T\backslash mathbf\{M\}\backslash mathbf\{A\}\; =\; -\backslash delta\; \backslash mathbf\{Q\}$

Recall that $\backslash delta$ is a small displacement in time. Letting $\backslash delta$ go to zero brings us closer and closer to having continuous dynamics—and in the limit we achieve them. It stands to reason that we should also recover the continuous-time Lyapunov equations in the limit as well. Dividing through by $\backslash delta$ on both sides, and then letting $\backslash delta\; \backslash to\; 0$ we find that:

$\backslash mathbf\{A\}^T\backslash mathbf\{M\}\; +\; \backslash mathbf\{M\}\backslash mathbf\{A\}\; =\; -\backslash mathbf\{Q\}$

which is the continuous-time Lyapunov equation, as desired.

• Sylvester equation, which generalizes the Lyapunov equation
• Algebraic Riccati equation
• Kalman filter

{{Reflist}}

{{Authority control}}

{{DEFAULTSORT:Lyapunov Equation}}

Category:Control theory