solve_discrete_are#
- scipy.linalg.solve_discrete_are(a, b, q, r, e=None, s=None, balanced=True)[源代码]#
求解离散时间代数Riccati方程 (DARE)。
DARE 定义为
\[A^HXA - X - (A^HXB) (R + B^HXB)^{-1} (B^HXA) + Q = 0\]存在解的限制条件是
\(A\) 的所有特征值在单位圆之外,应该是可控的。
相关的辛铅笔(参见注释)的特征值应充分远离单位圆。
此外,如果
e
和s
并非都精确地为None
,则广义 DARE\[A^HXA - E^HXE - (A^HXB+S) (R+B^HXB)^{-1} (B^HXA+S^H) + Q = 0\]被求解。当省略时,
e
假定为单位矩阵,而s
假定为零矩阵。- 参数:
- a(M, M) 数组型
方阵
- b(M, N) 数组型
输入
- q(M, M) 数组型
输入
- r(N, N) 数组型
方阵
- e(M, M) 数组型,可选
非奇异方阵
- s(M, N) 数组型,可选
输入
- balancedbool
指示是否对数据执行平衡步骤的布尔值。默认设置为 True。
- 返回:
- x(M, M) ndarray
离散代数Riccati方程的解。
- 引发:
- LinAlgError
对于无法隔离铅笔的稳定子空间的情况。有关详细信息,请参见注释部分和参考资料。
另请参阅
solve_continuous_are
求解连续代数Riccati方程
注释
该方程通过形成扩展的辛矩阵铅笔来求解,如 [1] 中所述,\(H - \lambda J\) 由分块矩阵给出
[ A 0 B ] [ E 0 B ] [ -Q E^H -S ] - \lambda * [ 0 A^H 0 ] [ S^H 0 R ] [ 0 -B^H 0 ]
并使用 QZ 分解方法。
在此算法中,失败条件与乘积 \(U_2 U_1^{-1}\) 的对称性和 \(U_1\) 的条件数有关。在这里,\(U\) 是 2m×m 矩阵,它保存跨越稳定子空间的特征向量,具有 2-m 行并被分成两个 m 行矩阵。有关更多详细信息,请参见 [1] 和 [2]。
为了提高 QZ 分解精度,铅笔会经历一个平衡步骤,其中 \(H\) 和 \(J\) 行/列的绝对值之和(去除对角线项后)按照 [3] 中给出的方法进行平衡。如果数据具有小的数值噪声,平衡可能会放大其影响,因此需要进行一些清理。
在 0.11.0 版本中添加。
参考资料
[1] (1,2)P. van Dooren , “求解 Riccati 方程的广义特征值方法”, SIAM Journal on Scientific and Statistical Computing, Vol.2(2), DOI:10.1137/0902010
[2]A.J. Laub, “求解代数Riccati方程的Schur方法”, 马萨诸塞理工学院。信息与决策系统实验室。 LIDS-R ; 859。在线可用:http://hdl.handle.net/1721.1/1301
[3]P. Benner, “哈密顿矩阵的辛平衡”,2001 年,SIAM J. Sci. Comput., 2001, Vol.22(5), DOI:10.1137/S1064827500367993
示例
给定 a、b、q 和 r,求解 x
>>> import numpy as np >>> from scipy import linalg as la >>> a = np.array([[0, 1], [0, -1]]) >>> b = np.array([[1, 0], [2, 1]]) >>> q = np.array([[-4, -4], [-4, 7]]) >>> r = np.array([[9, 3], [3, 1]]) >>> x = la.solve_discrete_are(a, b, q, r) >>> x array([[-4., -4.], [-4., 7.]]) >>> R = la.solve(r + b.T.dot(x).dot(b), b.T.dot(x).dot(a)) >>> np.allclose(a.T.dot(x).dot(a) - x - a.T.dot(x).dot(b).dot(R), -q) True