Menu
▾
▴
[python-control-discuss] SF.net SVN: python-control:[60] branches/control-0.4a/src
|
From: <kk...@us...> - 2011-02-08 22:14:44
|
Revision: 60
http://python-control.svn.sourceforge.net/python-control/?rev=60&view=rev
Author: kkchen
Date: 2011-02-08 22:14:38 +0000 (Tue, 08 Feb 2011)
Log Message:
-----------
Implemented StateSpace.evalfr and StateSpace.freqresp in statesp.py; added
appropriate tests in TestStateSp.py.
Also made miscellaneous fixes in xferfcn.py and TestXferFcn.py.
Kevin K. Chen <kk...@pr...>
Modified Paths:
--------------
branches/control-0.4a/src/TestStateSp.py
branches/control-0.4a/src/TestXferFcn.py
branches/control-0.4a/src/statesp.py
branches/control-0.4a/src/xferfcn.py
Modified: branches/control-0.4a/src/TestStateSp.py
===================================================================
--- branches/control-0.4a/src/TestStateSp.py 2011-02-08 22:14:32 UTC (rev 59)
+++ branches/control-0.4a/src/TestStateSp.py 2011-02-08 22:14:38 UTC (rev 60)
@@ -2,8 +2,53 @@
import numpy as np
import matlab
+from statesp import StateSpace
import unittest
+class TestStateSpace(unittest.TestCase):
+ """Tests for the StateSpace class."""
+
+ def testEvalFr(self):
+ """Evaluate the frequency response at one frequency."""
+
+ A = [[-2, 0.5], [0.5, -0.3]]
+ B = [[0.3, -1.3], [0.1, 0.]]
+ C = [[0., 0.1], [-0.3, -0.2]]
+ D = [[0., -0.8], [-0.3, 0.]]
+ sys = StateSpace(A, B, C, D)
+
+ resp = [[4.37636761487965e-05 - 0.0152297592997812j,
+ -0.792603938730853 + 0.0261706783369803j],
+ [-0.331544857768052 + 0.0576105032822757j,
+ 0.128919037199125 - 0.143824945295405j]]
+
+ np.testing.assert_almost_equal(sys.evalfr(1.), resp)
+
+ def testFreqResp(self):
+ """Evaluate the frequency response at multiple frequencies."""
+
+ A = [[-2, 0.5], [0.5, -0.3]]
+ B = [[0.3, -1.3], [0.1, 0.]]
+ C = [[0., 0.1], [-0.3, -0.2]]
+ D = [[0., -0.8], [-0.3, 0.]]
+ sys = StateSpace(A, B, C, D)
+
+ truemag = [[[0.0852992637230322, 0.00103596611395218],
+ [0.935374692849736, 0.799380720864549]],
+ [[0.55656854563842, 0.301542699860857],
+ [0.609178071542849, 0.0382108097985257]]]
+ truephase = [[[-0.566195599644593, -1.68063565332582],
+ [3.0465958317514, 3.14141384339534]],
+ [[2.90457947657161, 3.10601268291914],
+ [-0.438157380501337, -1.40720969147217]]]
+ trueomega = [0.1, 10.]
+
+ mag, phase, omega = sys.freqresp(trueomega)
+
+ np.testing.assert_almost_equal(mag, truemag)
+ np.testing.assert_almost_equal(phase, truephase)
+ np.testing.assert_equal(omega, trueomega)
+
class TestRss(unittest.TestCase):
"""These are tests for the proper functionality of statesp.rss."""
@@ -73,4 +118,4 @@
self.assertTrue(abs(z) < 1)
if __name__ == "__main__":
- unittest.main()
\ No newline at end of file
+ unittest.main()
Modified: branches/control-0.4a/src/TestXferFcn.py
===================================================================
--- branches/control-0.4a/src/TestXferFcn.py 2011-02-08 22:14:32 UTC (rev 59)
+++ branches/control-0.4a/src/TestXferFcn.py 2011-02-08 22:14:38 UTC (rev 60)
@@ -377,7 +377,7 @@
np.testing.assert_array_almost_equal(mag, truemag)
np.testing.assert_array_almost_equal(phase, truephase)
- np.testing.assert_array_almost_equal(omega, trueomega)
+ np.testing.assert_array_equal(omega, trueomega)
if __name__ == "__main__":
Modified: branches/control-0.4a/src/statesp.py
===================================================================
--- branches/control-0.4a/src/statesp.py 2011-02-08 22:14:32 UTC (rev 59)
+++ branches/control-0.4a/src/statesp.py 2011-02-08 22:14:38 UTC (rev 60)
@@ -43,6 +43,7 @@
import scipy as sp
from scipy import concatenate, zeros
+from numpy.linalg import solve
import xferfcn
from lti2 import Lti2
@@ -92,34 +93,31 @@
str += "D = " + self.D.__str__() + "\n"
return str
+ def evalfr(self, freq):
+ """Method for evaluating a system at one frequency."""
+
+ fresp = self.C * solve(freq * 1.j * sp.eye(self.states) - self.A,
+ self.B) + self.D
+ return fresp
+
# Method for generating the frequency response of the system
def freqresp(self, omega=None):
"""Compute the response of a system to a list of frequencies."""
- # Generate and save a transfer function matrix
- #! TODO: This is currently limited to SISO systems
- #nout, nin = self.D.shape
+ # Preallocate outputs.
+ numfreq = len(omega)
+ mag = sp.empty((self.outputs, self.inputs, numfreq))
+ phase = sp.empty((self.outputs, self.inputs, numfreq))
+ fresp = sp.empty((self.outputs, self.inputs, numfreq), dtype=complex)
- # Compute the denominator from the A matrix
- den = sp.poly1d(sp.poly(self.A))
+ for k in range(numfreq):
+ fresp[:, :, k] = self.evalfr(omega[k])
- # Compute the numerator based on zeros
- #! TODO: This is currently limited to SISO systems
- num = sp.poly1d(\
- sp.poly(self.A - sp.dot(self.B, self.C)) + (self.D[0] - 1) * den)
-
- # Generate the frequency response at each frequency
- fresp = map(lambda w: num(w*1j) / den(w*1j), omega)
- mag = sp.sqrt(sp.multiply(fresp, sp.conjugate(fresp)))
+ mag = abs(fresp)
phase = sp.angle(fresp)
return mag, phase, omega
- # Method for evaluating a system at one frequency
- def evalfr(self, freq):
- #! TODO: Not implemented
- return None
-
# Compute poles and zeros
def poles(self):
return sp.roots(sp.poly(self.A))
Modified: branches/control-0.4a/src/xferfcn.py
===================================================================
--- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:14:32 UTC (rev 59)
+++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:14:38 UTC (rev 60)
@@ -145,8 +145,6 @@
if mimo:
outstr += "\nInput %i to output %i:" % (i + 1, j + 1)
- lablen = 0
-
# Convert the numerator and denominator polynomials to strings.
numstr = _tfpolyToString(self.num[j][i]);
denstr = _tfpolyToString(self.den[j][i]);
@@ -158,11 +156,11 @@
# Center the numerator or denominator
if (len(numstr) < dashcount):
numstr = ' ' * \
- int(round((dashcount - len(numstr))/2) + lablen) + \
+ int(round((dashcount - len(numstr))/2)) + \
numstr
if (len(denstr) < dashcount):
denstr = ' ' * \
- int(round((dashcount - len(denstr))/2) + lablen) + \
+ int(round((dashcount - len(denstr))/2)) + \
denstr
outstr += "\n" + numstr + "\n" + dashes + "\n" + denstr + "\n"
@@ -229,7 +227,7 @@
return xTransferFunction(num, den)
def __radd__(self, other):
- """Add two transfer functions (parallel connection)"""
+ """Add two transfer functions (parallel connection)."""
return self + other;
@@ -244,7 +242,7 @@
return other + (-self)
def __mul__(self, other):
- """Multiply two transfer functions (serial connection)"""
+ """Multiply two transfer functions (serial connection)."""
# Convert the second argument to a transfer function.
if not isinstance(other, xTransferFunction):
@@ -278,13 +276,13 @@
return xTransferFunction(num, den)
def __rmul__(self, other):
- """Multiply two transfer functions (serial connection)"""
+ """Multiply two transfer functions (serial connection)."""
return self * other
# TODO: Division of MIMO transfer function objects is quite difficult.
def __div__(self, other):
- """Divide two transfer functions"""
+ """Divide two transfer functions."""
if self.inputs > 1 or self.outputs > 1 or \
other.inputs > 1 or other.outputs > 1:
@@ -302,7 +300,7 @@
# TODO: Division of MIMO transfer function objects is quite difficult.
def __rdiv__(self, other):
- """Reverse divide two transfer functions"""
+ """Reverse divide two transfer functions."""
if self.inputs > 1 or self.outputs > 1 or \
other.inputs > 1 or other.outputs > 1:
@@ -312,7 +310,7 @@
return other / self
def evalfr(self, freq):
- """Evaluate a transfer function at a single frequency"""
+ """Evaluate a transfer function at a single frequency."""
# Preallocate the output.
out = sp.empty((self.outputs, self.inputs), dtype=complex)
@@ -325,23 +323,22 @@
return out
# Method for generating the frequency response of the system
- def freqresp(self, omega):
- """Evaluate a transfer function at a list of frequencies"""
+ def freqresp(self, omega=None):
+ """Evaluate a transfer function at a list of frequencies."""
+ # Preallocate outputs.
numfreq = len(omega)
+ mag = sp.empty((self.outputs, self.inputs, numfreq))
+ phase = sp.empty((self.outputs, self.inputs, numfreq))
- # Preallocate outputs.
- mag = sp.empty((self.outputs, self.inputs, numfreq), dtype=complex)
- phase = sp.empty((self.outputs, self.inputs, numfreq), dtype=complex)
-
for i in range(self.outputs):
for j in range(self.inputs):
fresp = map(lambda w: sp.polyval(self.num[i][j], w * 1.j) / \
sp.polyval(self.den[i][j], w * 1.j), omega)
fresp = sp.array(fresp)
- mag[i][j] = abs(fresp)
- phase[i][j] = sp.angle(fresp)
+ mag[i, j] = abs(fresp)
+ phase[i, j] = sp.angle(fresp)
return mag, phase, omega
@@ -356,7 +353,7 @@
pass
def feedback(sys1, sys2, sign=-1):
- """Feedback interconnection between two transfer functions"""
+ """Feedback interconnection between two transfer functions."""
pass
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
