Skip to content

Fix test sg03ad #82

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Closed
wants to merge 5 commits into from
Closed

Fix test sg03ad #82

Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
fd830ee
sg03ad: arange instead of range
Nov 8, 2019
d173984
remove unused module, clean pep8 warning
Nov 8, 2019
df3b944
fix sg03ad test for continuos time and add discrete time test
bnavigator Nov 15, 2019
caecd0b
docstrings for sg03ad tests
bnavigator Nov 15, 2019
7f393f2
Merge branch 'master' into fix_test_sg03ad_static
bnavigator Apr 3, 2020
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
63 changes: 39 additions & 24 deletions slycot/tests/test_sg03ad.py
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,59 +1,74 @@
#!/usr/bin/env python
#
# test_sg03ad.py - test suit for stability margin commands
# test_sg03ad.py - test suite for stability margin commands
# RvP, 15 Jun 2017
from __future__ import print_function

import unittest
from slycot import synthesis
import numpy as np

from numpy.testing import assert_raises, assert_almost_equal
from numpy.testing import assert_almost_equal

# test cases from
# Penzl T., Numerical Solution of Generalized Lyapunov Equations
# http://www.qucosa.de/fileadmin/data/qucosa/documents/4168/data/b002.pdf


class test_sg03ad(unittest.TestCase):

def test_sg03ad_a(self):
# Example 1
def test_sg03ad_ex1c(self):
""" Example 1 continuous case"""
n = 100
Xref = np.ones((n,n))
Xref = np.ones((n, n))
U = np.tril(Xref)
for t in range(0, 50, 10):
A = 2.0**(-t) - np.eye(n) + np.diag(range(1,n+1)) + U.T
E = np.eye(n) + 2**(-t)*U
A = (2**(-t) - 1) * np.eye(n) + np.diag(np.arange(1., n+1.)) + U.T
E = np.eye(n) + 2**(-t) * U
Y = A.T.dot(Xref).dot(E) + E.T.dot(Xref).dot(A)
Q = np.zeros((n,n))
Z = np.zeros((n,n))
Q = np.zeros((n, n))
Z = np.zeros((n, n))
A, E, Q, Z, X, scale, sep, ferr, alphar, alphai, beta = \
synthesis.sg03ad('C', 'B', 'N', 'N', 'L', n, A, E, Q, Z, Y)
assert_almost_equal(Xref, X)

def test_sg03ad_3(self):
assert_almost_equal(X, Xref)

def test_sg03ad_ex1d(self):
""" Example 1 discrete case"""
n = 100
Xref = np.ones((n, n))
U = np.tril(Xref)
for t in range(0, 50, 10):
A = 2**(-t) * np.eye(n) + np.diag(np.arange(1., n+1.)) + U.T
E = np.eye(n) + 2**(-t) * U
Y = A.T.dot(Xref).dot(A) - E.T.dot(Xref).dot(E)
Q = np.zeros((n, n))
Z = np.zeros((n, n))
A, E, Q, Z, X, scale, sep, ferr, alphar, alphai, beta = \
synthesis.sg03ad('D', 'B', 'N', 'N', 'L', n, A, E, Q, Z, Y)
assert_almost_equal(X, Xref)

def test_sg03ad_b1(self):
""" SLICOT doc example / Penzl B.1 """
n = 3
A = np.array([[3.0, 1.0, 1.0],
[1.0, 3.0, 0.0],
[1.0, 0.0, 2.0]])
E = np.array([[1.0, 3.0, 0.0],
[3.0, 2.0, 1.0],
[1.0, 0.0, 1.0]])
[3.0, 2.0, 1.0],
[1.0, 0.0, 1.0]])
Y = np.array([[64.0, 73.0, 28.0],
[73.0, 70.0, 25.0],
[28.0, 25.0, 18.0]])
[73.0, 70.0, 25.0],
[28.0, 25.0, 18.0]])
Xref = np.array([[-2.0000, -1.0000, 0.0000],
[-1.0000, -3.0000, -1.0000],
[0.0000, -1.0000, -3.0000]])
Q = np.zeros((3,3))
Z = np.zeros((3,3))
[-1.0000, -3.0000, -1.0000],
[0.0000, -1.0000, -3.0000]])
Q = np.zeros((3, 3))
Z = np.zeros((3, 3))
A, E, Q, Z, X, scale, sep, ferr, alphar, alphai, beta = \
synthesis.sg03ad('C', 'B', 'N', 'N', 'L', n, A, E, Q, Z, -Y)
#print(A, E, Q, Z, X, scale, sep)
# print(A, E, Q, Z, X, scale, sep)
assert_almost_equal(X, Xref)

def suite():
return unittest.TestLoader().loadTestsFromTestCase(TestConvert)


if __name__ == "__main__":
unittest.main()