matrix_eqn¶
Filename: matrix_eqn.py
This files holds several functions that are used to solve matrix equations. Currently has functionality to solve:
- Lyapunov Equations
- Ricatti Equations
- TODO: 1. See issue 47 on github repository, should add support for
- Sylvester equations 2. Fix warnings from checking conditioning of matrices
-
quantecon.matrix_eqn.solve_discrete_lyapunov(A, B, max_it=50, method='doubling')[source]¶ Computes the solution to the discrete lyapunov equation
\[AXA' - X + B = 0\]X is computed by using a doubling algorithm. In particular, we iterate to convergence on X_j with the following recursions for j = 1, 2,... starting from X_0 = B, a_0 = A:
\[a_j = a_{j-1} a_{j-1}\]\[X_j = X_{j-1} + a_{j-1} X_{j-1} a_{j-1}'\]Parameters: A : array_like(float, ndim=2)
An n x n matrix as described above. We assume in order for convergence that the eigenvalues of A have moduli bounded by unity
B : array_like(float, ndim=2)
An n x n matrix as described above. We assume in order for convergence that the eigenvalues of A have moduli bounded by unity
max_it : scalar(int), optional(default=50)
The maximum number of iterations
method : string, optional(default=”doubling”)
Describes the solution method to use. If it is “doubling” then uses the doubling algorithm to solve, if it is “bartels-stewart” then it uses scipy’s implementation of the Bartels-Stewart approach.
Returns: gamma1: array_like(float, ndim=2)
Represents the value V
-
quantecon.matrix_eqn.solve_discrete_riccati(A, B, Q, R, N=None, tolerance=1e-10, max_iter=500)[source]¶ Solves the discrete-time algebraic Riccati equation
X = A’XA - (N + B’XA)’(B’XB + R)^{-1}(N + B’XA) + Qvia a modified structured doubling algorithm. An explanation of the algorithm can be found in the reference below.
Note that SciPy also has a discrete riccati equation solver. However it cannot handle the case where R is not invertible, or when N is nonzero. Both of these cases can be handled in the algorithm implemented below.
Parameters: A : array_like(float, ndim=2)
k x k array.
B : array_like(float, ndim=2)
k x n array
Q : array_like(float, ndim=2)
k x k, should be symmetric and non-negative definite
R : array_like(float, ndim=2)
n x n, should be symmetric and positive definite
N : array_like(float, ndim=2)
n x k array
tolerance : scalar(float), optional(default=1e-10)
The tolerance level for convergence
max_iter : scalar(int), optional(default=500)
The maximum number of iterations allowed
Returns: X : array_like(float, ndim=2)
The fixed point of the Riccati equation; a k x k array representing the approximate solution
References
Chiang, Chun-Yueh, Hung-Yuan Fan, and Wen-Wei Lin. “STRUCTURED DOUBLING ALGORITHM FOR DISCRETE-TIME ALGEBRAIC RICCATI EQUATIONS WITH SINGULAR CONTROL WEIGHTING MATRICES.” Taiwanese Journal of Mathematics 14, no. 3A (2010): pp-935.