Skip to content

Allow minreal to work with static gain StateSpace objects #108

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.

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
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
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
80 changes: 39 additions & 41 deletions control/statesp.py
Original file line number Diff line number Diff line change
Expand Up @@ -122,34 +122,35 @@ def __init__(self, *args):
else:
raise ValueError("Needs 1 or 4 arguments; received %i." % len(args))

# Here we're going to convert inputs to matrices, if the user gave a
# non-matrix type.
#! TODO: [A, B, C, D] = map(matrix, [A, B, C, D])?
matrices = [A, B, C, D]
for i in range(len(matrices)):
# Convert to matrix first, if necessary.
matrices[i] = matrix(matrices[i])
[A, B, C, D] = matrices

LTI.__init__(self, B.shape[1], C.shape[0], dt)
A, B, C, D = [matrix(M) for M in (A, B, C, D)]

# TODO: use super here?
LTI.__init__(self, inputs=D.shape[1], outputs=D.shape[0], dt=dt)
self.A = A
self.B = B
self.C = C
self.D = D

self.states = A.shape[0]
self.states = A.shape[1]

if 0 == self.states:
# static gain
# matrix's default "empty" shape is 1x0
A.shape = (0,0)
B.shape = (0,self.inputs)
C.shape = (self.outputs,0)

# Check that the matrix sizes are consistent.
if self.states != A.shape[1]:
if self.states != A.shape[0]:
raise ValueError("A must be square.")
if self.states != B.shape[0]:
raise ValueError("B must have the same row size as A.")
raise ValueError("A and B must have the same number of rows.")
if self.states != C.shape[1]:
raise ValueError("C must have the same column size as A.")
if self.inputs != D.shape[1]:
raise ValueError("D must have the same column size as B.")
if self.outputs != D.shape[0]:
raise ValueError("D must have the same row size as C.")
raise ValueError("A and C must have the same number of columns.")
if self.inputs != B.shape[1]:
raise ValueError("B and D must have the same number of columns.")
if self.outputs != C.shape[0]:
raise ValueError("C and D must have the same number of rows.")

# Check for states that don't do anything, and remove them.
self._remove_useless_states()
Expand Down Expand Up @@ -179,17 +180,10 @@ def _remove_useless_states(self):
useless.append(i)

# Remove the useless states.
if all(useless == range(self.states)):
# All the states were useless.
self.A = zeros((1, 1))
self.B = zeros((1, self.inputs))
self.C = zeros((self.outputs, 1))
else:
# A more typical scenario.
self.A = delete(self.A, useless, 0)
self.A = delete(self.A, useless, 1)
self.B = delete(self.B, useless, 0)
self.C = delete(self.C, useless, 1)
self.A = delete(self.A, useless, 0)
self.A = delete(self.A, useless, 1)
self.B = delete(self.B, useless, 0)
self.C = delete(self.C, useless, 1)

self.states = self.A.shape[0]
self.inputs = self.B.shape[1]
Expand Down Expand Up @@ -477,18 +471,22 @@ def feedback(self, other=1, sign=-1):
def minreal(self, tol=0.0):
"""Calculate a minimal realization, removes unobservable and
uncontrollable states"""
try:
from slycot import tb01pd
B = empty((self.states, max(self.inputs, self.outputs)))
B[:,:self.inputs] = self.B
C = empty((max(self.outputs, self.inputs), self.states))
C[:self.outputs,:] = self.C
A, B, C, nr = tb01pd(self.states, self.inputs, self.outputs,
self.A, B, C, tol=tol)
return StateSpace(A[:nr,:nr], B[:nr,:self.inputs],
C[:self.outputs,:nr], self.D)
except ImportError:
raise TypeError("minreal requires slycot tb01pd")
if self.states:
try:
from slycot import tb01pd
B = empty((self.states, max(self.inputs, self.outputs)))
B[:,:self.inputs] = self.B
C = empty((max(self.outputs, self.inputs), self.states))
C[:self.outputs,:] = self.C
A, B, C, nr = tb01pd(self.states, self.inputs, self.outputs,
self.A, B, C, tol=tol)
return StateSpace(A[:nr,:nr], B[:nr,:self.inputs],
C[:self.outputs,:nr], self.D)
except ImportError:
raise TypeError("minreal requires slycot tb01pd")
else:
return StateSpace(self)


