Stability analysis of a time dependent Lyapunov equation

17 views (last 30 days)
Anil on 29 Feb 2024
Commented: Anil on 2 Mar 2024
I want solution of coupled linear equations for which I am using time dependnet Lyapunov equation:
dP/dt=A*P + P*A'+ D,
where A (A' is transpose) is a system matrix, P is covariance matrix to be found and D is diffusion matrix. I solve it staright-forward as coupled equations extracted from a matrix in terms of a col. vector but the solution blows up i.e. "instabilities", even though the real eignevalues of matrix A are negative. Whie seraching literatures for this, an options is given, using vectorization: A*P+P*A'=(AI+IA)*P,. using this, giving me an error of incompatibel sizes. I am not sure if this vectorization is enough to make sure the stability of the system but I dont know other options, please suggets if any. For that vectorization, I am using a function:
% showing just a stucture of the function file
function dpdt=(t,p0)
dpdt=zero(m,1)% m is number of rows
[m,n]=size(A);
P1=reshape(P0, m,1);
D1=reshape(D,m,1);% D is some matrix
b=eye(m);
a=real(eig(A))
if all(a<0)
dprdt=(tensorprod(A.',b)+tensorprod(b,A.')).*P1 + D1;% error appears here;
dpdt(1)=dprdt(1);
.
.
.
.
dpdt(m)=dprdt(m)
end
end
is the above right way of doing the vectorization? Any other suggestions for the stability of this time depedent Lyapunov equation, I shall be grateful. Thanks.

Christine Tobler on 29 Feb 2024
The tensorprod function doesn't do what you are looking for here, that would be kron (short for Kronecker product).
It also seems like P being a column vector of size mx1 would have the wrong size, are you looking for it to be a m x m matrix that you reshape to be m*m x 1?
In terms of performance, you should expect A*P + P*A' (and then possibly reshaping the result of this) to be much more efficient than forming the (A' kron I + I kron A') matrix explicitly.
3 CommentsShow 1 older commentHide 1 older comment
Christine Tobler on 1 Mar 2024
The equation A*B + P*A.' implies that P is an m x m matrix. That means that if it is reshaped, it must still have m^2 elements, so mx1 would not be a possible result of that reshaping.
You might start out to check your code is doing the right thing by using diagonal matrices A - then eigenvalues being negative should be enough. For non-symmetric A, it's possible for the values to grow before starting to stabilize.
Anil on 2 Mar 2024
Dear @Christine Tobler, Yes, that's true, reshaped P is m*m x1. I missunderstood that. Thanks for the later part, I am checking it.