# Solving lyapunov equation, Matlab has different solution, why?

I need to solve the lyapunov equation i.e. $A^TP + PA = -Q$. With $A = \begin{bmatrix} -2 & 1 \\ -1 & 0 \end{bmatrix}$ and $Q = \begin{bmatrix}1 & 0 \\ 0 & 1\end{bmatrix}$.

Hence...

$$\begin{bmatrix} -2 & -1 \\ 1 & 0 \end{bmatrix} \begin{bmatrix} P_{11} & P_{12} \\ P_{12} & P_{22} \end{bmatrix} + \begin{bmatrix} P_{11} & P_{12} \\ P_{12} & P_{22} \end{bmatrix} \begin{bmatrix} -2 & 1 \\ -1 & 0 \end{bmatrix} + \begin{bmatrix} 1 & 0 \\ 0 & 1 \end{bmatrix} = 0$$

$$\begin{bmatrix} -4P_{11} - 2P_{12} + 1 & P_{11} -2P_{12} - P_{22} \\ P_{11} - 2P_{12} - P_{22} & 2P_{12} + 1 \end{bmatrix} = 0$$

$$\begin{bmatrix} -4 & -2 & 0 \\ 1 & -2 & -1 \\ 1 & -2 & -1 \\ 0 & 2 & 0 \end{bmatrix} \begin{bmatrix} P_{11} \\ P_{12} \\ P_{22} \end{bmatrix} = \begin{bmatrix} -1 \\ 0 \\ 0 \\ -1 \end{bmatrix} \Rightarrow \begin{bmatrix} -4 & -2 & 0 \\ 1 & -2 & -1 \\ 0 & 2 & 0 \end{bmatrix} \begin{bmatrix} P_{11} \\ P_{12} \\ P_{22} \end{bmatrix} = \begin{bmatrix} -1 \\ 0 \\ -1 \end{bmatrix}$$

And such we get that $P = \begin{bmatrix} 1/2 & -1/2 \\ -1/2 & 3/2\end{bmatrix}$.

This is the same solution as given by my professor. I wanted to check however if I can also find the solution using Matlab. I entered the following:

A = [-2 1; -1 0];
Q = [1 0; 0 1];
P = lyap(A,Q)


This however tells me that $P = \begin{bmatrix} 1/2 & 1/2 \\ 1/2 & 3/2\end{bmatrix}$.

What is going on here? Is Matlab correct or wrong? Or is my solution wrong? Or are we both correct?

-
Taking a look at [the documentation][1], you can see, that matlab actually solves $$AX+XA^T+Q=0$$ [1]: mathworks.de/de/help/control/ref/lyap.html – k1next Jan 24 '13 at 14:50

The matlab definition must be different from yours. With Matlab's $P$, one has $AP+PA^T=-Q$. Doing things your way with matlab's $P$ gives $A^TP+PA=[[-3,-2],[-2,1]]$. (sorry, to lazy to set up as usual matrix format).

Maybe look in matlab's documentation to see what equation their lyap(A,Q) is solving; I'd guess it had the transpose switched.

EDIT: It seems Mathematica9 has the other $AP+PA^T=-Q$ definition, or close to it. Their command LyapunovSolve[a,c] gives solution to $ax+xa^T=c$, which has the transposed matrix mentioned last, opposite to your version. Reference page:

http://reference.wolfram.com/mathematica/ref/LyapunovSolve.html

-

X = lyap(A,B,-C) solves the continuous-time Sylvester equation AX + XB = C

and

X = lyap(A’,Q) solves the continuous-time Lyapunov equation ATP + PA + Q = 0

so, you can solve the lyapunov function.

A = [-2 1; -1 0]; Q = [1 0; 0 1]; P = lyap(A',Q)

P = [0.5000 -0.5000] [-0.5000 1.5000]

-