# TODO: add discrete time check
def returnScipySignalLTI(self):
Expand Down
77 changes: 75 additions & 2 deletions control/tests/statesp_test.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,9 +5,10 @@

import unittest
import numpy as np
from scipy.linalg import eigvals
from numpy.linalg import solve
from scipy.linalg import eigvals, block_diag
from control import matlab
from control.statesp import StateSpace, _convertToStateSpace
from control.statesp import StateSpace, _convertToStateSpace,tf2ss
from control.xferfcn import TransferFunction

class TestStateSpace(unittest.TestCase):
Expand Down Expand Up @@ -235,6 +236,78 @@ def test_dcgain(self):
sys3 = StateSpace(0., 1., 1., 0.)
np.testing.assert_equal(sys3.dcgain(), np.nan)


def test_scalarStaticGain(self):
"""Regression: can we create a scalar static gain?"""
g1=StateSpace([],[],[],[2])
g2=StateSpace([],[],[],[3])

# make sure StateSpace internals, specifically ABC matrix
# sizes, are OK for LTI operations
g3 = g1*g2
self.assertEqual(6, g3.D[0,0])
g4 = g1+g2
self.assertEqual(5, g4.D[0,0])
g5 = g1.feedback(g2)
self.assertAlmostEqual(2./7, g5.D[0,0])
g6 = g1.append(g2)
np.testing.assert_array_equal(np.diag([2,3]),g6.D)

def test_matrixStaticGain(self):
"""Regression: can we create matrix static gains?"""
d1 = np.matrix([[1,2,3],[4,5,6]])
d2 = np.matrix([[7,8],[9,10],[11,12]])
g1=StateSpace([],[],[],d1)

# _remove_useless_states was making A = [[0]]
self.assertEqual((0,0), g1.A.shape)

g2=StateSpace([],[],[],d2)
g3=StateSpace([],[],[],d2.T)

h1 = g1*g2
np.testing.assert_array_equal(d1*d2, h1.D)
h2 = g1+g3
np.testing.assert_array_equal(d1+d2.T, h2.D)
h3 = g1.feedback(g2)
np.testing.assert_array_almost_equal(solve(np.eye(2)+d1*d2,d1), h3.D)
h4 = g1.append(g2)
np.testing.assert_array_equal(block_diag(d1,d2),h4.D)


def test_remove_useless_states(self):
"""Regression: _remove_useless_states gives correct ABC sizes"""
g1 = StateSpace(np.zeros((3,3)),
np.zeros((3,4)),
np.zeros((5,3)),
np.zeros((5,4)))
self.assertEqual((0,0), g1.A.shape)
self.assertEqual((0,4), g1.B.shape)
self.assertEqual((5,0), g1.C.shape)
self.assertEqual((5,4), g1.D.shape)


def test_BadEmptyMatrices(self):
"""Mismatched ABCD matrices when some are empty"""
self.assertRaises(ValueError,StateSpace, [1], [], [], [1])
self.assertRaises(ValueError,StateSpace, [1], [1], [], [1])
self.assertRaises(ValueError,StateSpace, [1], [], [1], [1])
self.assertRaises(ValueError,StateSpace, [], [1], [], [1])
self.assertRaises(ValueError,StateSpace, [], [1], [1], [1])
self.assertRaises(ValueError,StateSpace, [], [], [1], [1])
self.assertRaises(ValueError,StateSpace, [1], [1], [1], [])


def test_minrealStaticGain(self):
"""Regression: minreal on static gain was failing"""
g1 = StateSpace([],[],[],[1])
g2 = g1.minreal()
np.testing.assert_array_equal(g1.A, g2.A)
np.testing.assert_array_equal(g1.B, g2.B)
np.testing.assert_array_equal(g1.C, g2.C)
np.testing.assert_array_equal(g1.D, g2.D)


class TestRss(unittest.TestCase):
"""These are tests for the proper functionality of statesp.rss."""

Expand Down