矩阵方程的计算求解（Matlab实现）

一、Lyapunov方程的计算求解

1、连续Lyapunov方程
连续Lyapunov方程可以表示为：

AX + XA* = -C  % 其中A*是A的转置

Lyapunov方程源于微分方程稳定性理论，其中要求-C为对称正定的nxn矩阵，从而可以证明解X亦为nxn对称矩阵。Lyapunov类的方程求解是很困难的，可以利用Matlab控制系统工具箱中提供的lyap函数求解，调用格式为

X = lyap(A, C)

matlab代码：

A = [1 2 3;4 5 6;7 8 0];  C = -[10 5 4;5 6 7;4 7 9];
X = lyap(A, C)
norm(A*X + X*A' + C)    % 验证解的情况

% 结果：
>> Matrix_equation

X =

   -3.9444    3.8889    0.3889
    3.8889   -2.7778    0.2222
    0.3889    0.2222   -0.1111


ans =

   2.3211e-14

---

2、Lyapunov方程的解析解

---

3、离散 Lyapunov方程
离散Lyapunov方程可以表示为：

AXA* - X + Q = 0   % 其中A*是A的转置矩阵

该方程可以由MATLAB控制系统工具箱的dlyap函数直接求解。该函数的调用格式为：

X = dlyap(A, Q)

matlab代码为：

A = [1 2 3;4 5 6;7 8 0];  Q = -[10 5 4;5 6 7;4 7 9];
X = dlyap(A, Q)
norm(A*X*A'- X + Q)   % 精度验证

% 结果：

X =

   -2.8439    3.2500   -3.0501
    3.2500   -3.3780    2.8107
   -3.0501    2.8107   -0.5462


ans =

   7.6172e-14

二、Sylvester方程的计算求解

Sylvester方程的一般形式为：

AX + XB = -C 

其中，A为nxn矩阵，B为mxm矩阵，C和X均为nxm矩阵。该方程又称为广义的Lyapunov方程。仍然可以用Matlab中控制系统工具箱中的lyap函数直接求解该方程。函数的一般调用格式为：

X = lyap(A,B,C)

该函数采用的是Schur分解的数值解法求解方程。
matlab代码：

A = [1 2 3;2 4 1;4 6 1];
B = [2 3 5;2 7 5;5 4 3];
C = -[3 4 2;3 2 4;2 0 9];
X = lyap(A, B, C)
norm(A*X + X*B + C )

% 结果：

X =

   -9.5651   10.3207   -4.3218
    1.4515   -1.7102    1.3843
    9.9199   -9.7210    4.2467


ans =

   3.9005e-14

三、Riccati方程的计算求解

Riccati方程是一类很著名的二次型矩阵方程式，其一般形式为：

A*X + XA - XBX + C = 0    % A*是A的转置矩阵

由于含有未知矩阵X的二次项， 所以Riccati方程的求解数学上要比Lyapunov方程更难。Matlab的控制系统工具箱提供了现成函数are，调用形式如下：

X = are(A, B, C)

matlab代码：

A = [1 2 3;2 4 1;4 6 1];
B = [-2 3 5;2 7 5;5 4 3];
C = [3 4 2;3 2 4;2 0 9];
X = are(A, B, C)
norm(A'*X + X*A - X*B*X + C) 

% 结果：

X =

   -0.1180    1.4662   -0.6059
    0.4316    1.4014    0.0150
    0.9982   -0.4684    2.0600


ans =

   3.2899e-14