Merge pull request #108 from ChadFulton/dlyap

FIX: Discrete Lyapunov with complex input
QuantEcon · Jan 24, 2015 · 9b60437 · 9b60437
2 parents 8046e17 + 8956e1e
commit 9b60437
Show file tree

Hide file tree

Showing 2 changed files with 12 additions and 1 deletion.
diff --git a/quantecon/matrix_eqn.py b/quantecon/matrix_eqn.py
@@ -75,7 +75,7 @@ def solve_discrete_lyapunov(A, B, max_it=50, method="doubling"):
         while diff > 1e-15:
 
             alpha1 = alpha0.dot(alpha0)
-            gamma1 = gamma0 + np.dot(alpha0.dot(gamma0), alpha0.T)
+            gamma1 = gamma0 + np.dot(alpha0.dot(gamma0), alpha0.conjugate().T)
 
             diff = np.max(np.abs(gamma1 - gamma0))
             alpha0 = alpha1

diff --git a/quantecon/tests/test_matrix_eqn.py b/quantecon/tests/test_matrix_eqn.py
@@ -31,3 +31,14 @@ def test_solve_discrete_lyapunov_B():
 
     assert_allclose(B, X)
 
+def test_solve_discrete_lyapunov_complex():
+    'Complex test, A is companion matrix'
+    A = np.array([[0.5 + 0.3j, 0.1 + 0.1j],
+                  [         1,          0]])
+    B = np.eye(2)
+
+    X = qme.solve_discrete_lyapunov(A, B)
+
+    assert_allclose(np.dot(np.dot(A, X), A.conj().transpose()) - X, -B,
+                    atol=1e-15)
+