You can subscribe to this list here.
2010 |
Jan
|
Feb
|
Mar
|
Apr
|
May
|
Jun
(19) |
Jul
|
Aug
|
Sep
|
Oct
|
Nov
(11) |
Dec
(5) |
---|---|---|---|---|---|---|---|---|---|---|---|---|
2011 |
Jan
|
Feb
(113) |
Mar
(12) |
Apr
(27) |
May
(2) |
Jun
(16) |
Jul
(6) |
Aug
(6) |
Sep
|
Oct
(3) |
Nov
(9) |
Dec
(2) |
2012 |
Jan
(5) |
Feb
(11) |
Mar
|
Apr
(3) |
May
|
Jun
|
Jul
|
Aug
(3) |
Sep
(7) |
Oct
(18) |
Nov
(18) |
Dec
|
2013 |
Jan
(4) |
Feb
(1) |
Mar
(3) |
Apr
(1) |
May
|
Jun
(33) |
Jul
(2) |
Aug
(5) |
Sep
|
Oct
|
Nov
|
Dec
(1) |
2014 |
Jan
(1) |
Feb
|
Mar
(8) |
Apr
|
May
(3) |
Jun
(3) |
Jul
(9) |
Aug
(5) |
Sep
(6) |
Oct
|
Nov
|
Dec
|
2015 |
Jan
|
Feb
|
Mar
|
Apr
|
May
|
Jun
|
Jul
|
Aug
|
Sep
|
Oct
|
Nov
|
Dec
(4) |
2017 |
Jan
|
Feb
|
Mar
|
Apr
|
May
|
Jun
|
Jul
|
Aug
|
Sep
|
Oct
(5) |
Nov
|
Dec
|
2018 |
Jan
|
Feb
|
Mar
|
Apr
|
May
|
Jun
(1) |
Jul
|
Aug
|
Sep
|
Oct
|
Nov
|
Dec
|
2019 |
Jan
|
Feb
|
Mar
|
Apr
|
May
|
Jun
|
Jul
|
Aug
(6) |
Sep
|
Oct
|
Nov
(2) |
Dec
|
2020 |
Jan
(1) |
Feb
(1) |
Mar
|
Apr
|
May
|
Jun
|
Jul
|
Aug
|
Sep
|
Oct
|
Nov
|
Dec
|
2021 |
Jan
|
Feb
|
Mar
|
Apr
|
May
|
Jun
|
Jul
|
Aug
|
Sep
(4) |
Oct
|
Nov
|
Dec
|
From: <kk...@us...> - 2011-02-08 22:17:19
|
Revision: 91 http://python-control.svn.sourceforge.net/python-control/?rev=91&view=rev Author: kkchen Date: 2011-02-08 22:17:13 +0000 (Tue, 08 Feb 2011) Log Message: ----------- freqresp now sorts omega, to be more like MATLAB. Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/matlab.py branches/control-0.4a/src/statesp.py branches/control-0.4a/src/xferfcn.py Modified: branches/control-0.4a/src/matlab.py =================================================================== --- branches/control-0.4a/src/matlab.py 2011-02-08 22:17:07 UTC (rev 90) +++ branches/control-0.4a/src/matlab.py 2011-02-08 22:17:13 UTC (rev 91) @@ -646,7 +646,7 @@ Parameters ---------- sys: StateSpace or TransferFunction object - omega: list, tuple, or ndarray + omega: list or ndarray Returns ------- @@ -662,7 +662,8 @@ Notes ----- This function is a wrapper for StateSpace.freqresp and - TransferFunction.freqresp. + TransferFunction.freqresp. The output omega is a sorted version of the + input omega. Examples -------- Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:17:07 UTC (rev 90) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:17:13 UTC (rev 91) @@ -301,7 +301,8 @@ reports the value of the magnitude, phase, and angular frequency of the system's transfer function matrix evaluated at s = i * omega, where - omega is a list of angular frequencies. + omega is a list of angular frequencies, and is a sorted version of the + input omega. """ @@ -311,6 +312,8 @@ phase = empty((self.outputs, self.inputs, numfreq)) fresp = empty((self.outputs, self.inputs, numfreq), dtype=complex) + omega.sort() + for k in range(numfreq): fresp[:, :, k] = self.evalfr(omega[k]) Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:17:07 UTC (rev 90) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:17:13 UTC (rev 91) @@ -386,7 +386,7 @@ reports the value of the magnitude, phase, and angular frequency of the transfer function matrix evaluated at s = i * omega, where omega is a - list of angular frequencies. + list of angular frequencies, and is a sorted version of the input omega. """ @@ -395,6 +395,8 @@ mag = empty((self.outputs, self.inputs, numfreq)) phase = empty((self.outputs, self.inputs, numfreq)) + omega.sort() + for i in range(self.outputs): for j in range(self.inputs): fresp = map(lambda w: (polyval(self.num[i][j], w * 1.j) / This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:17:13
|
Revision: 90 http://python-control.svn.sourceforge.net/python-control/?rev=90&view=rev Author: kkchen Date: 2011-02-08 22:17:07 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Deleted a print statement used while debugging. Lauren Padilla <lpa...@pr...> Modified Paths: -------------- branches/control-0.4a/src/TestConvert.py branches/control-0.4a/src/statesp.py Modified: branches/control-0.4a/src/TestConvert.py =================================================================== --- branches/control-0.4a/src/TestConvert.py 2011-02-08 22:17:03 UTC (rev 89) +++ branches/control-0.4a/src/TestConvert.py 2011-02-08 22:17:07 UTC (rev 90) @@ -27,10 +27,10 @@ print "sys1:\n", sys sys2 = matlab.tf(sys) print "sys2:\n", sys2 - sys3 = matlab.ss(sys2) - print "sys3:\n", sys3 - sys4 = matlab.tf(sys3) - print "sys4:\n", sys4 + #sys3 = matlab.ss(sys2) + #print "sys3:\n", sys3 + #sys4 = matlab.tf(sys3) + #print "sys4:\n", sys4 if __name__ == "__main__": unittest.main() Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:17:03 UTC (rev 89) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:17:07 UTC (rev 90) @@ -437,7 +437,6 @@ # Repeat the common denominator along the rows. den = array([den for i in range(sys.outputs)]) - print "outputs = %g\n" % sys.outputs ssout = td04ad(sys.inputs, sys.outputs, index, den, num) states = ssout[0] This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:17:09
|
Revision: 89 http://python-control.svn.sourceforge.net/python-control/?rev=89&view=rev Author: kkchen Date: 2011-02-08 22:17:03 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Bug fixes in convertToStateSpace; new TestConvert.py Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/TestStateSp.py branches/control-0.4a/src/statesp.py Added Paths: ----------- branches/control-0.4a/src/TestConvert.py Added: branches/control-0.4a/src/TestConvert.py =================================================================== --- branches/control-0.4a/src/TestConvert.py (rev 0) +++ branches/control-0.4a/src/TestConvert.py 2011-02-08 22:17:03 UTC (rev 89) @@ -0,0 +1,36 @@ +#!/usr/bin/env python + +import numpy as np +import matlab +import unittest + +class TestConvert(unittest.TestCase): + """Test state space and transfer function conversions.""" + + def setUp(self): + """Set up testing parameters.""" + + # Number of times to run each of the randomized tests. + self.numTests = 10 + # Maximum number of states to test + 1 + self.maxStates = 5 + # Maximum number of inputs and outputs to test + 1 + self.maxIO = 5 + + def testConvert(self): + """Test state space to transfer function conversion.""" + + for states in range(1, self.maxStates): + for inputs in range(1, self.maxIO): + for outputs in range(1, self.maxIO): + sys = matlab.rss(states, inputs, outputs) + print "sys1:\n", sys + sys2 = matlab.tf(sys) + print "sys2:\n", sys2 + sys3 = matlab.ss(sys2) + print "sys3:\n", sys3 + sys4 = matlab.tf(sys3) + print "sys4:\n", sys4 + +if __name__ == "__main__": + unittest.main() Property changes on: branches/control-0.4a/src/TestConvert.py ___________________________________________________________________ Added: svn:executable + * Modified: branches/control-0.4a/src/TestStateSp.py =================================================================== --- branches/control-0.4a/src/TestStateSp.py 2011-02-08 22:16:57 UTC (rev 88) +++ branches/control-0.4a/src/TestStateSp.py 2011-02-08 22:17:03 UTC (rev 89) @@ -1,9 +1,9 @@ #!/usr/bin/env python import numpy as np +import unittest import matlab from statesp import StateSpace -import unittest class TestStateSpace(unittest.TestCase): """Tests for the StateSpace class.""" Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:16:57 UTC (rev 88) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:17:03 UTC (rev 89) @@ -433,13 +433,18 @@ # function matrix has a common denominator. num, den = sys._common_den() # Make a list of the orders of the denominator polynomials. - index = [len(den) for i in range(sys.outputs)] + index = [len(den) - 1 for i in range(sys.outputs)] # Repeat the common denominator along the rows. den = array([den for i in range(sys.outputs)]) + print "outputs = %g\n" % sys.outputs ssout = td04ad(sys.inputs, sys.outputs, index, den, num) - return StateSpace(ssout[1], ssout[2], ssout[3], ssout[4]) + states = ssout[0] + return StateSpace(ssout[1][:states, :states], + ssout[2][:states, :sys.inputs], + ssout[3][:sys.outputs, :states], + ssout[4]) elif isinstance(sys, (int, long, float, complex)): if "inputs" in kw: inputs = kw["inputs"] This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:17:03
|
Revision: 88 http://python-control.svn.sourceforge.net/python-control/?rev=88&view=rev Author: kkchen Date: 2011-02-08 22:16:57 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Major bug fix in convertToStateSpace. Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/statesp.py Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:16:53 UTC (rev 87) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:16:57 UTC (rev 88) @@ -437,7 +437,7 @@ # Repeat the common denominator along the rows. den = array([den for i in range(sys.outputs)]) - ssout = td04ad(sys.inputs, sys.outputs, index, num, den) + ssout = td04ad(sys.inputs, sys.outputs, index, den, num) return StateSpace(ssout[1], ssout[2], ssout[3], ssout[4]) elif isinstance(sys, (int, long, float, complex)): This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:59
|
Revision: 87 http://python-control.svn.sourceforge.net/python-control/?rev=87&view=rev Author: kkchen Date: 2011-02-08 22:16:53 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Modified a few of the docstrings in files so that Sphinx autodocumenation runs. Lauren Padilla <lpa...@pr...> Modified Paths: -------------- branches/control-0.4a/src/ctrlutil.py branches/control-0.4a/src/matlab.py branches/control-0.4a/src/modelsimp.py Modified: branches/control-0.4a/src/ctrlutil.py =================================================================== --- branches/control-0.4a/src/ctrlutil.py 2011-02-08 22:16:47 UTC (rev 86) +++ branches/control-0.4a/src/ctrlutil.py 2011-02-08 22:16:53 UTC (rev 87) @@ -43,23 +43,24 @@ # Packages that we need access to import scipy as sp import statesp -import xferfcn +import xferfcn # Specific functions that we use from scipy import pi # Utility function to unwrap an angle measurement + def unwrap(angle, period=2*pi): """Unwrap a phase angle to give a continuous curve - Usage: Y = unwrap(X, period=2*pi) + Usage: Y = unwrap(X, period=2``*``pi) Parameters ---------- X : array_like Input array period : number - Input period (usually either 2*pi or 360) + Input period (usually either 2``*``pi or 360) Returns ------- Modified: branches/control-0.4a/src/matlab.py =================================================================== --- branches/control-0.4a/src/matlab.py 2011-02-08 22:16:47 UTC (rev 86) +++ branches/control-0.4a/src/matlab.py 2011-02-08 22:16:53 UTC (rev 87) @@ -8,7 +8,7 @@ (python-control) for people who are familiar with the MATLAB Control Systems Toolbox (tm). Most of the functions are just calls to python-control functions defined elsewhere. Use 'from -control.matlab import *' in python to include all of the functions +control.matlab import \*' in python to include all of the functions defined here. Functions that are defined in other libraries that have the same names as their MATLAB equivalents are automatically imported here. @@ -67,7 +67,7 @@ import freqplot from statesp import StateSpace, rss_generate, convertToStateSpace from xferfcn import TransferFunction, convertToTransferFunction -from exception import * +from exception import ControlArgument # Import MATLAB-like functions that can be used as-is from ctrlutil import unwrap @@ -80,13 +80,13 @@ __doc__ += """ The control.matlab module defines functions that are roughly the equivalents of those in the MATLAB Control Toolbox. Items marked by a -'*' are currently implemented; those marked with a '-' are not planned +\* are currently implemented; those marked with a ``-`` are not planned for implementation. Creating linear models. -* tf - create transfer function (TF) models +\* tf - create transfer function (TF) models zpk - create zero/pole/gain (ZPK) models. -* ss - create state-space (SS) models +\* ss - create state-space (SS) models dss - create descriptor state-space models delayss - create state-space models with delayed terms frd - create frequency response data (FRD) models @@ -105,26 +105,26 @@ ss/getDelayModel - access internal delay model (state space only) Conversions -* tf - conversion to transfer function +\* tf - conversion to transfer function zpk - conversion to zero/pole/gain -* ss - conversion to state space +\* ss - conversion to state space frd - conversion to frequency data c2d - continuous to discrete conversion d2c - discrete to continuous conversion d2d - resample discrete-time model upsample - upsample discrete-time LTI systems -* ss2tf - state space to transfer function +\* ss2tf - state space to transfer function ss2zpk - transfer function to zero-pole-gain -* tf2ss - transfer function to state space +\* tf2ss - transfer function to state space tf2zpk - transfer function to zero-pole-gain zpk2ss - zero-pole-gain to state space zpk2tf - zero-pole-gain to transfer function System interconnections append - group LTI models by appending inputs and outputs -* parallel - connect LTI models in parallel (see also overloaded +) -* series - connect LTI models in series (see also overloaded *) -* feedback - connect lti models with a feedback loop +\* parallel - connect LTI models in parallel (see also overloaded +) +\* series - connect LTI models in series (see also overloaded \*) +\* feedback - connect lti models with a feedback loop lti/lft - generalized feedback interconnection lti/connect - arbitrary interconnection of lti models sumblk - specify summing junction (for use with connect) @@ -134,10 +134,10 @@ dcgain - steady-state (D.C.) gain lti/bandwidth - system bandwidth lti/norm - h2 and Hinfinity norms of LTI models -* lti/pole - system poles +\* lti/pole - system poles lti/zero - system (transmission) zeros lti/order - model order (number of states) -* pzmap - pole-zero map +\* pzmap - pole-zero map lti/iopzmap - input/output pole-zero map damp - natural frequency and damping of system poles esort - sort continuous poles by real part @@ -146,25 +146,25 @@ lti/modsep - region-based modal decomposition Time-domain analysis -* step - step response +\* step - step response stepinfo - step response characteristics (rise time, ...) -* impulse - impulse response +\* impulse - impulse response initial - free response with initial conditions -* lsim - response to user-defined input signal +\* lsim - response to user-defined input signal lsiminfo - linear response characteristics gensig - generate input signal for LSIM covar - covariance of response to white noise Frequency-domain analysis -* bode - Bode plot of the frequency response +\* bode - Bode plot of the frequency response lti/bodemag - Bode magnitude diagram only sigma - singular value frequency plot -* nyquist - Nyquist plot -* nichols - Nichols plot +\* nyquist - Nyquist plot +\* nichols - Nichols plot margin - gain and phase margins lti/allmargin - all crossover frequencies and related gain/phase margins -* lti/freqresp - frequency response over a frequency grid -* lti/evalfr - evaluate frequency response at given frequency +\* lti/freqresp - frequency response over a frequency grid +\* lti/evalfr - evaluate frequency response at given frequency Model simplification minreal - minimal realization and pole/zero cancellation @@ -192,13 +192,13 @@ augstate - augment output by appending states State-space (SS) models -* rss - random stable continuous-time state-space models -* drss - random stable discrete-time state-space models +\* rss - random stable continuous-time state-space models +\* drss - random stable discrete-time state-space models ss2ss - state coordinate transformation canon - canonical forms of state-space models -* ctrb - controllability matrix -* obsv - observability matrix -* gram - controllability and observability gramians +\* ctrb - controllability matrix +\* obsv - observability matrix +\* gram - controllability and observability gramians ss/prescale - optimal scaling of state-space models. balreal - gramian-based input/output balancing ss/xperm - reorder states. @@ -219,7 +219,7 @@ lti/hasdelay - true for models with time delays lti/totaldelay - total delay between each input/output pair lti/delay2z - replace delays by poles at z=0 or FRD phase shift -* pade - pade approximation of time delays +\* pade - pade approximation of time delays Model dimensions and characteristics class - model type ('tf', 'zpk', 'ss', or 'frd') @@ -235,14 +235,14 @@ lti/reshape - reshape array of linear models Overloaded arithmetic operations -* + and - - add and subtract systems (parallel connection) -* * - multiply systems (series connection) -* / - left divide -- sys1\sys2 means inv(sys1)*sys2 -* / - right divide -- sys1/sys2 means sys1*inv(sys2) +\* + and - - add and subtract systems (parallel connection) +\* \* - multiply systems (series connection) +\* / - left divide -- sys1\sys2 means inv(sys1)\*sys2 +\* / - right divide -- sys1/sys2 means sys1\*inv(sys2) ^ - powers of a given system ' - pertransposition .' - transposition of input/output map - .* - element-by-element multiplication + .\* - element-by-element multiplication [..] - concatenate models along inputs or outputs lti/stack - stack models/arrays along some array dimension lti/inv - inverse of an LTI system @@ -256,10 +256,10 @@ bdschur - block diagonalization of a square matrix Additional functions -* gangof4 - generate the Gang of 4 sensitivity plots -* linspace - generate a set of numbers that are linearly spaced -* logspace - generate a set of numbers that are logarithmically spaced -* unwrap - unwrap a phase angle to give a continuous curve +\* gangof4 - generate the Gang of 4 sensitivity plots +\* linspace - generate a set of numbers that are linearly spaced +\* logspace - generate a set of numbers that are logarithmically spaced +\* unwrap - unwrap a phase angle to give a continuous curve """ Modified: branches/control-0.4a/src/modelsimp.py =================================================================== --- branches/control-0.4a/src/modelsimp.py 2011-02-08 22:16:47 UTC (rev 86) +++ branches/control-0.4a/src/modelsimp.py 2011-02-08 22:16:53 UTC (rev 87) @@ -42,7 +42,7 @@ # External packages and modules import numpy as np import ctrlutil -from control.exception import * +from exception import * from statefbk import * from statesp import StateSpace This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:53
|
Revision: 86 http://python-control.svn.sourceforge.net/python-control/?rev=86&view=rev Author: kkchen Date: 2011-02-08 22:16:47 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Modified convertTo* to use keywords. Tests for TransferFunction.pole. The optional arguments inputs, outputs for convertTo{StateSpace,TransferFunction} are now keyword arguments. This lets us make sure the user isn't passing the number of inputs and outputs when it shouldn't. A test for TransferFunction.pole passes. Kevin K. Chen <kk...@pr...> Modified Paths: -------------- 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/TestXferFcn.py =================================================================== --- branches/control-0.4a/src/TestXferFcn.py 2011-02-08 22:16:42 UTC (rev 85) +++ branches/control-0.4a/src/TestXferFcn.py 2011-02-08 22:16:47 UTC (rev 86) @@ -381,6 +381,17 @@ np.testing.assert_array_almost_equal(phase, truephase) np.testing.assert_array_equal(omega, trueomega) + # Tests for TransferFunction.pole and TransferFunction.zero. + + def testPoleMIMO(self): + """Test for correct MIMO poles.""" + + sys = TransferFunction([[[1.], [1.]], [[1.], [1.]]], + [[[1., 2.], [1., 3.]], [[1., 4., 4.], [1., 9., 14.]]]) + p = sys.pole() + + np.testing.assert_array_almost_equal(p, [-7., -3., -2., -2.]) + # Tests for TransferFunction.feedback. def testFeedbackSISO(self): Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:16:42 UTC (rev 85) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:16:47 UTC (rev 86) @@ -401,20 +401,34 @@ return out -def convertToStateSpace(sys, inputs=1, outputs=1): +def convertToStateSpace(sys, **kw): """Convert a system to state space form (if needed). If sys is already a state space, then it is returned. If sys is a transfer function object, then it is converted to a state space and returned. If sys is a scalar, then the number of inputs and outputs can be specified - manually. + manually, as in: + + >>> sys = convertToStateSpace(3.) # Assumes inputs = outputs = 1 + >>> sys = convertToStateSpace(1., inputs=3, outputs=2) + + In the latter example, A = B = C = 0 and D = [[1., 1., 1.] + [1., 1., 1.]]. """ if isinstance(sys, StateSpace): + if len(kw): + raise TypeError("If sys is a StateSpace, convertToStateSpace \ +cannot take keywords.") + # Already a state space system; just return it return sys elif isinstance(sys, xferfcn.TransferFunction): + if len(kw): + raise TypeError("If sys is a TransferFunction, convertToStateSpace \ +cannot take keywords.") + # Change the numerator and denominator arrays so that the transfer # function matrix has a common denominator. num, den = sys._common_den() @@ -427,6 +441,15 @@ return StateSpace(ssout[1], ssout[2], ssout[3], ssout[4]) elif isinstance(sys, (int, long, float, complex)): + if "inputs" in kw: + inputs = kw["inputs"] + else: + inputs = 1 + if "outputs" in kw: + outputs = kw["outputs"] + else: + outputs = 1 + # Generate a simple state space system of the desired dimension # The following Doesn't work due to inconsistencies in ltisys: # return StateSpace([[]], [[]], [[]], eye(outputs, inputs)) Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:42 UTC (rev 85) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:47 UTC (rev 86) @@ -253,7 +253,8 @@ # Convert the second argument to a transfer function. if not isinstance(other, TransferFunction): - other = convertToTransferFunction(other, self.inputs, self.outputs) + other = convertToTransferFunction(other, inputs=self.inputs, + outputs=self.outputs) # Check that the input-output sizes are consistent. if self.inputs != other.inputs: @@ -294,7 +295,8 @@ # Convert the second argument to a transfer function. if isinstance(other, (int, float, long, complex)): - other = convertToTransferFunction(other, self.inputs, self.inputs) + other = convertToTransferFunction(other, inputs=self.inputs, + outputs=self.inputs) else: other = convertToTransferFunction(other) @@ -413,7 +415,7 @@ def zero(self): """Compute the zeros of a transfer function.""" - if (self.inputs > 1 or self.outputs > 1): + if self.inputs > 1 or self.outputs > 1: raise NotImplementedError("TransferFunction.zero is currently \ only implemented for SISO systems.") else: @@ -638,19 +640,33 @@ return num, den -def convertToTransferFunction(sys, inputs=1, outputs=1): +def convertToTransferFunction(sys, **kw): """Convert a system to transfer function form (if needed). If sys is already a transfer function, then it is returned. If sys is a state space object, then it is converted to a transfer function and returned. If sys is a scalar, then the number of inputs and outputs can be - specified manually. + specified manually, as in: + + >>> sys = convertToTransferFunction(3.) # Assumes inputs = outputs = 1 + >>> sys = convertToTransferFunction(1., inputs=3, outputs=2) + + In the latter example, sys's matrix transfer function is [[1., 1., 1.] + [1., 1., 1.]]. """ if isinstance(sys, TransferFunction): + if len(kw): + raise TypeError("If sys is a TransferFunction, \ +convertToTransferFunction cannot take keywords.") + return sys elif isinstance(sys, statesp.StateSpace): + if len(kw): + raise TypeError("If sys is a StateSpace, convertToTransferFunction \ +cannot take keywords.") + # Use Slycot to make the transformation. tfout = tb04ad(sys.states, sys.inputs, sys.outputs, sys.A, sys.B, sys.C, sys.D, sys.outputs, sys.outputs, sys.inputs) @@ -667,6 +683,15 @@ return TransferFunction(num, den) elif isinstance(sys, (int, long, float, complex)): + if "inputs" in kw: + inputs = kw["inputs"] + else: + inputs = 1 + if "outputs" in kw: + outputs = kw["outputs"] + else: + outputs = 1 + num = [[[sys] for j in range(inputs)] for i in range(outputs)] den = [[[1] for j in range(inputs)] for i in range(outputs)] This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:49
|
Revision: 85 http://python-control.svn.sourceforge.net/python-control/?rev=85&view=rev Author: kkchen Date: 2011-02-08 22:16:42 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Bode for SISO systems now compatible with new statespace and transferfunction classes. Added a script for testing frequency response functions. Lauren Padilla <lpa...@pr...> Modified Paths: -------------- branches/control-0.4a/src/ctrlutil.py branches/control-0.4a/src/freqplot.py branches/control-0.4a/src/xferfcn.py Added Paths: ----------- branches/control-0.4a/src/TestFreqRsp.py Added: branches/control-0.4a/src/TestFreqRsp.py =================================================================== --- branches/control-0.4a/src/TestFreqRsp.py (rev 0) +++ branches/control-0.4a/src/TestFreqRsp.py 2011-02-08 22:16:42 UTC (rev 85) @@ -0,0 +1,52 @@ +#!/usr/bin/env python + +# Script to test frequency response and frequency response plots like bode, nyquist and gang of 4. +# Especially need to ensure that nothing SISO is broken and that MIMO at least handles exceptions and has some default to SISO in place. + +from statesp import StateSpace +from matlab import ss, tf, bode +import numpy as np +import matplotlib.pyplot as plt + +# SISO +plt.close('all') + +A = np.matrix('1,1;0,1') +B = np.matrix('0;1') +C = np.matrix('1,0') +D = 0 +sys = StateSpace(A,B,C,D) +#or try +#sys = ss(A,B,C,D) + +# test frequency response +omega = np.linspace(10e-2,10e2,1000) +frq=sys.freqresp(omega) + +# MIMO +B = np.matrix('1,0;0,1') +D = np.matrix('0,0') +sysMIMO = ss(A,B,C,D) +frqMIMO = sysMIMO.freqresp(omega) + +plt.figure(1) +bode(sys) + +systf = tf(sys) +tfMIMO = tf(sysMIMO) + +print systf.pole() +#print tfMIMO.pole() # - should throw not implemented exception +#print tfMIMO.zero() # - should throw not implemented exception + +plt.figure(2) +bode(systf) + +#bode(sysMIMO) # - should throw not implemented exception +#bode(tfMIMO) # - should throw not implemented exception + +#plt.figure(3) +#plt.semilogx(omega,20*np.log10(np.squeeze(frq[0]))) + +#plt.figure(4) +#bode(sysMIMO,omega) \ No newline at end of file Modified: branches/control-0.4a/src/ctrlutil.py =================================================================== --- branches/control-0.4a/src/ctrlutil.py 2011-02-08 22:16:36 UTC (rev 84) +++ branches/control-0.4a/src/ctrlutil.py 2011-02-08 22:16:42 UTC (rev 85) @@ -68,10 +68,11 @@ """ wrap = 0; last = angle[0]; + for k in range(len(angle)): # See if we need to account for angle wrapping if (angle[k] - last > period/2): - wrap = -period + wrap = -period elif (last - angle[k] > period/2): wrap = period Modified: branches/control-0.4a/src/freqplot.py =================================================================== --- branches/control-0.4a/src/freqplot.py 2011-02-08 22:16:36 UTC (rev 84) +++ branches/control-0.4a/src/freqplot.py 2011-02-08 22:16:42 UTC (rev 85) @@ -88,50 +88,56 @@ if (not getattr(syslist, '__iter__', False)): syslist = (syslist,) - # Select a default range if none is provided - if (omega == None): - omega = default_frequency_range(syslist) - for sys in syslist: - # Get the magnitude and phase of the system - mag, phase, omega = sys.freqresp(omega) - if Hz: omega = omega/(2*sp.pi) - if dB: mag = 20*sp.log10(mag) - phase = unwrap(phase*180/sp.pi, 360) + if (sys.inputs > 1 or sys.outputs > 1): + #TODO: Add MIMO bode plots. + raise NotImplementedError("Bode is currently only implemented for SISO systems.") + else: + # Select a default range if none is provided + if (omega == None): + omega = default_frequency_range(syslist) - # Get the dimensions of the current axis, which we will divide up - #! TODO: Not current implemented; just use subplot for now + # Get the magnitude and phase of the system + mag_tmp, phase_tmp, omega = sys.freqresp(omega) + mag = np.squeeze(mag_tmp) + phase = np.squeeze(phase_tmp) + if Hz: omega = omega/(2*sp.pi) + if dB: mag = 20*sp.log10(mag) + phase = unwrap(phase*180/sp.pi, 360) - # Magnitude plot - plt.subplot(211); - if dB: - plt.semilogx(omega, mag) - plt.ylabel("Magnitude (dB)") - else: - plt.loglog(omega, mag) - plt.ylabel("Magnitude") + # Get the dimensions of the current axis, which we will divide up + #! TODO: Not current implemented; just use subplot for now - # Add a grid to the plot - plt.grid(True) - plt.grid(True, which='minor') - plt.hold(True); + # Magnitude plot + plt.subplot(211); + if dB: + plt.semilogx(omega, mag) + plt.ylabel("Magnitude (dB)") + else: + plt.loglog(omega, mag) + plt.ylabel("Magnitude") - # Phase plot - plt.subplot(212); - plt.semilogx(omega, phase) - plt.hold(True) + # Add a grid to the plot + plt.grid(True) + plt.grid(True, which='minor') + plt.hold(True); - # Add a grid to the plot - plt.grid(True) - plt.grid(True, which='minor') - plt.ylabel("Phase (deg)") + # Phase plot + plt.subplot(212); + plt.semilogx(omega, phase) + plt.hold(True) - # Label the frequency axis - if Hz: - plt.xlabel("Frequency (Hz)") - else: - plt.xlabel("Frequency (rad/sec)") + # Add a grid to the plot + plt.grid(True) + plt.grid(True, which='minor') + plt.ylabel("Phase (deg)") + # Label the frequency axis + if Hz: + plt.xlabel("Frequency (Hz)") + else: + plt.xlabel("Frequency (rad/sec)") + return (211, 212) # Nyquist plot @@ -316,8 +322,8 @@ features = np.array(()) for sys in syslist: # Add new features to the list - features = np.concatenate((features, np.abs(sys.poles))) - features = np.concatenate((features, np.abs(sys.zeros))) + features = np.concatenate((features, np.abs(sys.pole()))) + features = np.concatenate((features, np.abs(sys.zero()))) # Get rid of poles and zeros at the origin features = features[features != 0]; Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:36 UTC (rev 84) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:42 UTC (rev 85) @@ -413,8 +413,12 @@ def zero(self): """Compute the zeros of a transfer function.""" - raise NotImplementedError("TransferFunction.zero is not implemented \ -yet.") + if (self.inputs > 1 or self.outputs > 1): + raise NotImplementedError("TransferFunction.zero is currently \ +only implemented for SISO systems.") + else: + #for now, just give zeros of a SISO tf + return roots(self.num[0][0]) def feedback(self, other, sign=-1): """Feedback interconnection between two LTI objects.""" This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:43
|
Revision: 84 http://python-control.svn.sourceforge.net/python-control/?rev=84&view=rev Author: kkchen Date: 2011-02-08 22:16:36 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Finished writing TransferFunction._common_den. tf2ss should now work. In initial tests, however, a "nan" appears in some results. This needs to be investigated further. Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/statesp.py branches/control-0.4a/src/xferfcn.py Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:16:31 UTC (rev 83) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:16:36 UTC (rev 84) @@ -424,7 +424,7 @@ den = array([den for i in range(sys.outputs)]) ssout = td04ad(sys.inputs, sys.outputs, index, num, den) - + return StateSpace(ssout[1], ssout[2], ssout[3], ssout[4]) elif isinstance(sys, (int, long, float, complex)): # Generate a simple state space system of the desired dimension Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:31 UTC (rev 83) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:36 UTC (rev 84) @@ -72,8 +72,8 @@ """ # External function declarations -from numpy import angle, array, empty, ndarray, ones, polyadd, polymul, \ - polyval, roots, zeros +from numpy import angle, array, empty, finfo, insert, ndarray, ones, polyadd, \ + polymul, polyval, roots, sort, zeros from scipy.signal import lti from copy import deepcopy from slycot import tb04ad @@ -470,33 +470,108 @@ >>> n, d = sys._common_den() computes the single denominator containing all the poles of sys.den, and - reports it as the array d. + reports it as the array d. The output numerator array n is modified to + use the common denominator. It is an sys.outputs-by-sys.inputs-by- + [something] array. - The output numerator array n is modified to use the common denominator. - It is an sys.outputs-by-sys.inputs-by-[something] array. - """ - # Preallocate some variables. Start by figuring out the maximum number - # of numerator coefficients. - numcoeffs = 0 - for i in range(self.outputs): - for j in range(self.inputs): - numcoeffs = max(numcoeffs, len(self.num[i][j])) - # The output 3-D adjusted numerator array. - num = empty((i, j, numcoeffs)) - # A list to keep track of roots found as we scan self.den. + # A sorted list to keep track of cumulative poles found as we scan + # self.den. poles = [] - # A 3-D list to keep track of common denominator roots not present in + # A 3-D list to keep track of common denominator poles not present in # the self.den[i][j]. missingpoles = [[[] for j in range(self.inputs)] for i in range(self.outputs)] - for i in range(sys.outputs): - for j in range(sys.inputs): - currentpoles = roots(self.den[i][j]) - #TODO: finish this + for i in range(self.outputs): + for j in range(self.inputs): + # A sorted array of the poles of this SISO denominator. + currentpoles = sort(roots(self.den[i][j])) + cp_ind = 0 # Index in currentpoles. + p_ind = 0 # Index in poles. + + # Crawl along the list of current poles and the list of + # cumulative poles, until one of them reaches the end. Keep in + # mind that both lists are always sorted. + while cp_ind < len(currentpoles) and p_ind < len(poles): + if abs(currentpoles[cp_ind] - poles[p_ind]) < (10 * + finfo(float).eps): + # If the current element of both lists match, then we're + # good. Move to the next pair of elements. + cp_ind += 1 + elif currentpoles[cp_ind] < poles[p_ind]: + # We found a pole in this transfer function that's not + # in the list of cumulative poles. Add it to the list. + poles.insert(p_ind, currentpoles[cp_ind]) + # Now mark this pole as "missing" in all previous + # denominators. + for k in range(i): + for m in range(self.inputs): + # All previous rows. + missingpoles[k][m].append(currentpoles[cp_ind]) + for m in range(j): + # This row only. + missingpoles[i][m].append(currentpoles[cp_ind]) + cp_ind += 1 + else: + # There is a pole in the cumulative list of poles that + # is not in our transfer function denominator. Mark + # this pole as "missing", and do not increment cp_ind. + missingpoles[i][j].append(poles[p_ind]) + p_ind += 1 + + if cp_ind == len(currentpoles) and p_ind < len(poles): + # If we finished scanning currentpoles first, then all the + # remaining cumulative poles are missing poles. + missingpoles[i][j].extend(poles[p_ind:]) + elif cp_ind < len(currentpoles) and p_ind == len(poles): + # If we finished scanning the cumulative poles first, then + # all the reamining currentpoles need to be added to poles. + poles.extend(currentpoles[cp_ind:]) + # Now mark these poles as "missing" in previous + # denominators. + for k in range(i): + for m in range(self.inputs): + # All previous rows. + missingpoles[k][m].extend(currentpoles[cp_ind:]) + for m in range(j): + # This row only. + missingpoles[i][m].extend(currentpoles[cp_ind:]) + + # Construct the common denominator. + den = 1. + for p in poles: + den = polymul(den, [1., -p]) + + # Modify the numerators so that they each take the common denominator. + num = deepcopy(self.num) + for i in range(self.outputs): + for j in range(self.inputs): + # The common denominator has leading coefficient 1. Scale out + # the existing denominator's leading coefficient. + assert self.den[i][j][0], "The i = %i, j = %i denominator has \ +a zero leading coefficient." % (i, j) + num[i][j] = num[i][j] / self.den[i][j][0] + # Multiply in the missing poles. + for p in missingpoles[i][j]: + num[i][j] = polymul(num[i][j], [1., -p]) + # Find the largest numerator polynomial size. + largest = 0 + for i in range(self.outputs): + for j in range(self.inputs): + largest = max(largest, len(num[i][j])) + # Pad all smaller numerator polynomials with zeros. + for i in range(self.outputs): + for j in range(self.inputs): + num[i][j] = insert(num[i][j], zeros(largest - len(num[i][j])), + zeros(largest - len(num[i][j]))) + # Finally, convert the numerator to a 3-D array. + num = array(num) + + return num, den + # Utility function to convert a transfer function polynomial to a string # Borrowed from poly1d library def _tfpolyToString(coeffs, var='s'): This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:37
|
Revision: 83 http://python-control.svn.sourceforge.net/python-control/?rev=83&view=rev Author: kkchen Date: 2011-02-08 22:16:31 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Begun work on tf2ss wrapping. Wrote StateSpace._remove_useless_states. TransferFunction.pole and convertToStateSpace have been edited to use TransferFunction._common_den, which returns a common denominator formulation of a MIMO transfer function. _common_den is not completed yet. Also, StateSpace._remove_useless_states was written to remove states marked by entire rows or columns of zeros in the A, B, or C matrix. Now, TestBDAlg.py passes. Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/statesp.py branches/control-0.4a/src/xferfcn.py Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:16:26 UTC (rev 82) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:16:31 UTC (rev 83) @@ -68,12 +68,13 @@ """ -from numpy import angle, any, array, concatenate, cos, dot, empty, exp, eye, \ - ones, pi, poly, poly1d, matrix, roots, sin, zeros +from numpy import all, angle, any, array, concatenate, cos, delete, dot, \ + empty, exp, eye, matrix, ones, pi, poly, poly1d, roots, sin, zeros from numpy.random import rand, randn from numpy.linalg import inv, det, solve from numpy.linalg.linalg import LinAlgError from scipy.signal import lti +from slycot import td04ad from lti import Lti import xferfcn @@ -121,6 +122,46 @@ if self.outputs != D.shape[0]: raise ValueError("D must have the same row size as C.") + # Check for states that don't do anything, and remove them. + self._remove_useless_states() + + def _remove_useless_states(self): + """Check for states that don't do anything, and remove them. + + Scan the A, B, and C matrices for rows or columns of zeros. If the + zeros are such that a particular state has no effect on the input-output + dynamics, then remove that state from the A, B, and C matrices. + + """ + + # Indices of useless states. + useless = [] + + # Search for useless states. + for i in range(self.states): + if (all(self.A[i, :] == zeros((1, self.states))) and + all(self.B[i, :] == zeros((1, self.inputs)))): + useless.append(i) + # To avoid duplucate indices in useless, jump to the next + # iteration. + continue + if (all(self.A[:, i] == zeros((self.states, 1))) and + all(self.C[:, i] == zeros((self.outputs, 1)))): + useless.append(i) + + # Remove the useless states. + if all(useless == range(self.states)): + # All the states were useless. + self.A = 0 + 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) + def __str__(self): """String representation of the state space.""" @@ -374,15 +415,23 @@ # Already a state space system; just return it return sys elif isinstance(sys, xferfcn.TransferFunction): - # TODO: Wrap SLICOT to do transfer function to state space conversion. - raise NotImplementedError("Transfer function to state space conversion \ -is not implemented yet.") - elif (isinstance(sys, (int, long, float, complex))): + # Change the numerator and denominator arrays so that the transfer + # function matrix has a common denominator. + num, den = sys._common_den() + # Make a list of the orders of the denominator polynomials. + index = [len(den) for i in range(sys.outputs)] + # Repeat the common denominator along the rows. + den = array([den for i in range(sys.outputs)]) + + ssout = td04ad(sys.inputs, sys.outputs, index, num, den) + + return StateSpace(ssout[1], ssout[2], ssout[3], ssout[4]) + elif isinstance(sys, (int, long, float, complex)): # Generate a simple state space system of the desired dimension # The following Doesn't work due to inconsistencies in ltisys: # return StateSpace([[]], [[]], [[]], eye(outputs, inputs)) return StateSpace(0., zeros((1, inputs)), zeros((outputs, 1)), - sys * ones(outputs, inputs)) + sys * ones((outputs, inputs))) else: raise TypeError("Can't convert given type to StateSpace system.") Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:26 UTC (rev 82) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:31 UTC (rev 83) @@ -26,6 +26,7 @@ TransferFunction.zero TransferFunction.feedback TransferFunction.returnScipySignalLti +TransferFunction._common_den _tfpolyToString _addSISO convertToTransferFunction @@ -72,7 +73,7 @@ # External function declarations from numpy import angle, array, empty, ndarray, ones, polyadd, polymul, \ - polyval, zeros + polyval, roots, zeros from scipy.signal import lti from copy import deepcopy from slycot import tb04ad @@ -406,9 +407,9 @@ def pole(self): """Compute the poles of a transfer function.""" - raise NotImplementedError("TransferFunction.pole is not implemented \ -yet.") - + num, den = self._common_den() + return roots(den) + def zero(self): """Compute the zeros of a transfer function.""" @@ -463,6 +464,39 @@ return out + def _common_den(self): + """Compute MIMO common denominator; return it and an adjusted numerator. + + >>> n, d = sys._common_den() + + computes the single denominator containing all the poles of sys.den, and + reports it as the array d. + + The output numerator array n is modified to use the common denominator. + It is an sys.outputs-by-sys.inputs-by-[something] array. + + """ + + # Preallocate some variables. Start by figuring out the maximum number + # of numerator coefficients. + numcoeffs = 0 + for i in range(self.outputs): + for j in range(self.inputs): + numcoeffs = max(numcoeffs, len(self.num[i][j])) + # The output 3-D adjusted numerator array. + num = empty((i, j, numcoeffs)) + # A list to keep track of roots found as we scan self.den. + poles = [] + # A 3-D list to keep track of common denominator roots not present in + # the self.den[i][j]. + missingpoles = [[[] for j in range(self.inputs)] for i in + range(self.outputs)] + + for i in range(sys.outputs): + for j in range(sys.inputs): + currentpoles = roots(self.den[i][j]) + #TODO: finish this + # Utility function to convert a transfer function polynomial to a string # Borrowed from poly1d library def _tfpolyToString(coeffs, var='s'): This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:32
|
Revision: 82 http://python-control.svn.sourceforge.net/python-control/?rev=82&view=rev Author: kkchen Date: 2011-02-08 22:16:26 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Implemented and tested convertToTransferFunction (which ss2tf uses). Various fixes: - changed control.exception to exception in the imports of statefbk.py - changed sys + (const) to behave more like MATLAB Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/TestXferFcn.py branches/control-0.4a/src/statefbk.py branches/control-0.4a/src/statesp.py branches/control-0.4a/src/xferfcn.py Modified: branches/control-0.4a/src/TestXferFcn.py =================================================================== --- branches/control-0.4a/src/TestXferFcn.py 2011-02-08 22:16:21 UTC (rev 81) +++ branches/control-0.4a/src/TestXferFcn.py 2011-02-08 22:16:26 UTC (rev 82) @@ -1,7 +1,8 @@ #!/usr/bin/env python import numpy as np -from xferfcn import TransferFunction +from statesp import StateSpace +from xferfcn import TransferFunction, convertToTransferFunction import unittest class TestXferFcn(unittest.TestCase): @@ -383,7 +384,8 @@ # Tests for TransferFunction.feedback. def testFeedbackSISO(self): - + """Test for correct SISO transfer function feedback.""" + sys1 = TransferFunction([-1., 4.], [1., 3., 5.]) sys2 = TransferFunction([2., 3., 0.], [1., -3., 4., 0]) @@ -395,5 +397,27 @@ np.testing.assert_array_equal(sys4.num, [[[-1., 7., -16., 16., 0.]]]) np.testing.assert_array_equal(sys4.den, [[[1., 0., 2., -8., 8., 0.]]]) + def testConvertToTransferFunction(self): + """Test for correct state space to transfer function conversion.""" + + A = [[1., -2.], [-3., 4.]] + B = [[6., 5.], [4., 3.]] + C = [[1., -2.], [3., -4.], [5., -6.]] + D = [[1., 0.], [0., 1.], [1., 0.]] + sys = StateSpace(A, B, C, D) + + tfsys = convertToTransferFunction(sys) + + num = [[np.array([1., -7., 10.]), np.array([-1., 10.])], + [np.array([2., -8.]), np.array([1., -2., -8.])], + [np.array([1., 1., -30.]), np.array([7., -22.])]] + den = [[np.array([1., -5., -2.]) for j in range(sys.inputs)] + for i in range(sys.outputs)] + + for i in range(sys.outputs): + for j in range(sys.inputs): + np.testing.assert_array_almost_equal(tfsys.num[i][j], num[i][j]) + np.testing.assert_array_almost_equal(tfsys.den[i][j], den[i][j]) + if __name__ == "__main__": unittest.main() Modified: branches/control-0.4a/src/statefbk.py =================================================================== --- branches/control-0.4a/src/statefbk.py 2011-02-08 22:16:21 UTC (rev 81) +++ branches/control-0.4a/src/statefbk.py 2011-02-08 22:16:26 UTC (rev 82) @@ -42,7 +42,7 @@ # External packages and modules import numpy as np import ctrlutil -from control.exception import * +from exception import * # Pole placement def place(A, B, p): Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:16:21 UTC (rev 81) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:16:26 UTC (rev 82) @@ -69,7 +69,7 @@ """ from numpy import angle, any, array, concatenate, cos, dot, empty, exp, eye, \ - pi, poly, poly1d, matrix, roots, sin, zeros + ones, pi, poly, poly1d, matrix, roots, sin, zeros from numpy.random import rand, randn from numpy.linalg import inv, det, solve from numpy.linalg.linalg import LinAlgError @@ -382,7 +382,7 @@ # The following Doesn't work due to inconsistencies in ltisys: # return StateSpace([[]], [[]], [[]], eye(outputs, inputs)) return StateSpace(0., zeros((1, inputs)), zeros((outputs, 1)), - sys * eye(outputs, inputs)) + sys * ones(outputs, inputs)) else: raise TypeError("Can't convert given type to StateSpace system.") Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:21 UTC (rev 81) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:26 UTC (rev 82) @@ -75,6 +75,7 @@ polyval, zeros from scipy.signal import lti from copy import deepcopy +from slycot import tb04ad from lti import Lti import statesp @@ -397,8 +398,8 @@ polyval(self.den[i][j], w * 1.j)), omega) fresp = array(fresp) - mag[i, j] = abs(fresp) - phase[i, j] = angle(fresp) + mag[i, j, :] = abs(fresp) + phase[i, j, :] = angle(fresp) return mag, phase, omega @@ -537,15 +538,24 @@ if isinstance(sys, TransferFunction): return sys elif isinstance(sys, statesp.StateSpace): - # TODO: Wrap SLICOT to do state space to transfer function conversion. - raise NotImplementedError("State space to transfer function conversion \ -is not implemented yet.") + # Use Slycot to make the transformation. + tfout = tb04ad(sys.states, sys.inputs, sys.outputs, sys.A, sys.B, sys.C, + sys.D, sys.outputs, sys.outputs, sys.inputs) + + # Preallocate outputs. + num = [[[] for j in range(sys.inputs)] for i in range(sys.outputs)] + den = [[[] for j in range(sys.inputs)] for i in range(sys.outputs)] + + for i in range(sys.outputs): + for j in range(sys.inputs): + num[i][j] = list(tfout[6][i, j, :]) + # Each transfer function matrix row has a common denominator. + den[i][j] = list(tfout[5][i, :]) + + return TransferFunction(num, den) elif isinstance(sys, (int, long, float, complex)): - # Make an identity system. - num = [[[0] for j in range(inputs)] for i in range(outputs)] + num = [[[sys] for j in range(inputs)] for i in range(outputs)] den = [[[1] for j in range(inputs)] for i in range(outputs)] - for i in range(min(inputs, outputs)): - num[i][i] = [sys] return TransferFunction(num, den) else: This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:27
|
Revision: 81 http://python-control.svn.sourceforge.net/python-control/?rev=81&view=rev Author: kkchen Date: 2011-02-08 22:16:21 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Added function hinfsyn to robust.py, based on slycot wrapper sb10ad. Steven Brunton <sbr...@pr...> Modified Paths: -------------- branches/control-0.4a/src/robust.py Modified: branches/control-0.4a/src/robust.py =================================================================== --- branches/control-0.4a/src/robust.py 2011-02-08 22:16:17 UTC (rev 80) +++ branches/control-0.4a/src/robust.py 2011-02-08 22:16:21 UTC (rev 81) @@ -47,15 +47,15 @@ from statesp import StateSpace def h2syn(P,nmeas,ncon): - """H_2 control synthesis for plant. + """H_2 control synthesis for plant P. Usage ===== - K, rcond, info = h2syn(plant,nmeas,ncon) + K = h2syn(P,nmeas,ncon) Inputs ====== - plant : partitioned lti plant + P : partitioned lti plant nmeas : number of measurements (input to controller) ncon : number of control inputs (output from controller) @@ -74,12 +74,6 @@ # else: dico = 'C' - #Check system is stable - D,V = np.linalg.eig(P.A) - for e in D: - if e.real >= 0: - raise ValueError, "Oops, the system is unstable!" - try: from slycot import sb10hd except ImportError: @@ -97,3 +91,63 @@ K = StateSpace(Ak, Bk, Ck, Dk) return K + +def hinfsyn(P,nmeas,ncon): + """H_{inf} control synthesis for plant P. + + Usage + ===== + K, CL, gam, info = hinfsyn(P,nmeas,ncon) + + Inputs + ====== + P : partitioned lti plant + nmeas : number of measurements (input to controller) + ncon : number of control inputs (output from controller) + + Outputs + ======= + K : controller to stabilize P + CL : closed loop system + gam : infinity norm of closed loop system + info : info returned from siycot routine + """ + + #Check for ss system object, need a utility for this? + + #TODO: Check for continous or discrete, only continuous supported right now + # if isCont(): + # dico = 'C' + # elif isDisc(): + # dico = 'D' + # else: + dico = 'C' + + try: + from slycot import sb10ad + except ImportError: + raise ControlSlycot("can't find slycot subroutine sb10ad") + + job = 3 + n = np.size(P.A,0) + m = np.size(P.B,1) + np = np.size(P.C,0) + gamma = 1.e100 + out = sb10ad(job,n,m,np,ncon,nmeas,gamma,P.A,P.B,P.C,P.D) + gam = out[0] + Ak = out[1] + Bk = out[2] + Ck = out[3] + Dk = out[4] + Ac = out[5] + Bc = out[6] + Cc = out[7] + Dc = out[8] + rcond = out[9] + info = out[10] + + K = StateSpace(Ak, Bk, Ck, Dk) + CL = StateSpace(Ac, Bc, Cc, Dc) + + return K, CL, gam, info + This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:23
|
Revision: 80 http://python-control.svn.sourceforge.net/python-control/?rev=80&view=rev Author: kkchen Date: 2011-02-08 22:16:17 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Added robust.py for robust control tools. Function h2syn added, based on slycot wrapper sb10hd. No tests in place, probably buggy. Steven Brunton <sbr...@pr...> Added Paths: ----------- branches/control-0.4a/src/robust.py Added: branches/control-0.4a/src/robust.py =================================================================== --- branches/control-0.4a/src/robust.py (rev 0) +++ branches/control-0.4a/src/robust.py 2011-02-08 22:16:17 UTC (rev 80) @@ -0,0 +1,99 @@ +# robust.py - tools for robust control +# +# Author: Steve Brunton, Kevin Chen, Lauren Padilla +# Date: 24 Dec 2010 +# +# This file contains routines for obtaining reduced order models +# +# Copyright (c) 2010 by California Institute of Technology +# All rights reserved. +# +# Redistribution and use in source and binary forms, with or without +# modification, are permitted provided that the following conditions +# are met: +# +# 1. Redistributions of source code must retain the above copyright +# notice, this list of conditions and the following disclaimer. +# +# 2. Redistributions in binary form must reproduce the above copyright +# notice, this list of conditions and the following disclaimer in the +# documentation and/or other materials provided with the distribution. +# +# 3. Neither the name of the California Institute of Technology nor +# the names of its contributors may be used to endorse or promote +# products derived from this software without specific prior +# written permission. +# +# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +# "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT +# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS +# FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL CALTECH +# OR THE CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, +# SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT +# LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF +# USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND +# ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, +# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT +# OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF +# SUCH DAMAGE. +# +# $Id$ + +# External packages and modules +import numpy as np +import ctrlutil +from control.exception import * +from statefbk import * +from statesp import StateSpace + +def h2syn(P,nmeas,ncon): + """H_2 control synthesis for plant. + + Usage + ===== + K, rcond, info = h2syn(plant,nmeas,ncon) + + Inputs + ====== + plant : partitioned lti plant + nmeas : number of measurements (input to controller) + ncon : number of control inputs (output from controller) + + Outputs + ======= + K : controller to stabilize P + """ + + #Check for ss system object, need a utility for this? + + #TODO: Check for continous or discrete, only continuous supported right now + # if isCont(): + # dico = 'C' + # elif isDisc(): + # dico = 'D' + # else: + dico = 'C' + + #Check system is stable + D,V = np.linalg.eig(P.A) + for e in D: + if e.real >= 0: + raise ValueError, "Oops, the system is unstable!" + + try: + from slycot import sb10hd + except ImportError: + raise ControlSlycot("can't find slycot subroutine sb10hd") + + n = np.size(P.A,0) + m = np.size(P.B,1) + np = np.size(P.C,0) + out = sb10hd(n,m,np,ncon,nmeas,P.A,P.B,P.C,P.D) + Ak = out[0] + Bk = out[1] + Ck = out[2] + Dk = out[3] + + K = StateSpace(Ak, Bk, Ck, Dk) + + return K This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:20
|
Revision: 79 http://python-control.svn.sourceforge.net/python-control/?rev=79&view=rev Author: kkchen Date: 2011-02-08 22:16:12 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Modified matlab functions impulse and step for keywords. Tests pass. Steven Brunton <sbr...@pr...> Modified Paths: -------------- branches/control-0.4a/src/TestMatlab.py branches/control-0.4a/src/matlab.py Modified: branches/control-0.4a/src/TestMatlab.py =================================================================== --- branches/control-0.4a/src/TestMatlab.py 2011-02-08 22:16:06 UTC (rev 78) +++ branches/control-0.4a/src/TestMatlab.py 2011-02-08 22:16:12 UTC (rev 79) @@ -12,7 +12,7 @@ D = np.matrix("9.") sys = ss(A,B,C,D) t = np.linspace(0, 1, 10) - t, yout = step(sys, t) + t, yout = step(sys, T=t) youttrue = np.matrix("9. 17.6457 24.7072 30.4855 35.2234 39.1165 42.3227 44.9694 47.1599 48.9776") np.testing.assert_array_almost_equal(yout, youttrue,decimal=4) @@ -23,7 +23,7 @@ D = np.matrix("9.") sys = ss(A,B,C,D) t = np.linspace(0, 1, 10) - t, yout = impulse(sys, t) + t, yout = impulse(sys, T=t) youttrue = np.matrix("86. 70.1808 57.3753 46.9975 38.5766 31.7344 26.1668 21.6292 17.9245 14.8945") np.testing.assert_array_almost_equal(yout, youttrue,decimal=4) @@ -35,7 +35,7 @@ sys = ss(A,B,C,D) t = np.linspace(0, 1, 10) x0 = np.matrix(".5; 1.") - t, yout = initial(sys, t, x0) + t, yout = initial(sys, T=t, X0=x0) youttrue = np.matrix("11. 8.1494 5.9361 4.2258 2.9118 1.9092 1.1508 0.5833 0.1645 -0.1391") np.testing.assert_array_almost_equal(yout, youttrue,decimal=4) Modified: branches/control-0.4a/src/matlab.py =================================================================== --- branches/control-0.4a/src/matlab.py 2011-02-08 22:16:06 UTC (rev 78) +++ branches/control-0.4a/src/matlab.py 2011-02-08 22:16:12 UTC (rev 79) @@ -791,7 +791,12 @@ ltiobjs = sys.returnScipySignalLti() ltiobj = ltiobjs[0][0] - return sp.signal.step(ltiobj, **keywords) + out = sp.signal.step(ltiobj, **keywords) + yout = [] + yout.append(np.mat(out[0])) + yout.append(out[1]) + yout = tuple(yout) + return yout # Redefine initial to use lsim2 #! Not yet implemented (uses step for now) @@ -815,7 +820,8 @@ ltiobjs = sys.returnScipySignalLti() ltiobj = ltiobjs[0][0] - return sp.signal.initial(ltiobj, **keywords) + yout = sp.signal.initial(ltiobj, **keywords) + return np.mat(yout) # Redefine impulse to use initial() #! Not yet implemented (uses impulse for now) @@ -839,5 +845,5 @@ ltiobjs = sys.returnScipySignalLti() ltiobj = ltiobjs[0][0] - return sp.signal.impulse(ltiobj, **keywords) - + yout = sp.signal.impulse(ltiobj, **keywords) + return np.mat(yout) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:13
|
Revision: 78 http://python-control.svn.sourceforge.net/python-control/?rev=78&view=rev Author: kkchen Date: 2011-02-08 22:16:06 +0000 (Tue, 08 Feb 2011) Log Message: ----------- More docstrings written in preparation for Sphinx documentation. Miscellaneous edits, including raising NotImplementedError, using parentheses instead of '\' for line continuation, and converting types. Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/bdalg.py branches/control-0.4a/src/matlab.py branches/control-0.4a/src/statesp.py branches/control-0.4a/src/xferfcn.py Modified: branches/control-0.4a/src/bdalg.py =================================================================== --- branches/control-0.4a/src/bdalg.py 2011-02-08 22:16:00 UTC (rev 77) +++ branches/control-0.4a/src/bdalg.py 2011-02-08 22:16:06 UTC (rev 78) @@ -68,7 +68,7 @@ Raises ------ ValueError - if sys2.inputs does not equal sys1.outputs + if `sys2.inputs` does not equal `sys1.outputs` See Also -------- @@ -78,8 +78,8 @@ Notes ----- This function is a wrapper for the __mul__ function in the StateSpace and - TransferFunction classes. The output type is usually the type of sys2. If - sys2 is a scalar, then the output type is the type of sys1. + TransferFunction classes. The output type is usually the type of `sys2`. + If `sys2` is a scalar, then the output type is the type of `sys1`. Examples -------- @@ -105,7 +105,7 @@ Raises ------ ValueError - if sys1 and sys2 do not have the same numbers of inputs and outputs + if `sys1` and `sys2` do not have the same numbers of inputs and outputs See Also -------- @@ -115,8 +115,8 @@ Notes ----- This function is a wrapper for the __add__ function in the StateSpace and - TransferFunction classes. The output type is usually the type of sys1. If - sys1 is a scalar, then the output type is the type of sys2. + TransferFunction classes. The output type is usually the type of `sys1`. If + `sys1` is a scalar, then the output type is the type of `sys2`. Examples ------- @@ -162,9 +162,9 @@ sys2: scalar, StateSpace, or TransferFunction The feedback plant (often a feedback controller). sign: scalar - The sign of feedback. sign = -1 indicates negative feedback, and sign = 1 - indicates positive feedback. sign is an optional argument; it assumes a - value of -1 if not specified. + The sign of feedback. `sign` = -1 indicates negative feedback, and + `sign` = 1 indicates positive feedback. `sign` is an optional argument; it + assumes a value of -1 if not specified. Returns ------- @@ -173,8 +173,8 @@ Raises ------ ValueError - if sys1 does not have as many inputs as sys2 has outputs, or if sys2 - does not have as many inputs as sys1 has outputs + if `sys1` does not have as many inputs as `sys2` has outputs, or if + `sys2` does not have as many inputs as `sys1` has outputs NotImplementedError if an attempt is made to perform a feedback on a MIMO TransferFunction object @@ -187,11 +187,11 @@ Notes ---- This function is a wrapper for the feedback function in the StateSpace and - TransferFunction classes. It calls TransferFunction.feedback if sys1 is a - TransferFunction object, and StateSpace.feedback if sys1 is a StateSpace - object. If sys1 is a scalar, then it is converted to sys2's type, and the - corresponding feedback function is used. If sys1 and sys2 are both scalars, - then TransferFunction.feedback is used. + TransferFunction classes. It calls TransferFunction.feedback if `sys1` is a + TransferFunction object, and StateSpace.feedback if `sys1` is a StateSpace + object. If `sys1` is a scalar, then it is converted to `sys2`'s type, and + the corresponding feedback function is used. If `sys1` and `sys2` are both + scalars, then TransferFunction.feedback is used. """ Modified: branches/control-0.4a/src/matlab.py =================================================================== --- branches/control-0.4a/src/matlab.py 2011-02-08 22:16:00 UTC (rev 77) +++ branches/control-0.4a/src/matlab.py 2011-02-08 22:16:06 UTC (rev 78) @@ -274,7 +274,7 @@ C: numpy matrix or matrix-like object D: numpy matrix or matrix-like object sys: StateSpace or TransferFunction object - ss accepts a set of A, B, C, D matrices or sys. + ss accepts a set of `A`, `B`, `C`, `D` matrices or `sys`. Returns ------- @@ -321,7 +321,7 @@ num: vector, or list of lists of vectors den: vector, or list of lists of vectors sys: StateSpace or TransferFunction object - tf accepts a num and den, or sys. + tf accepts a `num` and `den`, or `sys``. Returns ------- @@ -330,9 +330,9 @@ Raises ------ ValueError - if num and den have invalid or unequal dimensions + if `num` and `den` have invalid or unequal dimensions TypeError - if num or den are of incorrect type + if `num` or `den` are of incorrect type See Also -------- @@ -342,9 +342,9 @@ Notes -------- - num[i][j] is the vector of polynomial coefficients of the transfer function - numerator from the (j+1)st output to the (i+1)st input. den[i][j] works the - same way. + `num`[`i`][`j`] is the vector of polynomial coefficients of the transfer + function numerator from the (`j`+1)st output to the (`i`+1)st input. + `den`[`i`][`j`] works the same way. Examples -------- @@ -382,7 +382,7 @@ C: numpy matrix or matrix-like object D: numpy matrix or matrix-like object sys: StateSpace object - ss accepts a set of A, B, C, D matrices or a StateSpace object. + ss accepts a set of `A`, `B`, `C`, `D` matrices, or `sys`. Returns ------- @@ -394,7 +394,7 @@ if matrix sizes are not self-consistent, or if an invalid number of arguments is passed in TypeError - if sys is not a StateSpace object + if `sys` is not a StateSpace object See Also -------- @@ -432,7 +432,7 @@ num: vector, or list of lists of vectors den: vector, or list of lists of vectors sys: TransferFunction object - tf2ss accepts num and den, or sys. + tf2ss accepts `num` and `den`, or `sys`. Returns ------- @@ -441,11 +441,11 @@ Raises ------ ValueError - if num and den have invalid or unequal dimensions, or if an invalid + if `num` and `den` have invalid or unequal dimensions, or if an invalid number of arguments is passed in TypeError - if num or den are of incorrect type, or if sys is not a TransferFunction - object + if `num` or `den` are of incorrect type, or if sys is not a + TransferFunction object See Also -------- @@ -455,9 +455,9 @@ Notes -------- - num[i][j] is the vector of polynomial coefficients of the transfer function - numerator from the (j+1)st output to the (i+1)st input. den[i][j] works the - same way. + `num`[`i`][`j`] is the vector of polynomial coefficients of the transfer + function numerator from the (`j`+1)st output to the (`i`+1)st input. + `den`[`i`][`j`] works the same way. Examples -------- @@ -481,35 +481,203 @@ raise ValueError("Needs 1 or 2 arguments; received %i." % len(args)) def rss(states=1, inputs=1, outputs=1): - """Create a stable continuous random state space object.""" + """ + Create a stable continuous random state space object. + Parameters + ---------- + states: integer + inputs: integer + outputs: integer + + Returns + ------- + sys: StateSpace object + + Raises + ------ + ValueError + if any input is not a positive integer + + See Also + -------- + drss + + Notes + ----- + If the number of states, inputs, or outputs is not specified, then the + missing numbers are assumed to be 1. The poles of the returned system will + always have a negative real part. + + """ + return rss_generate(states, inputs, outputs, 'c') def drss(states=1, inputs=1, outputs=1): - """Create a stable discrete random state space object.""" + """ + Create a stable discrete random state space object. + Parameters + ---------- + states: integer + inputs: integer + outputs: integer + + Returns + ------- + sys: StateSpace object + + Raises + ------ + ValueError + if any input is not a positive integer + + See Also + -------- + rss + + Notes + ----- + If the number of states, inputs, or outputs is not specified, then the + missing numbers are assumed to be 1. The poles of the returned system will + always have a magnitude less than 1. + + """ + return rss_generate(states, inputs, outputs, 'd') def pole(sys): - """Return system poles.""" + """ + Return system poles. + Parameters + ---------- + sys: StateSpace or TransferFunction object + + Returns + ------- + poles: ndarray + + Raises + ------ + NotImplementedError + when called on a TransferFunction object + + See Also + -------- + zero + + Notes + ----- + This function is a wrapper for StateSpace.pole and TransferFunction.pole. + + """ + return sys.pole() def zero(sys): - """Return system zeros.""" + """ + Return system zeros. + Parameters + ---------- + sys: StateSpace or TransferFunction object + + Returns + ------- + zeros: ndarray + + Raises + ------ + NotImplementedError + when called on a TransferFunction object or a MIMO StateSpace object + + See Also + -------- + pole + + Notes + ----- + This function is a wrapper for StateSpace.zero and TransferFunction.zero. + + """ + return sys.zero() def evalfr(sys, omega): - """Evaluate the transfer function of an LTI system at a single frequency - omega.""" + """ + Evaluate the transfer function of an LTI system at an angular frequency. + Parameters + ---------- + sys: StateSpace or TransferFunction object + omega: scalar + + Returns + ------- + fresp: ndarray + + See Also + -------- + freqresp + bode + + Notes + ----- + This function is a wrapper for StateSpace.evalfr and + TransferFunction.evalfr. + + Examples + -------- + >>> sys = rss(3, 2, 2) + >>> evalfr(sys, 1.) + array([[ 4.09376126 -6.2171555j , 23.71332080-35.24245284j], + [ 0.83405186 -1.82896006j, 8.10962251-12.66640309j]]) + This is the transfer function matrix evaluated at s = i. + + """ + return sys.evalfr(omega) def freqresp(sys, omega): - """Return the frequency response for an LTI object at a list of frequencies - omega.""" + """ + Frequency response of an LTI system at multiple angular frequencies. + Parameters + ---------- + sys: StateSpace or TransferFunction object + omega: list, tuple, or ndarray + + Returns + ------- + mag: ndarray + phase: ndarray + omega: list, tuple, or ndarray + + See Also + -------- + evalfr + bode + + Notes + ----- + This function is a wrapper for StateSpace.freqresp and + TransferFunction.freqresp. + + Examples + -------- + >>> sys = rss(3, 2, 2) + >>> mag, phase, omega = freqresp(sys, [0.1, 1., 10.]) + >>> mag[0, 1, :] + array([ 55.43747231, 42.47766549, 1.97225895]) + >>> phase[1, 0, :] + array([-0.12611087, -1.14294316, 2.5764547 ]) + This is the magnitude of the frequency response from the 2nd input to the + 1st output, and the phase (in radians) of the frequency response from the + 1st input to the 2nd output, for s = 0.1i, i, 10i. + + """ + return sys.freqresp(omega) # Bode plots Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:16:00 UTC (rev 77) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:16:06 UTC (rev 78) @@ -68,8 +68,8 @@ """ -from numpy import angle, any, concatenate, cos, dot, empty, exp, eye, pi, \ - poly, poly1d, matrix, roots, sin, zeros +from numpy import angle, any, array, concatenate, cos, dot, empty, exp, eye, \ + pi, poly, poly1d, matrix, roots, sin, zeros from numpy.random import rand, randn from numpy.linalg import inv, det, solve from numpy.linalg.linalg import LinAlgError @@ -250,10 +250,10 @@ fresp = self.C * solve(omega * 1.j * eye(self.states) - self.A, self.B) + self.D - return fresp + return array(fresp) # Method for generating the frequency response of the system - def freqresp(self, omega=None): + def freqresp(self, omega): """Evaluate the system's transfer func. at a list of ang. frequencies. mag, phase, omega = self.freqresp(omega) @@ -294,10 +294,10 @@ den = poly1d(poly(self.A)) # Compute the numerator based on zeros #! TODO: This is currently limited to SISO systems - num = poly1d(\ - poly(self.A - dot(self.B, self.C)) + (self.D[0, 0] - 1) * den) + num = poly1d(poly(self.A - dot(self.B, self.C)) + ((self.D[0, 0] - 1) * + den)) - return (roots(num)) + return roots(num) # Feedback around a state space system def feedback(self, other, sign=-1): @@ -407,6 +407,17 @@ # Probability that D = 0. pDzero = 0.5 + # Check for valid input arguments. + if states < 1 or states % 1: + raise ValueError(("states must be a positive integer. states = %g." % + states)) + if inputs < 1 or inputs % 1: + raise ValueError(("inputs must be a positive integer. inputs = %g." % + inputs)) + if outputs < 1 or outputs % 1: + raise ValueError(("outputs must be a positive integer. outputs = %g." % + outputs)) + # Make some poles for A. Preallocate a complex array. poles = zeros(states) + zeros(states) * 0.j i = 0 Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:00 UTC (rev 77) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:16:06 UTC (rev 78) @@ -106,14 +106,14 @@ if isinstance(data[i], (int, float, long, complex)): # Convert scalar to list of list of array. data[i] = [[array([data[i]])]] - elif isinstance(data[i], (list, tuple, ndarray)) and \ - isinstance(data[i][0], (int, float, long, complex)): + elif (isinstance(data[i], (list, tuple, ndarray)) and + isinstance(data[i][0], (int, float, long, complex))): # Convert array to list of list of array. data[i] = [[array(data[i])]] - elif isinstance(data[i], list) and \ - isinstance(data[i][0], list) and \ - isinstance(data[i][0][0], (list, tuple, ndarray)) and \ - isinstance(data[i][0][0][0], (int, float, long, complex)): + elif (isinstance(data[i], list) and + isinstance(data[i][0], list) and + isinstance(data[i][0][0], (list, tuple, ndarray)) and + isinstance(data[i][0][0][0], (int, float, long, complex))): # We might already have the right format. Convert the # coefficient vectors to arrays, if necessary. for j in range(len(data[i])): @@ -226,14 +226,12 @@ dashes = '-' * dashcount # Center the numerator or denominator - if (len(numstr) < dashcount): - numstr = ' ' * \ - int(round((dashcount - len(numstr))/2)) + \ - numstr - if (len(denstr) < dashcount): - denstr = ' ' * \ - int(round((dashcount - len(denstr))/2)) + \ - denstr + if len(numstr) < dashcount: + numstr = (' ' * int(round((dashcount - len(numstr))/2)) + + numstr) + if len(denstr) < dashcount: + denstr = (' ' * int(round((dashcount - len(denstr))/2)) + + denstr) outstr += "\n" + numstr + "\n" + dashes + "\n" + denstr + "\n" return outstr @@ -334,8 +332,8 @@ def __div__(self, other): """Divide two LTI objects.""" - if self.inputs > 1 or self.outputs > 1 or \ - other.inputs > 1 or other.outputs > 1: + if (self.inputs > 1 or self.outputs > 1 or + other.inputs > 1 or other.outputs > 1): raise NotImplementedError("TransferFunction.__div__ is currently \ implemented only for SISO systems.") @@ -351,8 +349,8 @@ def __rdiv__(self, other): """Reverse divide two LTI objects.""" - if self.inputs > 1 or self.outputs > 1 or \ - other.inputs > 1 or other.outputs > 1: + if (self.inputs > 1 or self.outputs > 1 or + other.inputs > 1 or other.outputs > 1): raise NotImplementedError("TransferFunction.__rdiv__ is currently \ implemented only for SISO systems.") @@ -371,13 +369,13 @@ for i in range(self.outputs): for j in range(self.inputs): - out[i][j] = polyval(self.num[i][j], omega * 1.j) / \ - polyval(self.den[i][j], omega * 1.j) + out[i][j] = (polyval(self.num[i][j], omega * 1.j) / + polyval(self.den[i][j], omega * 1.j)) return out # Method for generating the frequency response of the system - def freqresp(self, omega=None): + def freqresp(self, omega): """Evaluate a transfer function at a list of angular frequencies. mag, phase, omega = self.freqresp(omega) @@ -395,8 +393,8 @@ for i in range(self.outputs): for j in range(self.inputs): - fresp = map(lambda w: polyval(self.num[i][j], w * 1.j) / \ - polyval(self.den[i][j], w * 1.j), omega) + fresp = map(lambda w: (polyval(self.num[i][j], w * 1.j) / + polyval(self.den[i][j], w * 1.j)), omega) fresp = array(fresp) mag[i, j] = abs(fresp) @@ -407,20 +405,22 @@ def pole(self): """Compute the poles of a transfer function.""" - pass + raise NotImplementedError("TransferFunction.pole is not implemented \ +yet.") def zero(self): """Compute the zeros of a transfer function.""" - pass + raise NotImplementedError("TransferFunction.zero is not implemented \ +yet.") def feedback(self, other, sign=-1): """Feedback interconnection between two LTI objects.""" other = convertToTransferFunction(other) - if self.inputs > 1 or self.outputs > 1 or \ - other.inputs > 1 or other.outputs > 1: + if (self.inputs > 1 or self.outputs > 1 or + other.inputs > 1 or other.outputs > 1): # TODO: MIMO feedback raise NotImplementedError("TransferFunction.feedback is currently \ only implemented for SISO functions.") This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:06
|
Revision: 77 http://python-control.svn.sourceforge.net/python-control/?rev=77&view=rev Author: kkchen Date: 2011-02-08 22:16:00 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Added TestMatlab to test impulse, step and initial. Tests fail. There is no scipy.signal.lti initial() function, so we will have to use lsim. Steven Brunton <sbr...@pr...> Modified Paths: -------------- branches/control-0.4a/src/matlab.py Added Paths: ----------- branches/control-0.4a/src/TestMatlab.py Added: branches/control-0.4a/src/TestMatlab.py =================================================================== --- branches/control-0.4a/src/TestMatlab.py (rev 0) +++ branches/control-0.4a/src/TestMatlab.py 2011-02-08 22:16:00 UTC (rev 77) @@ -0,0 +1,44 @@ +#!/usr/bin/env python + +from matlab import * +import numpy as np +import unittest + +class TestMatlab(unittest.TestCase): + def testStep(self): + A = np.matrix("1. -2.; 3. -4.") + B = np.matrix("5.; 7.") + C = np.matrix("6. 8.") + D = np.matrix("9.") + sys = ss(A,B,C,D) + t = np.linspace(0, 1, 10) + t, yout = step(sys, t) + youttrue = np.matrix("9. 17.6457 24.7072 30.4855 35.2234 39.1165 42.3227 44.9694 47.1599 48.9776") + np.testing.assert_array_almost_equal(yout, youttrue,decimal=4) + + def testImpulse(self): + A = np.matrix("1. -2.; 3. -4.") + B = np.matrix("5.; 7.") + C = np.matrix("6. 8.") + D = np.matrix("9.") + sys = ss(A,B,C,D) + t = np.linspace(0, 1, 10) + t, yout = impulse(sys, t) + youttrue = np.matrix("86. 70.1808 57.3753 46.9975 38.5766 31.7344 26.1668 21.6292 17.9245 14.8945") + np.testing.assert_array_almost_equal(yout, youttrue,decimal=4) + + def testInitial(self): + A = np.matrix("1. -2.; 3. -4.") + B = np.matrix("5.; 7.") + C = np.matrix("6. 8.") + D = np.matrix("9.") + sys = ss(A,B,C,D) + t = np.linspace(0, 1, 10) + x0 = np.matrix(".5; 1.") + t, yout = initial(sys, t, x0) + youttrue = np.matrix("11. 8.1494 5.9361 4.2258 2.9118 1.9092 1.1508 0.5833 0.1645 -0.1391") + np.testing.assert_array_almost_equal(yout, youttrue,decimal=4) + + +if __name__ == '__main__': + unittest.main() Property changes on: branches/control-0.4a/src/TestMatlab.py ___________________________________________________________________ Added: svn:executable + * Modified: branches/control-0.4a/src/matlab.py =================================================================== --- branches/control-0.4a/src/matlab.py 2011-02-08 22:15:55 UTC (rev 76) +++ branches/control-0.4a/src/matlab.py 2011-02-08 22:16:00 UTC (rev 77) @@ -599,7 +599,7 @@ ltiobjs = sys.returnScipySignalLti() ltiobj = ltiobjs[0][0] - return sp.signal.lsim2(*args, **keywords) + return sp.signal.lsim2(ltiobj, **keywords) #! Redefine step to use lsim2 #! Not yet implemented @@ -622,15 +622,8 @@ sys = args[0] ltiobjs = sys.returnScipySignalLti() ltiobj = ltiobjs[0][0] - newargs = [] - newargs.append(ltiobj) - for i in range(1, len(args)): - newargs.append(args[i]) - newargs = tuple(newargs) - print len(args) - print len(newargs) - return sp.signal.step(*newargs, **keywords) + return sp.signal.step(ltiobj, **keywords) # Redefine initial to use lsim2 #! Not yet implemented (uses step for now) @@ -654,7 +647,7 @@ ltiobjs = sys.returnScipySignalLti() ltiobj = ltiobjs[0][0] - return sp.signal.initial(*args, **keywords) + return sp.signal.initial(ltiobj, **keywords) # Redefine impulse to use initial() #! Not yet implemented (uses impulse for now) This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:16:01
|
Revision: 76 http://python-control.svn.sourceforge.net/python-control/?rev=76&view=rev Author: kkchen Date: 2011-02-08 22:15:55 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Merge of revisions 72 to 74. Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/matlab.py Modified: branches/control-0.4a/src/matlab.py =================================================================== --- branches/control-0.4a/src/matlab.py 2011-02-08 22:15:51 UTC (rev 75) +++ branches/control-0.4a/src/matlab.py 2011-02-08 22:15:55 UTC (rev 76) @@ -595,6 +595,10 @@ yout response of the system xout time evolution of the state vector """ + sys = args[0] + ltiobjs = sys.returnScipySignalLti() + ltiobj = ltiobjs[0][0] + return sp.signal.lsim2(*args, **keywords) #! Redefine step to use lsim2 @@ -615,8 +619,19 @@ T time values of the output yout response of the system """ - return sp.signal.step(*args, **keywords) + sys = args[0] + ltiobjs = sys.returnScipySignalLti() + ltiobj = ltiobjs[0][0] + newargs = [] + newargs.append(ltiobj) + for i in range(1, len(args)): + newargs.append(args[i]) + newargs = tuple(newargs) + print len(args) + print len(newargs) + return sp.signal.step(*newargs, **keywords) + # Redefine initial to use lsim2 #! Not yet implemented (uses step for now) def initial(*args, **keywords): @@ -635,6 +650,10 @@ T time values of the output yout response of the system """ + sys = args[0] + ltiobjs = sys.returnScipySignalLti() + ltiobj = ltiobjs[0][0] + return sp.signal.initial(*args, **keywords) # Redefine impulse to use initial() @@ -655,5 +674,9 @@ T time values of the output yout response of the system """ - return sp.signal.impulse(*args, **keywords) + sys = args[0] + ltiobjs = sys.returnScipySignalLti() + ltiobj = ltiobjs[0][0] + return sp.signal.impulse(ltiobj, **keywords) + This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:15:58
|
Revision: 75 http://python-control.svn.sourceforge.net/python-control/?rev=75&view=rev Author: kkchen Date: 2011-02-08 22:15:51 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Numerous docstring edits to comply with Sphinx and Python styles. Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/bdalg.py branches/control-0.4a/src/lti.py branches/control-0.4a/src/matlab.py branches/control-0.4a/src/statesp.py branches/control-0.4a/src/xferfcn.py Modified: branches/control-0.4a/src/bdalg.py =================================================================== --- branches/control-0.4a/src/bdalg.py 2011-02-08 22:15:45 UTC (rev 74) +++ branches/control-0.4a/src/bdalg.py 2011-02-08 22:15:51 UTC (rev 75) @@ -1,98 +1,200 @@ -# bdalg.py - functions for implmeenting block diagram algebra -# -# Author: Richard M. Murray -# Date: 24 May 09 -# -# This file contains some standard block diagram algebra. If all -# arguments are SISO transfer functions, the results are a transfer -# function. Otherwise, the computation is done in state space. -# -#! State space operations are not currently implemented. -# -# Copyright (c) 2010 by California Institute of Technology -# All rights reserved. -# -# Redistribution and use in source and binary forms, with or without -# modification, are permitted provided that the following conditions -# are met: -# -# 1. Redistributions of source code must retain the above copyright -# notice, this list of conditions and the following disclaimer. -# -# 2. Redistributions in binary form must reproduce the above copyright -# notice, this list of conditions and the following disclaimer in the -# documentation and/or other materials provided with the distribution. -# -# 3. Neither the name of the California Institute of Technology nor -# the names of its contributors may be used to endorse or promote -# products derived from this software without specific prior -# written permission. -# -# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS -# "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT -# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS -# FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL CALTECH -# OR THE CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, -# SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT -# LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF -# USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND -# ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, -# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT -# OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF -# SUCH DAMAGE. -# -# $Id: bdalg.py 17 2010-05-29 23:50:52Z murrayrm $ +"""bdalg.py +This file contains some standard block diagram algebra. + +Routines in this module: + +series +parallel +negate +feedback + +Copyright (c) 2010 by California Institute of Technology +All rights reserved. + +Redistribution and use in source and binary forms, with or without +modification, are permitted provided that the following conditions +are met: + +1. Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + +2. Redistributions in binary form must reproduce the above copyright + notice, this list of conditions and the following disclaimer in the + documentation and/or other materials provided with the distribution. + +3. Neither the name of the California Institute of Technology nor + the names of its contributors may be used to endorse or promote + products derived from this software without specific prior + written permission. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT +LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS +FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL CALTECH +OR THE CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, +SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT +LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF +USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND +ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, +OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT +OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF +SUCH DAMAGE. + +Author: Richard M. Murray +Date: 24 May 09 +Revised: Kevin K. Chen, Dec 10 + +$Id: bdalg.py 17 2010-05-29 23:50:52Z murrayrm $ + +""" + import scipy as sp import xferfcn as tf import statesp as ss def series(sys1, sys2): - """Return the series connection sys1 * sys2 for --> sys2 --> sys1 -->.""" + """Return the series connection sys2 * sys1 for --> sys1 --> sys2 -->. + + Parameters + --------- + sys1: scalar, StateSpace, or TransferFunction + sys2: scalar, StateSpace, or TransferFunction + + Returns + ------- + out: scalar, StateSpace, or TransferFunction + + Raises + ------ + ValueError + if sys2.inputs does not equal sys1.outputs + + See Also + -------- + parallel + feedback + + Notes + ----- + This function is a wrapper for the __mul__ function in the StateSpace and + TransferFunction classes. The output type is usually the type of sys2. If + sys2 is a scalar, then the output type is the type of sys1. + + Examples + -------- + >>> sys3 = series(sys1, sys2) # Same as sys3 = sys2 * sys1. + + """ - return sys1 * sys2 + return sys2 * sys1 def parallel(sys1, sys2): - """Return the parallel connection sys1 + sys2.""" + """ + Return the parallel connection sys1 + sys2. + + Parameters + --------- + sys1: scalar, StateSpace, or TransferFunction + sys2: scalar, StateSpace, or TransferFunction + + Returns + ------- + out: scalar, StateSpace, or TransferFunction + + Raises + ------ + ValueError + if sys1 and sys2 do not have the same numbers of inputs and outputs + + See Also + -------- + series + feedback + Notes + ----- + This function is a wrapper for the __add__ function in the StateSpace and + TransferFunction classes. The output type is usually the type of sys1. If + sys1 is a scalar, then the output type is the type of sys2. + + Examples + ------- + >>> sys3 = parallel(sys1, sys2) # Same as sys3 = sys1 + sys2. + + """ + return sys1 + sys2 def negate(sys): - """Return the negative of a system.""" + """ + Return the negative of a system. + + Parameters + ---------- + sys: StateSpace or TransferFunction + + Returns + ------- + out: StateSpace or TransferFunction + + Notes + ----- + This function is a wrapper for the __neg__ function in the StateSpace and + TransferFunction classes. The output type is the same as the input type. + + Examples + -------- + >>> sys2 = negate(sys1) # Same as sys2 = -sys1. + + """ return -sys; def feedback(sys1, sys2, sign=-1): - """Feedback interconnection between two I/O systems. + """ + Feedback interconnection between two I/O systems. - Usage - ===== - sys = feedback(sys1, sys2) - sys = feedback(sys1, sys2, sign) - - Compute the system corresponding to a feedback interconnection between - sys1 and sys2. When sign is not specified, it assumes a value of -1 - (negative feedback). A sign of 1 indicates positive feedback. - Parameters ---------- - sys1, sys2: linsys - Linear input/output systems + sys1: scalar, StateSpace, or TransferFunction + The primary plant. + sys2: scalar, StateSpace, or TransferFunction + The feedback plant (often a feedback controller). sign: scalar - Feedback sign. + The sign of feedback. sign = -1 indicates negative feedback, and sign = 1 + indicates positive feedback. sign is an optional argument; it assumes a + value of -1 if not specified. - Return values - ------------- - sys: linsys + Returns + ------- + out: StateSpace or TransferFunction + Raises + ------ + ValueError + if sys1 does not have as many inputs as sys2 has outputs, or if sys2 + does not have as many inputs as sys1 has outputs + NotImplementedError + if an attempt is made to perform a feedback on a MIMO TransferFunction + object + + See Also + -------- + series + parallel + Notes - ----- - 1. This function calls calls xferfcn.feedback if sys1 is a TransferFunction - object and statesp.feedback if sys1 is a StateSpace object. If sys1 is a - scalar, then it is converted to sys2's type, and the corresponding - feedback function is used. If sys1 and sys2 are both scalars, then use - xferfcn.feedback.""" + ---- + This function is a wrapper for the feedback function in the StateSpace and + TransferFunction classes. It calls TransferFunction.feedback if sys1 is a + TransferFunction object, and StateSpace.feedback if sys1 is a StateSpace + object. If sys1 is a scalar, then it is converted to sys2's type, and the + corresponding feedback function is used. If sys1 and sys2 are both scalars, + then TransferFunction.feedback is used. + """ + # Check for correct input types. if not isinstance(sys1, (int, long, float, complex, tf.TransferFunction, ss.StateSpace)): Modified: branches/control-0.4a/src/lti.py =================================================================== --- branches/control-0.4a/src/lti.py 2011-02-08 22:15:45 UTC (rev 74) +++ branches/control-0.4a/src/lti.py 2011-02-08 22:15:51 UTC (rev 75) @@ -1,9 +1,48 @@ +"""lti.py + +The lti module contains the Lti parent class to the child classes StateSpace and +TransferFunction. It is designed for use in the python-control library. + +Routines in this module: + +Lti.__init__ + +""" + class Lti: - """The Lti is a parent class to the StateSpace and TransferFunction child - classes. It only contains the number of inputs and outputs, but this can be - expanded in the future.""" + + """Lti is a parent class to linear time invariant control (LTI) objects. + Lti is the parent to the StateSpace and TransferFunction child classes. It + contains the number of inputs and outputs, but this can be expanded in the + future. + + The StateSpace and TransferFunction child classes contain several common + "virtual" functions. These are: + + __init__ + __str__ + __neg__ + __add__ + __radd__ + __sub__ + __rsub__ + __mul__ + __rmul__ + __div__ + __rdiv__ + evalfr + freqresp + pole + zero + feedback + returnScipySignalLti + + """ + def __init__(self, inputs=1, outputs=1): + """Assign the LTI object's numbers of inputs and ouputs.""" + # Data members common to StateSpace and TransferFunction. self.inputs = inputs self.outputs = outputs Modified: branches/control-0.4a/src/matlab.py =================================================================== --- branches/control-0.4a/src/matlab.py 2011-02-08 22:15:45 UTC (rev 74) +++ branches/control-0.4a/src/matlab.py 2011-02-08 22:15:51 UTC (rev 75) @@ -1,53 +1,58 @@ -# matlab.py - MATLAB emulation functions -# -# Author: Richard M. Murray -# Date: 29 May 09 -# -# This file contains a number of functions that emulate some of the -# functionality of MATLAB. The intent of these functions is to -# provide a simple interface to the python control systems library -# (python-control) for people who are familiar with the MATLAB Control -# Systems Toolbox (tm). Most of the functions are just calls to -# python-control functions defined elsewhere. Use 'from -# control.matlab import *' in python to include all of the functions -# defined here. Functions that are defined in other libraries that -# have the same names as their MATLAB equivalents are automatically -# imported here. -# -# Copyright (c) 2009 by California Institute of Technology -# All rights reserved. -# -# Redistribution and use in source and binary forms, with or without -# modification, are permitted provided that the following conditions -# are met: -# -# 1. Redistributions of source code must retain the above copyright -# notice, this list of conditions and the following disclaimer. -# -# 2. Redistributions in binary form must reproduce the above copyright -# notice, this list of conditions and the following disclaimer in the -# documentation and/or other materials provided with the distribution. -# -# 3. Neither the name of the California Institute of Technology nor -# the names of its contributors may be used to endorse or promote -# products derived from this software without specific prior -# written permission. -# -# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS -# "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT -# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS -# FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL CALTECH -# OR THE CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, -# SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT -# LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF -# USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND -# ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, -# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT -# OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF -# SUCH DAMAGE. -# -# $Id: matlab.py 33 2010-11-26 21:59:57Z murrayrm $ +"""matlab.py +MATLAB emulation functions. + +This file contains a number of functions that emulate some of the +functionality of MATLAB. The intent of these functions is to +provide a simple interface to the python control systems library +(python-control) for people who are familiar with the MATLAB Control +Systems Toolbox (tm). Most of the functions are just calls to +python-control functions defined elsewhere. Use 'from +control.matlab import *' in python to include all of the functions +defined here. Functions that are defined in other libraries that +have the same names as their MATLAB equivalents are automatically +imported here. + +Copyright (c) 2009 by California Institute of Technology +All rights reserved. + +Redistribution and use in source and binary forms, with or without +modification, are permitted provided that the following conditions +are met: + +1. Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + +2. Redistributions in binary form must reproduce the above copyright + notice, this list of conditions and the following disclaimer in the + documentation and/or other materials provided with the distribution. + +3. Neither the name of the California Institute of Technology nor + the names of its contributors may be used to endorse or promote + products derived from this software without specific prior + written permission. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT +LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS +FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL CALTECH +OR THE CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, +SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT +LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF +USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND +ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, +OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT +OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF +SUCH DAMAGE. + +Author: Richard M. Murray +Date: 29 May 09 +Revised: Kevin K. Chen, Dec 10 + +$Id: matlab.py 33 2010-11-26 21:59:57Z murrayrm $ + +""" + # Libraries that we make use of import scipy as sp # SciPy library (used all over) import numpy as np # NumPy library @@ -72,7 +77,7 @@ from statefbk import ctrb, obsv, gram, place, lqr from delay import pade -__doc__ = """ +__doc__ += """ The control.matlab module defines functions that are roughly the equivalents of those in the MATLAB Control Toolbox. Items marked by a '*' are currently implemented; those marked with a '-' are not planned @@ -255,16 +260,44 @@ * linspace - generate a set of numbers that are linearly spaced * logspace - generate a set of numbers that are logarithmically spaced * unwrap - unwrap a phase angle to give a continuous curve + """ def ss(*args): - """Create a state space system. + """ + Create a state space system. - Usage - ===== - ss(A, B, C, D) - ss(sys)""" - + Parameters + ---------- + A: numpy matrix or matrix-like object + B: numpy matrix or matrix-like object + C: numpy matrix or matrix-like object + D: numpy matrix or matrix-like object + sys: StateSpace or TransferFunction object + ss accepts a set of A, B, C, D matrices or sys. + + Returns + ------- + out: StateSpace object + + Raises + ------ + ValueError + if matrix sizes are not self-consistent + + See Also + -------- + tf + ss2tf + tf2ss + + Examples + -------- + >>> sys = ss(A, B, C, D) # Create a StateSpace object from these matrices. + >>> sys = ss(sys1) # Convert a TransferFunction to a StateSpace object. + + """ + if len(args) == 4: return StateSpace(args[0], args[1], args[2], args[3]) elif len(args) == 1: @@ -275,21 +308,55 @@ return tf2ss(sys) else: raise TypeError("ss(sys): sys must be a StateSpace or \ -TransferFunction object.") +TransferFunction object. It is %s." % type(sys)) else: raise ValueError("Needs 1 or 4 arguments; received %i." % len(args)) def tf(*args): - """Create a transfer function system. + """ + Create a transfer function system. - Usage - ===== - tf(num, den) - tf(sys) - - num and den can be scalars or vectors for SISO systems, or lists of lists of - vectors for MIMO systems.""" + Parameters + ---------- + num: vector, or list of lists of vectors + den: vector, or list of lists of vectors + sys: StateSpace or TransferFunction object + tf accepts a num and den, or sys. + Returns + ------- + out: TransferFunction object + + Raises + ------ + ValueError + if num and den have invalid or unequal dimensions + TypeError + if num or den are of incorrect type + + See Also + -------- + ss + ss2tf + tf2ss + + Notes + -------- + num[i][j] is the vector of polynomial coefficients of the transfer function + numerator from the (j+1)st output to the (i+1)st input. den[i][j] works the + same way. + + Examples + -------- + >>> num = [[[1., 2.], [3., 4.]], [[5., 6.], [7., 8.]]] + >>> den = [[[9., 8., 7.], [6., 5., 4.]], [[3., 2., 1.], [-1., -2., -3.]]] + >>> sys = tf(num, den) + The transfer function from the 2nd input to the 1st output is + (3s + 4) / (6s^2 + 5s + 4). + >>> sys = tf(sys1) # Convert a StateSpace to a TransferFunction object. + + """ + if len(args) == 2: return TransferFunction(args[0], args[1]) elif len(args) == 1: @@ -300,38 +367,107 @@ return sys else: raise TypeError("tf(sys): sys must be a StateSpace or \ -TransferFunction object.") +TransferFunction object. It is %s." % type(sys)) else: raise ValueError("Needs 1 or 2 arguments; received %i." % len(args)) def ss2tf(*args): - """Transform a state space system to a transfer function. + """ + Transform a state space system to a transfer function. - Usage - ===== - ss2tf(A, B, C, D) - ss2tf(sys) - sys should have attributes A, B, C, D""" - + Parameters + ---------- + A: numpy matrix or matrix-like object + B: numpy matrix or matrix-like object + C: numpy matrix or matrix-like object + D: numpy matrix or matrix-like object + sys: StateSpace object + ss accepts a set of A, B, C, D matrices or a StateSpace object. + + Returns + ------- + out: TransferFunction object + + Raises + ------ + ValueError + if matrix sizes are not self-consistent, or if an invalid number of + arguments is passed in + TypeError + if sys is not a StateSpace object + + See Also + -------- + tf + ss + tf2ss + + Examples + -------- + >>> sys = ss2tf(A, B, C, D) + >>> sys = ss2tf(sys1) # Convert a StateSpace to a TransferFunction object. + + """ + if len(args) == 4: # Assume we were given the A, B, C, D matrix return convertToTransferFunction(StateSpace(args[0], args[1], args[2], args[3])) elif len(args) == 1: sys = args[0] - if not isinstance(sys, StateSpace): - raise TypeError("ss2tf(sys): sys must be a StateSpace object.") - return convertToTransferFunction(sys) + if isinstance(sys, StateSpace): + return convertToTransferFunction(sys) + else: + raise TypeError("ss2tf(sys): sys must be a StateSpace object. It \ +is %s." % type(sys)) else: raise ValueError("Needs 1 or 4 arguments; received %i." % len(args)) def tf2ss(*args): - """Transform a transfer function to a state space system. - - Usage - ===== - tf2ss(num, den) - tf2ss(sys) - sys should be a system object (lti or TransferFunction)""" + """ + Transform a transfer function to a state space system. + Parameters + ---------- + num: vector, or list of lists of vectors + den: vector, or list of lists of vectors + sys: TransferFunction object + tf2ss accepts num and den, or sys. + + Returns + ------- + out: StateSpace object + + Raises + ------ + ValueError + if num and den have invalid or unequal dimensions, or if an invalid + number of arguments is passed in + TypeError + if num or den are of incorrect type, or if sys is not a TransferFunction + object + + See Also + -------- + ss + tf + ss2tf + + Notes + -------- + num[i][j] is the vector of polynomial coefficients of the transfer function + numerator from the (j+1)st output to the (i+1)st input. den[i][j] works the + same way. + + Examples + -------- + >>> num = [[[1., 2.], [3., 4.]], [[5., 6.], [7., 8.]]] + >>> den = [[[9., 8., 7.], [6., 5., 4.]], [[3., 2., 1.], [-1., -2., -3.]]] + >>> sys = tf2ss(num, den) + >>> sys = tf2ss(sys1) # Convert a TransferFunction to a StateSpace object. + + """ + if len(args) == 2: # Assume we were given the num, den return convertToStateSpace(TransferFunction(args[0], args[1])) @@ -459,10 +595,6 @@ yout response of the system xout time evolution of the state vector """ - sys = args[0] - ltiobjs = sys.returnScipySignalLti() - ltiobj = ltiobjs[0][0] - return sp.signal.lsim2(*args, **keywords) #! Redefine step to use lsim2 @@ -483,19 +615,8 @@ T time values of the output yout response of the system """ - sys = args[0] - ltiobjs = sys.returnScipySignalLti() - ltiobj = ltiobjs[0][0] - newargs = [] - newargs.append(ltiobj) - for i in range(1, len(args)): - newargs.append(args[i]) - newargs = tuple(newargs) - print len(args) - print len(newargs) + return sp.signal.step(*args, **keywords) - return sp.signal.step(*newargs, **keywords) - # Redefine initial to use lsim2 #! Not yet implemented (uses step for now) def initial(*args, **keywords): @@ -514,10 +635,6 @@ T time values of the output yout response of the system """ - sys = args[0] - ltiobjs = sys.returnScipySignalLti() - ltiobj = ltiobjs[0][0] - return sp.signal.initial(*args, **keywords) # Redefine impulse to use initial() @@ -538,9 +655,5 @@ T time values of the output yout response of the system """ - sys = args[0] - ltiobjs = sys.returnScipySignalLti() - ltiobj = ltiobjs[0][0] + return sp.signal.impulse(*args, **keywords) - return sp.signal.impulse(ltiobj, **keywords) - Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:15:45 UTC (rev 74) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:15:51 UTC (rev 75) @@ -1,67 +1,104 @@ -# stateSpace.py - state space class for control systems library -# -# Author: Richard M. Murray -# Date: 24 May 09 -# -# This file contains the StateSpace class, which is used to represent -# linear systems in state space. This is the primary representation -# for the control system library. -# -# Copyright (c) 2010 by California Institute of Technology -# All rights reserved. -# -# Redistribution and use in source and binary forms, with or without -# modification, are permitted provided that the following conditions -# are met: -# -# 1. Redistributions of source code must retain the above copyright -# notice, this list of conditions and the following disclaimer. -# -# 2. Redistributions in binary form must reproduce the above copyright -# notice, this list of conditions and the following disclaimer in the -# documentation and/or other materials provided with the distribution. -# -# 3. Neither the name of the California Institute of Technology nor -# the names of its contributors may be used to endorse or promote -# products derived from this software without specific prior -# written permission. -# -# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS -# "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT -# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS -# FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL CALTECH -# OR THE CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, -# SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT -# LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF -# USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND -# ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, -# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT -# OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF -# SUCH DAMAGE. -# -# $Id: statesp.py 21 2010-06-06 17:29:42Z murrayrm $ +"""stateSpace.py -import scipy as sp -from scipy import concatenate, zeros -from numpy.linalg import solve +State space representation and functions. + +This file contains the StateSpace class, which is used to represent +linear systems in state space. This is the primary representation +for the python-control library. + +Routines in this module: + +StateSpace.__init__ +StateSpace.__str__ +StateSpace.__neg__ +StateSpace.__add__ +StateSpace.__radd__ +StateSpace.__sub__ +StateSpace.__rsub__ +StateSpace.__mul__ +StateSpace.__rmul__ +StateSpace.__div__ +StateSpace.__rdiv__ +StateSpace.evalfr +StateSpace.freqresp +StateSpace.pole +StateSpace.zero +StateSpace.feedback +StateSpace.returnScipySignalLti +convertToStateSpace +rss_generate + +Copyright (c) 2010 by California Institute of Technology +All rights reserved. + +Redistribution and use in source and binary forms, with or without +modification, are permitted provided that the following conditions +are met: + +1. Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + +2. Redistributions in binary form must reproduce the above copyright + notice, this list of conditions and the following disclaimer in the + documentation and/or other materials provided with the distribution. + +3. Neither the name of the California Institute of Technology nor + the names of its contributors may be used to endorse or promote + products derived from this software without specific prior + written permission. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT +LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS +FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL CALTECH +OR THE CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, +SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT +LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF +USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND +ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, +OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT +OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF +SUCH DAMAGE. + +Author: Richard M. Murray +Date: 24 May 09 +Revised: Kevin K. Chen, Dec 10 + +$Id: statepy 21 2010-06-06 17:29:42Z murrayrm $ + +""" + +from numpy import angle, any, concatenate, cos, dot, empty, exp, eye, pi, \ + poly, poly1d, matrix, roots, sin, zeros +from numpy.random import rand, randn +from numpy.linalg import inv, det, solve +from numpy.linalg.linalg import LinAlgError +from scipy.signal import lti +from lti import Lti import xferfcn -from lti import Lti class StateSpace(Lti): - """The StateSpace class is used throughout the python-control library to - represent systems in state space form. This class is derived from the Lti2 - base class.""" + """The StateSpace class represents state space instances and functions. + + The StateSpace class is used throughout the python-control library to + represent systems in state space form. This class is derived from the Lti + base class. + + The main data members are the A, B, C, and D matrices. The class also keeps + track of the number of states (i.e., the size of A). + + """ + def __init__(self, A=0, B=0, C=0, D=1): - """StateSpace constructor. The default constructor is the unit gain - direct feedthrough system.""" + """Construct a state space object. The default is unit static gain.""" # Here we're going to convert inputs to matrices, if the user gave a # non-matrix type. matrices = [A, B, C, D] for i in range(len(matrices)): # Convert to matrix first, if necessary. - matrices[i] = sp.matrix(matrices[i]) + matrices[i] = matrix(matrices[i]) [A, B, C, D] = matrices self.A = A @@ -85,7 +122,7 @@ raise ValueError("D must have the same row size as C.") def __str__(self): - """Style to use for printing.""" + """String representation of the state space.""" str = "A = " + self.A.__str__() + "\n\n" str += "B = " + self.B.__str__() + "\n\n" @@ -101,17 +138,18 @@ # Addition of two transfer functions (parallel interconnection) def __add__(self, other): - """Add two state space systems.""" + """Add two LTI systems (parallel connection).""" # Check for a couple of special cases if (isinstance(other, (int, long, float, complex))): # Just adding a scalar; put it in the D matrix A, B, C = self.A, self.B, self.C; D = self.D + other; + else: + other = convertToStateSpace(other) - else: # Check to make sure the dimensions are OK - if ((self.inputs != other.inputs) or \ + if ((self.inputs != other.inputs) or (self.outputs != other.outputs)): raise ValueError, "Systems have different shapes." @@ -125,52 +163,59 @@ B = concatenate((self.B, other.B), axis=0) C = concatenate((self.C, other.C), axis=1) D = self.D + other.D + return StateSpace(A, B, C, D) # Reverse addition - just switch the arguments def __radd__(self, other): - """Add two state space systems.""" + """Reverse add two LTI systems (parallel connection).""" - return self.__add__(other) + return self + other # Subtraction of two transfer functions (parallel interconnection) def __sub__(self, other): - """Subtract two state space systems.""" + """Subtract two LTI systems.""" - return self.__add__(-other) + return self + (-other) + def __rsub__(self, other): + """Reverse subtract two LTI systems.""" + + return other + (-self) + # Multiplication of two transfer functions (series interconnection) def __mul__(self, other): - """Serial interconnection between two state space systems.""" + """Multiply two LTI objects (serial connection).""" # Check for a couple of special cases if isinstance(other, (int, long, float, complex)): # Just multiplying by a scalar; change the output - A, B = self.A, self.B; - C = self.C * other; - D = self.D * other; + A, B = self.A, self.B + C = self.C * other + D = self.D * other else: - # Check to make sure the dimensions are OK - if (self.outputs != other.inputs): - raise ValueError, "Number of first's outputs must match number \ -of second's inputs." + other = convertToStateSpace(other) - # Concatenate the various arrays - A = concatenate( - (concatenate((other.A, zeros((other.A.shape[0], self.A.shape[1]))), - axis=1), - concatenate((self.B * other.C, self.A), axis=1)), - axis=0) - B = concatenate((other.B, self.B * other.D), axis=0) - C = concatenate((self.D * other.C, self.C),axis=1) - D = self.D * other.D + # Check to make sure the dimensions are OK + if self.inputs != other.outputs: + raise ValueError("C = A * B: A has %i column(s) (input(s)), \ +but B has %i row(s)\n(output(s))." % (self.inputs, other.outputs)) + # Concatenate the various arrays + A = concatenate( + (concatenate((other.A, zeros((other.A.shape[0], self.A.shape[1]))), + axis=1), + concatenate((self.B * other.C, self.A), axis=1)), axis=0) + B = concatenate((other.B, self.B * other.D), axis=0) + C = concatenate((self.D * other.C, self.C),axis=1) + D = self.D * other.D + return StateSpace(A, B, C, D) # Reverse multiplication of two transfer functions (series interconnection) # Just need to convert LH argument to a state space object def __rmul__(self, other): - """Serial interconnection between two state space systems""" + """Reverse multiply two LTI objects (serial connection).""" # Check for a couple of special cases if isinstance(other, (int, long, float, complex)): @@ -183,69 +228,88 @@ else: raise TypeError("can't interconnect systems") - # Compute poles and zeros - def pole(self): - """Compute the poles of a state space system.""" + # TODO: __div__ and __rdiv__ are not written yet. + def __div__(self, other): + """Divide two LTI systems.""" - return sp.roots(sp.poly(self.A)) + raise NotImplementedError("StateSpace.__div__ is not implemented yet.") - def zero(self): - """Compute the zeros of a state space system.""" + def __rdiv__(self, other): + """Reverse divide two LTI systems.""" - if self.inputs > 1 or self.outputs > 1: - raise NotImplementedError("StateSpace.zeros is currently \ -implemented only for SISO systems.") + raise NotImplementedError("StateSpace.__rdiv__ is not implemented yet.") - den = sp.poly1d(sp.poly(self.A)) - # 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, 0] - 1) * den) + def evalfr(self, omega): + """Evaluate a SS system's transfer function at a single frequency. - return (sp.roots(num)) + self.evalfr(omega) returns the value of the transfer function matrix with + input value s = i * omega. - 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, + fresp = self.C * solve(omega * 1.j * 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.""" + """Evaluate the system's transfer func. at a list of ang. frequencies. + + mag, phase, omega = self.freqresp(omega) + + reports the value of the magnitude, phase, and angular frequency of the + system's transfer function matrix evaluated at s = i * omega, where + omega is a list of angular frequencies. + + """ # 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) + mag = empty((self.outputs, self.inputs, numfreq)) + phase = empty((self.outputs, self.inputs, numfreq)) + fresp = empty((self.outputs, self.inputs, numfreq), dtype=complex) for k in range(numfreq): fresp[:, :, k] = self.evalfr(omega[k]) mag = abs(fresp) - phase = sp.angle(fresp) + phase = angle(fresp) return mag, phase, omega + # Compute poles and zeros + def pole(self): + """Compute the poles of a state space system.""" + + return roots(poly(self.A)) + + def zero(self): + """Compute the zeros of a state space system.""" + + if self.inputs > 1 or self.outputs > 1: + raise NotImplementedError("StateSpace.zeros is currently \ +implemented only for SISO systems.") + + den = poly1d(poly(self.A)) + # Compute the numerator based on zeros + #! TODO: This is currently limited to SISO systems + num = poly1d(\ + poly(self.A - dot(self.B, self.C)) + (self.D[0, 0] - 1) * den) + + return (roots(num)) + # Feedback around a state space system def feedback(self, other, sign=-1): - """Feedback interconnection between two state space systems.""" - - # Check for special cases - if (isinstance(other, (int, long, float, complex))): - # Scalar feedback, create state space system that is this case - other = StateSpace([[0]], [[0]], [[0]], [[ other ]]) + """Feedback interconnection between two LTI systems.""" + + other = convertToStateSpace(other) # Check to make sure the dimensions are OK if ((self.inputs != other.outputs) or (self.outputs != other.inputs)): raise ValueError, "State space systems don't have compatible \ inputs/outputs for feedback." - from numpy.linalg import inv, det - from numpy import eye - A1 = self.A B1 = self.B C1 = self.C @@ -276,17 +340,16 @@ return StateSpace(A, B, C, D) def returnScipySignalLti(self): - """Return a list of a list of scipy.signal.lti objects for a MIMO - system. For instance, + """Return a list of a list of scipy.signal.lti objects. + For instance, + >>> out = ssobject.returnScipySignalLti() >>> out[3][5] is a signal.scipy.lti object corresponding to the transfer function from the 6th input to the 4th output.""" - from scipy.signal import lti - # Preallocate the output. out = [[[] for j in range(self.inputs)] for i in range(self.outputs)] @@ -298,31 +361,39 @@ return out def convertToStateSpace(sys, inputs=1, outputs=1): - """Convert a system to state space form (if needed). If sys is a scalar, - then the number of inputs and outputs can be specified manually.""" + """Convert a system to state space form (if needed). + + If sys is already a state space, then it is returned. If sys is a transfer + function object, then it is converted to a state space and returned. If sys + is a scalar, then the number of inputs and outputs can be specified + manually. + """ + if isinstance(sys, StateSpace): # Already a state space system; just return it return sys elif isinstance(sys, xferfcn.TransferFunction): + # TODO: Wrap SLICOT to do transfer function to state space conversion. raise NotImplementedError("Transfer function to state space conversion \ is not implemented yet.") elif (isinstance(sys, (int, long, float, complex))): # Generate a simple state space system of the desired dimension # The following Doesn't work due to inconsistencies in ltisys: - # return StateSpace([[]], [[]], [[]], sp.eye(outputs, inputs)) + # return StateSpace([[]], [[]], [[]], eye(outputs, inputs)) return StateSpace(0., zeros((1, inputs)), zeros((outputs, 1)), - sys * sp.eye(outputs, inputs)) + sys * eye(outputs, inputs)) else: raise TypeError("Can't convert given type to StateSpace system.") def rss_generate(states, inputs, outputs, type): - """This does the actual random state space generation expected from rss and - drss. type is 'c' for continuous systems and 'd' for discrete systems.""" - - import numpy - from numpy.random import rand, randn + """Generate a random state space. + This does the actual random state space generation expected from rss and + drss. type is 'c' for continuous systems and 'd' for discrete systems. + + """ + # Probability of repeating a previous root. pRepeat = 0.05 # Probability of choosing a real root. Note that when choosing a complex @@ -337,7 +408,7 @@ pDzero = 0.5 # Make some poles for A. Preallocate a complex array. - poles = numpy.zeros(states) + numpy.zeros(states) * 0.j + poles = zeros(states) + zeros(states) * 0.j i = 0 while i < states: @@ -355,24 +426,24 @@ elif rand() < pReal or i == states - 1: # No-oscillation pole. if type == 'c': - poles[i] = -sp.exp(randn()) + 0.j + poles[i] = -exp(randn()) + 0.j elif type == 'd': poles[i] = 2. * rand() - 1. i += 1 else: # Complex conjugate pair of oscillating poles. if type == 'c': - poles[i] = complex(-sp.exp(randn()), 3. * sp.exp(randn())) + poles[i] = complex(-exp(randn()), 3. * exp(randn())) elif type == 'd': mag = rand() - phase = 2. * numpy.pi * rand() - poles[i] = complex(mag * numpy.cos(phase), - mag * numpy.sin(phase)) + phase = 2. * pi * rand() + poles[i] = complex(mag * cos(phase), + mag * sin(phase)) poles[i+1] = complex(poles[i].real, -poles[i].imag) i += 2 # Now put the poles in A as real blocks on the diagonal. - A = numpy.zeros((states, states)) + A = zeros((states, states)) i = 0 while i < states: if poles[i].imag == 0: @@ -387,9 +458,9 @@ while True: T = randn(states, states) try: - A = numpy.dot(solve(T, A), T) # A = T \ A * T + A = dot(solve(T, A), T) # A = T \ A * T break - except numpy.linalg.linalg.LinAlgError: + except LinAlgError: # In the unlikely event that T is rank-deficient, iterate again. pass @@ -401,14 +472,14 @@ # Make masks to zero out some of the elements. while True: Bmask = rand(states, inputs) < pBCmask - if sp.any(Bmask): # Retry if we get all zeros. + if any(Bmask): # Retry if we get all zeros. break while True: Cmask = rand(outputs, states) < pBCmask - if sp.any(Cmask): # Retry if we get all zeros. + if any(Cmask): # Retry if we get all zeros. break if rand() < pDzero: - Dmask = numpy.zeros((outputs, inputs)) + Dmask = zeros((outputs, inputs)) else: Dmask = rand(outputs, inputs) < pDmask Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:15:45 UTC (rev 74) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:15:51 UTC (rev 75) @@ -1,89 +1,124 @@ -# xferfcn.py - transfer function class and functions -# -# Author: Richard M. Murray -# Date: 24 May 09 -# -# This file contains the TransferFunction class and also functions -# that operate on transfer functions. This class extends the -# signal.lti class by adding some additional useful functions like -# block diagram algebra. -# -# NOTE: Transfer function in this module are restricted to be SISO -# systems. To perform calcualtiosn on MIMO systems, you need to use -# the state space module (statesp.py). -# -# Copyright (c) 2010 by California Institute of Technology -# All rights reserved. -# -# Redistribution and use in source and binary forms, with or without -# modification, are permitted provided that the following conditions -# are met: -# -# 1. Redistributions of source code must retain the above copyright -# notice, this list of conditions and the following disclaimer. -# -# 2. Redistributions in binary form must reproduce the above copyright -# notice, this list of conditions and the following disclaimer in the -# documentation and/or other materials provided with the distribution. -# -# 3. Neither the name of the California Institute of Technology nor -# the names of its contributors may be used to endorse or promote -# products derived from this software without specific prior -# written permission. -# -# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS -# "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT -# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS -# FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL CALTECH -# OR THE CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, -# SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT -# LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF -# USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND -# ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, -# OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT -# OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF -# SUCH DAMAGE. -# -# $Id: xferfcn.py 21 2010-06-06 17:29:42Z murrayrm $ +"""xferfcn.py +Transfer function representation and functions. + +This file contains the TransferFunction class and also functions +that operate on transfer functions. This is the primary representation +for the python-control library. + +Routines in this module: + +TransferFunction.__init__ +TransferFunction._truncatecoeff +TransferFunction.__str__ +TransferFunction.__neg__ +TransferFunction.__add__ +TransferFunction.__radd__ +TransferFunction.__sub__ +TransferFunction.__rsub__ +TransferFunction.__mul__ +TransferFunction.__rmul__ +TransferFunction.__div__ +TransferFunction.__rdiv__ +TransferFunction.evalfr +TransferFunction.freqresp +TransferFunction.pole +TransferFunction.zero +TransferFunction.feedback +TransferFunction.returnScipySignalLti +_tfpolyToString +_addSISO +convertToTransferFunction + +Copyright (c) 2010 by California Institute of Technology +All rights reserved. + +Redistribution and use in source and binary forms, with or without +modification, are permitted provided that the following conditions +are met: + +1. Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + +2. Redistributions in binary form must reproduce the above copyright + notice, this list of conditions and the following disclaimer in the + documentation and/or other materials provided with the distribution. + +3. Neither the name of the California Institute of Technology nor + the names of its contributors may be used to endorse or promote + products derived from this software without specific prior + written permission. + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +"AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT +LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS +FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL CALTECH +OR THE CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, +SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT +LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF +USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND +ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, +OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT +OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF +SUCH DAMAGE. + +Author: Richard M. Murray +Date: 24 May 09 +Revised: Kevin K. Chewn, Dec 10 + +$Id: xferfcn.py 21 2010-06-06 17:29:42Z murrayrm $ + +""" + # External function declarations -import scipy as sp +from numpy import angle, array, empty, ndarray, ones, polyadd, polymul, \ + polyval, zeros +from scipy.signal import lti +from copy import deepcopy from lti import Lti +import statesp class TransferFunction(Lti): - """The TransferFunction class is derived from the Lti2 parent class. The - main data members are 'num' and 'den', which are 2-D lists of arrays + + """The TransferFunction class represents TF instances and functions. + + The TransferFunction class is derived from the Lti parent class. It is used + throught the python-control library to represent systems in transfer + function form. + + The main data members are 'num' and 'den', which are 2-D lists of arrays containing MIMO numerator and denominator coefficients. For example, >>> num[2][5] = numpy.array([1., 4., 8.]) means that the numerator of the transfer function from the 6th input to the - 3rd output is set to s^2 + 4s + 8.""" + 3rd output is set to s^2 + 4s + 8. + """ + def __init__(self, num=1, den=1): - """This is the constructor. The default transfer function is 1 (unit - gain direct feedthrough).""" - + """Construct a transfer function. The default is unit static gain.""" + # Make num and den into lists of lists of arrays, if necessary. Beware: # this is a shallow copy! This should be okay, but be careful. data = [num, den] for i in range(len(data)): if isinstance(data[i], (int, float, long, complex)): # Convert scalar to list of list of array. - data[i] = [[sp.array([data[i]])]] - elif isinstance(data[i], (list, tuple, sp.ndarray)) and \ + data[i] = [[array([data[i]])]] + elif isinstance(data[i], (list, tuple, ndarray)) and \ isinstance(data[i][0], (int, float, long, complex)): # Convert array to list of list of array. - data[i] = [[sp.array(data[i])]] + data[i] = [[array(data[i])]] elif isinstance(data[i], list) and \ isinstance(data[i][0], list) and \ - isinstance(data[i][0][0], (list, tuple, sp.ndarray)) and \ + isinstance(data[i][0][0], (list, tuple, ndarray)) and \ isinstance(data[i][0][0][0], (int, float, long, complex)): # We might already have the right format. Convert the # coefficient vectors to arrays, if necessary. for j in range(len(data[i])): for k in range(len(data[i][j])): - data[i][j][k] = sp.array(data[i][j][k]) + data[i][j][k] = array(data[i][j][k]) else: # If the user passed in anything else, then it's unclear what # the meaning is. @@ -134,7 +169,7 @@ zeronum = False break if zeronum: - den[i][j] = sp.ones(1) + den[i][j] = ones(1) self.num = num self.den = den @@ -142,6 +177,35 @@ self._truncatecoeff() + def _truncatecoeff(self): + """Remove extraneous zero coefficients from num and den. + + Check every element of the numerator and denominator matrices, and + truncate leading zeros. For instance, running self._truncatecoeff() + will reduce self.num = [[[0, 0, 1, 2]]] to [[[1, 2]]]. + + """ + + # Beware: this is a shallow copy. This should be okay. + data = [self.num, self.den] + for p in range(len(data)): + for i in range(self.outputs): + for j in range(self.inputs): + # Find the first nontrivial coefficient. + nonzero = None + for k in range(data[p][i][j].size): + if data[p][i][j][k]: + nonzero = k + break + + if nonzero is None: + # The array is all zeros. + data[p][i][j] = zeros(1) + else: + # Truncate the trivial coefficients. + data[p][i][j] = data[p][i][j][nonzero:] + [self.num, self.den] = data + def __str__(self): """String representation of the transfer function.""" @@ -174,36 +238,10 @@ outstr += "\n" + numstr + "\n" + dashes + "\n" + denstr + "\n" return outstr - def _truncatecoeff(self): - """Remove extraneous zero coefficients from polynomials in numerator and - denominator matrices.""" - - # Beware: this is a shallow copy. This should be okay. - data = [self.num, self.den] - for p in range(len(data)): - for i in range(self.outputs): - for j in range(self.inputs): - # Find the first nontrivial coefficient. - nonzero = None - for k in range(data[p][i][j].size): - if data[p][i][j][k]: - nonzero = k - break - - if nonzero is None: - # The array is all zeros. - data[p][i][j] = sp.zeros(1) - else: - # Truncate the trivial coefficients. - data[p][i][j] = data[p][i][j][nonzero:] - [self.num, self.den] = data - def __neg__(self): """Negate a transfer function.""" - import copy - - num = copy.deepcopy(self.num) + num = deepcopy(self.num) for i in range(self.outputs): for j in range(self.inputs): num[i][j] *= -1 @@ -211,7 +249,7 @@ return TransferFunction(num, self.den) def __add__(self, other): - """Add two transfer functions (parallel connection).""" + """Add two LTI objects (parallel connection).""" # Convert the second argument to a transfer function. if not isinstance(other, TransferFunction): @@ -237,26 +275,28 @@ return TransferFunction(num, den) def __radd__(self, other): - """Add two transfer functions (parallel connection).""" + """Reverse add two LTI objects (parallel connection).""" return self + other; def __sub__(self, other): - """Subtract two transfer functions.""" + """Subtract two LTI objects.""" return self + (-other) def __rsub__(self, other): - """Subtract two transfer functions.""" + """Reverse subtract two LTI objects.""" return other + (-self) def __mul__(self, other): - """Multiply two transfer functions (serial connection).""" + """Multiply two LTI objects (serial connection).""" # Convert the second argument to a transfer function. - if not isinstance(other, TransferFunction): + if isinstance(other, (int, float, long, complex)): other = convertToTransferFunction(other, self.inputs, self.inputs) + else: + other = convertToTransferFunction(other) # Check that the input-output sizes are consistent. if self.inputs != other.outputs: @@ -278,21 +318,21 @@ for i in range(outputs): # Iterate through rows of product. for j in range(inputs): # Iterate through columns of product. for k in range(self.inputs): # Multiply & add. - num_summand[k] = sp.polymul(self.num[i][k], other.num[k][j]) - den_summand[k] = sp.polymul(self.den[i][k], other.den[k][j]) + num_summand[k] = polymul(self.num[i][k], other.num[k][j]) + den_summand[k] = polymul(self.den[i][k], other.den[k][j]) num[i][j], den[i][j] = _addSISO(num[i][j], den[i][j], num_summand[k], den_summand[k]) return TransferFunction(num, den) def __rmul__(self, other): - """Multiply two transfer functions (serial connection).""" + """Reverse multiply two LTI objects (serial connection).""" return self * other - # TODO: Division of MIMO transfer function objects is quite difficult. + # TODO: Division of MIMO transfer function objects is not written yet. def __div__(self, other): - """Divide two transfer functions.""" + """Divide two LTI objects.""" if self.inputs > 1 or self.outputs > 1 or \ other.inputs > 1 or other.outputs > 1: @@ -300,17 +340,16 @@ implemented only for SISO systems.") # Convert the second argument to a transfer function. - if not isinstance(other, TransferFunction): - other = convertToTransferFunction(other, 1, 1) + other = convertToTransferFunction(other) - num = sp.polymul(self.num[0][0], other.den[0][0]) - den = sp.polymul(self.den[0][0], other.num[0][0]) + num = polymul(self.num[0][0], other.den[0][0]) + den = polymul(self.den[0][0], other.num[0][0]) return TransferFunction(num, den) - # TODO: Division of MIMO transfer function objects is quite difficult. + # TODO: Division of MIMO transfer function objects is not written yet. def __rdiv__(self, other): - """Reverse divide two transfer functions.""" + """Reverse divide two LTI objects.""" if self.inputs > 1 or self.outputs > 1 or \ other.inputs > 1 or other.outputs > 1: @@ -319,56 +358,70 @@ return other / self - def evalfr(self, freq): - """Evaluate a transfer function at a single frequency.""" + def evalfr(self, omega): + """Evaluate a transfer function... [truncated message content] |
From: <kk...@us...> - 2011-02-08 22:15:52
|
Revision: 74 http://python-control.svn.sourceforge.net/python-control/?rev=74&view=rev Author: kkchen Date: 2011-02-08 22:15:45 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Modified lsim and initial to use returnScipySignalLti... all functions using scipy.signal.lti still need modification for their arguments list and MIMO Steven Brunton <sbr...@pr...> Modified Paths: -------------- branches/control-0.4a/src/matlab.py Modified: branches/control-0.4a/src/matlab.py =================================================================== --- branches/control-0.4a/src/matlab.py 2011-02-08 22:15:42 UTC (rev 73) +++ branches/control-0.4a/src/matlab.py 2011-02-08 22:15:45 UTC (rev 74) @@ -459,6 +459,10 @@ yout response of the system xout time evolution of the state vector """ + sys = args[0] + ltiobjs = sys.returnScipySignalLti() + ltiobj = ltiobjs[0][0] + return sp.signal.lsim2(*args, **keywords) #! Redefine step to use lsim2 @@ -510,6 +514,10 @@ T time values of the output yout response of the system """ + sys = args[0] + ltiobjs = sys.returnScipySignalLti() + ltiobj = ltiobjs[0][0] + return sp.signal.initial(*args, **keywords) # Redefine impulse to use initial() This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:15:48
|
Revision: 73 http://python-control.svn.sourceforge.net/python-control/?rev=73&view=rev Author: kkchen Date: 2011-02-08 22:15:42 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Added changes to step and impulse from earlier... still needs work with argument list Steven Brunton <sbr...@pr...> Modified Paths: -------------- branches/control-0.4a/src/matlab.py Modified: branches/control-0.4a/src/matlab.py =================================================================== --- branches/control-0.4a/src/matlab.py 2011-02-08 22:15:38 UTC (rev 72) +++ branches/control-0.4a/src/matlab.py 2011-02-08 22:15:42 UTC (rev 73) @@ -479,8 +479,19 @@ T time values of the output yout response of the system """ - return sp.signal.step(*args, **keywords) + sys = args[0] + ltiobjs = sys.returnScipySignalLti() + ltiobj = ltiobjs[0][0] + newargs = [] + newargs.append(ltiobj) + for i in range(1, len(args)): + newargs.append(args[i]) + newargs = tuple(newargs) + print len(args) + print len(newargs) + return sp.signal.step(*newargs, **keywords) + # Redefine initial to use lsim2 #! Not yet implemented (uses step for now) def initial(*args, **keywords): @@ -519,5 +530,9 @@ T time values of the output yout response of the system """ - return sp.signal.impulse(*args, **keywords) + sys = args[0] + ltiobjs = sys.returnScipySignalLti() + ltiobj = ltiobjs[0][0] + return sp.signal.impulse(ltiobj, **keywords) + This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:15:44
|
Revision: 72 http://python-control.svn.sourceforge.net/python-control/?rev=72&view=rev Author: kkchen Date: 2011-02-08 22:15:38 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Added StateSpace binary operator tests to TestStateSp.py; bug fixes. - {StateSpace,XferFcn}.{poles(),zeros()} changed to pole(), zero() to be consistent with MATLAB. - Bug fix in StateSpace.zero(). - Bug fix in StateSpace.__sub__(). - Bug fix in StateSpace.__mul__(). - Rearranged the StateSpace member functions to a more sensible order. Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/TestStateSp.py branches/control-0.4a/src/matlab.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:15:32 UTC (rev 71) +++ branches/control-0.4a/src/TestStateSp.py 2011-02-08 22:15:38 UTC (rev 72) @@ -8,6 +8,89 @@ class TestStateSpace(unittest.TestCase): """Tests for the StateSpace class.""" + def setUp(self): + """Set up a MIMO system to test operations on.""" + + A = [[-3., 4., 2.], [-1., -3., 0.], [2., 5., 3.]] + B = [[1., 4.], [-3., -3.], [-2., 1.]] + C = [[4., 2., -3.], [1., 4., 3.]] + D = [[-2., 4.], [0., 1.]] + + a = [[4., 1.], [2., -3]] + b = [[5., 2.], [-3., -3.]] + c = [[2., -4], [0., 1.]] + d = [[3., 2.], [1., -1.]] + + self.sys1 = StateSpace(A, B, C, D) + self.sys2 = StateSpace(a, b, c, d) + + def testPole(self): + """Evaluate the poles of a MIMO system.""" + + p = self.sys1.pole() + + np.testing.assert_array_almost_equal(p, [3.34747678408874, + -3.17373839204437 + 1.47492908003839j, + -3.17373839204437 - 1.47492908003839j]) + + def testZero(self): + """Evaluate the zeros of a SISO system.""" + + sys = StateSpace(self.sys1.A, [[3.], [-2.], [4.]], [[-1., 3., 2.]], [[-4.]]) + z = sys.zero() + + np.testing.assert_array_almost_equal(z, [4.26864638637134, + -3.75932319318567 + 1.10087776649554j, + -3.75932319318567 - 1.10087776649554j]) + + def testAdd(self): + """Add two MIMO systems.""" + + A = [[-3., 4., 2., 0., 0.], [-1., -3., 0., 0., 0.], + [2., 5., 3., 0., 0.], [0., 0., 0., 4., 1.], [0., 0., 0., 2., -3.]] + B = [[1., 4.], [-3., -3.], [-2., 1.], [5., 2.], [-3., -3.]] + C = [[4., 2., -3., 2., -4.], [1., 4., 3., 0., 1.]] + D = [[1., 6.], [1., 0.]] + + sys = self.sys1 + self.sys2 + + np.testing.assert_array_almost_equal(sys.A, A) + np.testing.assert_array_almost_equal(sys.B, B) + np.testing.assert_array_almost_equal(sys.C, C) + np.testing.assert_array_almost_equal(sys.D, D) + + def testSub(self): + """Subtract two MIMO systems.""" + + A = [[-3., 4., 2., 0., 0.], [-1., -3., 0., 0., 0.], + [2., 5., 3., 0., 0.], [0., 0., 0., 4., 1.], [0., 0., 0., 2., -3.]] + B = [[1., 4.], [-3., -3.], [-2., 1.], [5., 2.], [-3., -3.]] + C = [[4., 2., -3., -2., 4.], [1., 4., 3., 0., -1.]] + D = [[-5., 2.], [-1., 2.]] + + sys = self.sys1 - self.sys2 + + np.testing.assert_array_almost_equal(sys.A, A) + np.testing.assert_array_almost_equal(sys.B, B) + np.testing.assert_array_almost_equal(sys.C, C) + np.testing.assert_array_almost_equal(sys.D, D) + + def testMul(self): + """Multiply two MIMO systems.""" + + A = [[4., 1., 0., 0., 0.], [2., -3., 0., 0., 0.], [2., 0., -3., 4., 2.], + [-6., 9., -1., -3., 0.], [-4., 9., 2., 5., 3.]] + B = [[5., 2.], [-3., -3.], [7., -2.], [-12., -3.], [-5., -5.]] + C = [[-4., 12., 4., 2., -3.], [0., 1., 1., 4., 3.]] + D = [[-2., -8.], [1., -1.]] + + sys = self.sys1 * self.sys2 + + np.testing.assert_array_almost_equal(sys.A, A) + np.testing.assert_array_almost_equal(sys.B, B) + np.testing.assert_array_almost_equal(sys.C, C) + np.testing.assert_array_almost_equal(sys.D, D) + def testEvalFr(self): """Evaluate the frequency response at one frequency.""" @@ -79,7 +162,7 @@ for inputs in range(1, self.maxIO): for outputs in range(1, self.maxIO): sys = matlab.rss(states, inputs, outputs) - p = sys.poles() + p = sys.pole() for z in p: self.assertTrue(z.real < 0) @@ -113,7 +196,7 @@ for inputs in range(1, self.maxIO): for outputs in range(1, self.maxIO): sys = matlab.drss(states, inputs, outputs) - p = sys.poles() + p = sys.pole() for z in p: self.assertTrue(abs(z) < 1) Modified: branches/control-0.4a/src/matlab.py =================================================================== --- branches/control-0.4a/src/matlab.py 2011-02-08 22:15:32 UTC (rev 71) +++ branches/control-0.4a/src/matlab.py 2011-02-08 22:15:38 UTC (rev 72) @@ -357,8 +357,13 @@ def pole(sys): """Return system poles.""" - return sys.poles() + return sys.pole() +def zero(sys): + """Return system zeros.""" + + return sys.zero() + def evalfr(sys, omega): """Evaluate the transfer function of an LTI system at a single frequency omega.""" Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:15:32 UTC (rev 71) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:15:38 UTC (rev 72) @@ -93,52 +93,6 @@ 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.""" - - # 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) - - for k in range(numfreq): - fresp[:, :, k] = self.evalfr(omega[k]) - - mag = abs(fresp) - phase = sp.angle(fresp) - - return mag, phase, omega - - # Compute poles and zeros - def poles(self): - """Compute the poles of a state space system.""" - - return sp.roots(sp.poly(self.A)) - - def zeros(self): - """Compute the zeros of a state space system.""" - - if self.inputs > 1 or self.outputs > 1: - raise NotImplementedError("StateSpace.zeros is currently \ -implemented only for SISO systems.") - - den = sp.poly1d(sp.poly(self.A)) - # 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) - - return (sp.roots(num)) - # Negation of a system def __neg__(self): """Negate a state space system.""" @@ -183,19 +137,18 @@ def __sub__(self, other): """Subtract two state space systems.""" - return __add__(self, other.__neg__()) + return self.__add__(-other) # Multiplication of two transfer functions (series interconnection) def __mul__(self, other): """Serial interconnection between two state space systems.""" # Check for a couple of special cases - if (isinstance(other, (int, long, float, complex))): + if isinstance(other, (int, long, float, complex)): # Just multiplying by a scalar; change the output A, B = self.A, self.B; C = self.C * other; D = self.D * other; - else: # Check to make sure the dimensions are OK if (self.outputs != other.inputs): @@ -203,14 +156,15 @@ of second's inputs." # Concatenate the various arrays - A = concatenate(( - concatenate(( self.A, zeros((self.A.shape[0], - other.A.shape[-1])) ),axis=1), - concatenate(( other.B*self.C, other.A ),axis=1), - ),axis=0) - B = concatenate( (self.B, other.B*self.D), axis=0 ) - C = concatenate( (other.D*self.C, other.C), axis=1 ) - D = other.D*self.D + A = concatenate( + (concatenate((other.A, zeros((other.A.shape[0], self.A.shape[1]))), + axis=1), + concatenate((self.B * other.C, self.A), axis=1)), + axis=0) + B = concatenate((other.B, self.B * other.D), axis=0) + C = concatenate((self.D * other.C, self.C),axis=1) + D = self.D * other.D + return StateSpace(A, B, C, D) # Reverse multiplication of two transfer functions (series interconnection) @@ -219,7 +173,7 @@ """Serial interconnection between two state space systems""" # Check for a couple of special cases - if (isinstance(other, (int, long, float, complex))): + if isinstance(other, (int, long, float, complex)): # Just multiplying by a scalar; change the input A, C = self.A, self.C; B = self.B * other; @@ -229,6 +183,52 @@ else: raise TypeError("can't interconnect systems") + # Compute poles and zeros + def pole(self): + """Compute the poles of a state space system.""" + + return sp.roots(sp.poly(self.A)) + + def zero(self): + """Compute the zeros of a state space system.""" + + if self.inputs > 1 or self.outputs > 1: + raise NotImplementedError("StateSpace.zeros is currently \ +implemented only for SISO systems.") + + den = sp.poly1d(sp.poly(self.A)) + # 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, 0] - 1) * den) + + return (sp.roots(num)) + + 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.""" + + # 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) + + for k in range(numfreq): + fresp[:, :, k] = self.evalfr(omega[k]) + + mag = abs(fresp) + phase = sp.angle(fresp) + + return mag, phase, omega + # Feedback around a state space system def feedback(self, other, sign=-1): """Feedback interconnection between two state space systems.""" Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:15:32 UTC (rev 71) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:15:38 UTC (rev 72) @@ -352,12 +352,12 @@ return mag, phase, omega - def poles(self): + def pole(self): """Compute poles of a transfer function.""" pass - def zeros(self): + def zero(self): """Compute zeros of a transfer function.""" pass This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:15:39
|
Revision: 71 http://python-control.svn.sourceforge.net/python-control/?rev=71&view=rev Author: kkchen Date: 2011-02-08 22:15:32 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Renamed Lti2 to Lti. Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/statesp.py branches/control-0.4a/src/xferfcn.py Added Paths: ----------- branches/control-0.4a/src/lti.py Removed Paths: ------------- branches/control-0.4a/src/lti2.py Added: branches/control-0.4a/src/lti.py =================================================================== --- branches/control-0.4a/src/lti.py (rev 0) +++ branches/control-0.4a/src/lti.py 2011-02-08 22:15:32 UTC (rev 71) @@ -0,0 +1,9 @@ +class Lti: + """The Lti is a parent class to the StateSpace and TransferFunction child + classes. It only contains the number of inputs and outputs, but this can be + expanded in the future.""" + + def __init__(self, inputs=1, outputs=1): + # Data members common to StateSpace and TransferFunction. + self.inputs = inputs + self.outputs = outputs Deleted: branches/control-0.4a/src/lti2.py =================================================================== --- branches/control-0.4a/src/lti2.py 2011-02-08 22:15:26 UTC (rev 70) +++ branches/control-0.4a/src/lti2.py 2011-02-08 22:15:32 UTC (rev 71) @@ -1,9 +0,0 @@ -class Lti2: - """The Lti2 is a parent class to the StateSpace and TransferFunction child - classes. It only contains the number of inputs and outputs, but this can be - expanded in the future.""" - - def __init__(self, inputs=1, outputs=1): - # Data members common to StateSpace and TransferFunction. - self.inputs = inputs - self.outputs = outputs \ No newline at end of file Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:15:26 UTC (rev 70) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:15:32 UTC (rev 71) @@ -45,9 +45,9 @@ from scipy import concatenate, zeros from numpy.linalg import solve import xferfcn -from lti2 import Lti2 +from lti import Lti -class StateSpace(Lti2): +class StateSpace(Lti): """The StateSpace class is used throughout the python-control library to represent systems in state space form. This class is derived from the Lti2 base class.""" @@ -70,7 +70,7 @@ self.D = D self.states = A.shape[0] - Lti2.__init__(self, B.shape[1], C.shape[0]) + Lti.__init__(self, B.shape[1], C.shape[0]) # Check that the matrix sizes are consistent. if self.states != A.shape[1]: Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:15:26 UTC (rev 70) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:15:32 UTC (rev 71) @@ -48,9 +48,9 @@ # External function declarations import scipy as sp -from lti2 import Lti2 +from lti import Lti -class TransferFunction(Lti2): +class TransferFunction(Lti): """The TransferFunction class is derived from the Lti2 parent class. The main data members are 'num' and 'den', which are 2-D lists of arrays containing MIMO numerator and denominator coefficients. For example, @@ -138,7 +138,7 @@ self.num = num self.den = den - Lti2.__init__(self, inputs, outputs) + Lti.__init__(self, inputs, outputs) self._truncatecoeff() This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:15:32
|
Revision: 70 http://python-control.svn.sourceforge.net/python-control/?rev=70&view=rev Author: kkchen Date: 2011-02-08 22:15:26 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Added modred to modelsimp. Tests pass for MatchDC and Truncate methods. Steven Brunton <sbr...@pr...> Modified Paths: -------------- branches/control-0.4a/src/TestModelsimp.py branches/control-0.4a/src/modelsimp.py Modified: branches/control-0.4a/src/TestModelsimp.py =================================================================== --- branches/control-0.4a/src/TestModelsimp.py 2011-02-08 22:15:21 UTC (rev 69) +++ branches/control-0.4a/src/TestModelsimp.py 2011-02-08 22:15:26 UTC (rev 70) @@ -35,15 +35,15 @@ C = np.matrix('-0.9057, -0.4068, 0.3263, -0.3474') D = np.matrix('0.') sys = ss(A,B,C,D) - rsys = modred(sys,np.matrix('3, 4'),'matchdc') + rsys = modred(sys,[2, 3],'matchdc') Artrue = np.matrix('-4.431, -4.552; -4.552, -5.361') Brtrue = np.matrix('-1.362; -1.031') Crtrue = np.matrix('-1.362, -1.031') Drtrue = np.matrix('-0.08384') - np.testing.assert_array_almost_equal(sys.A, Artrue) - np.testing.assert_array_almost_equal(sys.B, Brtrue) - np.testing.assert_array_almost_equal(sys.C, Crtrue) - np.testing.assert_array_almost_equal(sys.D, Drtrue) + np.testing.assert_array_almost_equal(rsys.A, Artrue,decimal=3) + np.testing.assert_array_almost_equal(rsys.B, Brtrue,decimal=3) + np.testing.assert_array_almost_equal(rsys.C, Crtrue,decimal=3) + np.testing.assert_array_almost_equal(rsys.D, Drtrue,decimal=2) def testModredTruncate(self): #balanced realization computed in matlab for the transfer function: @@ -56,15 +56,15 @@ C = np.matrix('-0.9057, -0.4068, 0.3263, -0.3474') D = np.matrix('0.') sys = ss(A,B,C,D) - rsys = modred(sys,np.matrix('3, 4'),'truncate') + rsys = modred(sys,[2, 3],'truncate') Artrue = np.matrix('-1.958, -1.194; -1.194, -0.8344') - Brtrue = np.matrix('-0.9057, -0.4068') + Brtrue = np.matrix('-0.9057; -0.4068') Crtrue = np.matrix('-0.9057, -0.4068') Drtrue = np.matrix('0.') - np.testing.assert_array_almost_equal(sys.A, Artrue) - np.testing.assert_array_almost_equal(sys.B, Brtrue) - np.testing.assert_array_almost_equal(sys.C, Crtrue) - np.testing.assert_array_almost_equal(sys.D, Drtrue) + np.testing.assert_array_almost_equal(rsys.A, Artrue) + np.testing.assert_array_almost_equal(rsys.B, Brtrue) + np.testing.assert_array_almost_equal(rsys.C, Crtrue) + np.testing.assert_array_almost_equal(rsys.D, Drtrue) def testBalredMatchDC(self): #controlable canonical realization computed in matlab for the transfer function: Modified: branches/control-0.4a/src/modelsimp.py =================================================================== --- branches/control-0.4a/src/modelsimp.py 2011-02-08 22:15:21 UTC (rev 69) +++ branches/control-0.4a/src/modelsimp.py 2011-02-08 22:15:26 UTC (rev 70) @@ -109,27 +109,57 @@ # else: dico = 'C' - #TODO: Check system is stable, perhaps a utility in ctrlutil.py - # or a method of the StateSpace class? + #Check system is stable D,V = np.linalg.eig(sys.A) for e in D: if e.real >= 0: raise ValueError, "Oops, the system is unstable!" + print ELIM + ELIM = np.sort(ELIM) + print ELIM + NELIM = [] + # Create list of elements not to eliminate (NELIM) + for i in range(0,len(sys.A)): + if i not in ELIM: + NELIM.append(i) + print NELIM + # A1 is a matrix of all columns of sys.A not to eliminate + A1 = sys.A[:,NELIM[0]] + for i in NELIM[1:]: + A1 = np.hstack((A1, sys.A[:,i])) + A11 = A1[NELIM,:] + A21 = A1[ELIM,:] + # A2 is a matrix of all columns of sys.A to eliminate + A2 = sys.A[:,ELIM[0]] + for i in ELIM[1:]: + A2 = np.hstack((A2, sys.A[:,i])) + A12 = A2[NELIM,:] + A22 = A2[ELIM,:] + C1 = sys.C[:,NELIM] + C2 = sys.C[:,ELIM] + B1 = sys.B[NELIM,:] + B2 = sys.B[ELIM,:] + + print np.shape(A22) + A22I = np.linalg.inv(A22) + if method=='matchdc': - raise ValueError, "Not supported yet!" + # if matchdc, residualize + Ar = A11 - A12*A22.I*A21 + Br = B1 - A12*A22.I*B2 + Cr = C1 - C2*A22.I*A21 + Dr = sys.D - C2*A22.I*B2 elif method=='truncate': - raise ValueError, "Not supported yet!" + # if truncate, simply discard state x2 + Ar = A11 + Br = B1 + Cr = C1 + Dr = sys.D else: raise ValueError, "Oops, method is not supported!" - #Compute rhs using the Slycot routine XXXXXX - #make sure Slycot is installed - #try: - # from slycot import XXXXXX - #except ImportError: - # raise ControlSlycot("can't find slycot module 'XXXXXX'") - rsys = 0. + rsys = StateSpace(Ar,Br,Cr,Dr) return rsys def balred(sys,orders,method='truncate'): @@ -162,8 +192,7 @@ # else: dico = 'C' - #TODO: Check system is stable, perhaps a utility in ctrlutil.py - # or a method of the StateSpace class? + #Check system is stable D,V = np.linalg.eig(sys.A) for e in D: if e.real >= 0: This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:15:28
|
Revision: 69 http://python-control.svn.sourceforge.net/python-control/?rev=69&view=rev Author: kkchen Date: 2011-02-08 22:15:21 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Fixed hsv output, changed balred inputs to be less redundant (no longer takes 'elimination' argument, method='truncate' instead). Steven Brunton <sbr...@pr...> Modified Paths: -------------- branches/control-0.4a/src/TestModelsimp.py branches/control-0.4a/src/modelsimp.py branches/control-0.4a/src/statefbk.py Modified: branches/control-0.4a/src/TestModelsimp.py =================================================================== --- branches/control-0.4a/src/TestModelsimp.py 2011-02-08 22:15:17 UTC (rev 68) +++ branches/control-0.4a/src/TestModelsimp.py 2011-02-08 22:15:21 UTC (rev 69) @@ -78,7 +78,7 @@ D = np.matrix('0.') sys = ss(A,B,C,D) orders = 1 - rsys = balred(sys,orders,'elimination','matchdc') + rsys = balred(sys,orders,method='matchdc') Artrue = np.matrix('-0.566') Brtrue = np.matrix('-0.414') Crtrue = np.matrix('-0.5728') @@ -100,7 +100,7 @@ D = np.matrix('0.') sys = ss(A,B,C,D) orders = 2 - rsys = balred(sys,orders,'elimination','truncate') + rsys = balred(sys,orders,method='truncate') Artrue = np.matrix('-1.958, -1.194; -1.194, -0.8344') Brtrue = np.matrix('0.9057; 0.4068') Crtrue = np.matrix('0.9057, 0.4068') Modified: branches/control-0.4a/src/modelsimp.py =================================================================== --- branches/control-0.4a/src/modelsimp.py 2011-02-08 22:15:17 UTC (rev 68) +++ branches/control-0.4a/src/modelsimp.py 2011-02-08 22:15:21 UTC (rev 69) @@ -75,6 +75,8 @@ hsv = np.sqrt(w) hsv = np.matrix(hsv) + hsv = np.sort(hsv) + hsv = np.fliplr(hsv) # Return the Hankel singular values return hsv @@ -115,9 +117,9 @@ raise ValueError, "Oops, the system is unstable!" if method=='matchdc': - print "matchdc" + raise ValueError, "Not supported yet!" elif method=='truncate': - print "truncate" + raise ValueError, "Not supported yet!" else: raise ValueError, "Oops, method is not supported!" @@ -130,7 +132,7 @@ rsys = 0. return rsys -def balred(sys,orders,elimination,method): +def balred(sys,orders,method='truncate'): """Balanced reduced order model of sys of a given order. States are eliminated based on Hankel singular value. Usage Modified: branches/control-0.4a/src/statefbk.py =================================================================== --- branches/control-0.4a/src/statefbk.py 2011-02-08 22:15:17 UTC (rev 68) +++ branches/control-0.4a/src/statefbk.py 2011-02-08 22:15:21 UTC (rev 69) @@ -265,11 +265,9 @@ raise ValueError, "Oops, the system is unstable!" if type=='c': - print "controllable" trana = 'T' C = -np.dot(sys.B,sys.B.transpose()) elif type=='o': - print "observable" trana = 'N' C = -np.dot(sys.C.transpose(),sys.C) else: This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:15:23
|
Revision: 68 http://python-control.svn.sourceforge.net/python-control/?rev=68&view=rev Author: kkchen Date: 2011-02-08 22:15:17 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Docstring fix for {TransferFunction,StateSpace}.returnScipySignalLti(). Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/statesp.py branches/control-0.4a/src/xferfcn.py Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:15:12 UTC (rev 67) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:15:17 UTC (rev 68) @@ -279,7 +279,7 @@ """Return a list of a list of scipy.signal.lti objects for a MIMO system. For instance, - >>> out = ssobject.returnSignalScipyLti() + >>> out = ssobject.returnScipySignalLti() >>> out[3][5] is a signal.scipy.lti object corresponding to the transfer function from Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:15:12 UTC (rev 67) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:15:17 UTC (rev 68) @@ -391,7 +391,7 @@ """Return a list of a list of scipy.signal.lti objects for a MIMO system. For instance, - >>> out = ssobject.returnSignalScipyLti() + >>> out = tfobject.returnScipySignalLti() >>> out[3][5] is a signal.scipy.lti object corresponding to the transfer function from This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |
From: <kk...@us...> - 2011-02-08 22:15:18
|
Revision: 67 http://python-control.svn.sourceforge.net/python-control/?rev=67&view=rev Author: kkchen Date: 2011-02-08 22:15:12 +0000 (Tue, 08 Feb 2011) Log Message: ----------- Created the "virtual" function returnScipySignalLti() in the TransferFunction and StateSpace classes. Also moved some import commands into member functions (so as not to be global) and edited some docstrings. Kevin K. Chen <kk...@pr...> Modified Paths: -------------- branches/control-0.4a/src/statesp.py branches/control-0.4a/src/xferfcn.py Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:15:06 UTC (rev 66) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:15:12 UTC (rev 67) @@ -275,15 +275,31 @@ return StateSpace(A, B, C, D) -# -# convertToStateSpace - create a state space system from another type -# -# To allow scalar constants to be used in a simple way (k*P, 1+L), this -# function allows the dimension of the input/output system to be specified -# in the case of a scalar system -# + def returnScipySignalLti(self): + """Return a list of a list of scipy.signal.lti objects for a MIMO + system. For instance, + + >>> out = ssobject.returnSignalScipyLti() + >>> out[3][5] + + is a signal.scipy.lti object corresponding to the transfer function from + the 6th input to the 4th output.""" + + from scipy.signal import lti + + # Preallocate the output. + out = [[[] for j in range(self.inputs)] for i in range(self.outputs)] + + for i in range(self.outputs): + for j in range(self.inputs): + out[i][j] = lti(self.A, self.B[:, j], self.C[i, :], + self.D[i, j]) + + return out + def convertToStateSpace(sys, inputs=1, outputs=1): - """Convert a system to state space form (if needed).""" + """Convert a system to state space form (if needed). If sys is a scalar, + then the number of inputs and outputs can be specified manually.""" if isinstance(sys, StateSpace): # Already a state space system; just return it @@ -371,7 +387,7 @@ while True: T = randn(states, states) try: - A = numpy.dot(numpy.linalg.solve(T, A), T) # A = T \ A * T + A = numpy.dot(solve(T, A), T) # A = T \ A * T break except numpy.linalg.linalg.LinAlgError: # In the unlikely event that T is rank-deficient, iterate again. Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:15:06 UTC (rev 66) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:15:12 UTC (rev 67) @@ -48,10 +48,6 @@ # External function declarations import scipy as sp -import scipy.signal as signal -import copy -import bdalg as bd -import statesp from lti2 import Lti2 class TransferFunction(Lti2): @@ -67,7 +63,7 @@ def __init__(self, num=1, den=1): """This is the constructor. The default transfer function is 1 (unit gain direct feedthrough).""" - + # Make num and den into lists of lists of arrays, if necessary. Beware: # this is a shallow copy! This should be okay, but be careful. data = [num, den] @@ -205,6 +201,8 @@ def __neg__(self): """Negate a transfer function.""" + import copy + num = copy.deepcopy(self.num) for i in range(self.outputs): for j in range(self.inputs): @@ -389,6 +387,27 @@ # But this does not work correctly because the state size will be too # large. + def returnScipySignalLti(self): + """Return a list of a list of scipy.signal.lti objects for a MIMO + system. For instance, + + >>> out = ssobject.returnSignalScipyLti() + >>> out[3][5] + + is a signal.scipy.lti object corresponding to the transfer function from + the 6th input to the 4th output.""" + + from scipy.signal import lti + + # Preallocate the output. + out = [[[] for j in range(self.inputs)] for i in range(self.outputs)] + + for i in range(self.outputs): + for j in range(self.inputs): + out[i][j] = lti(self.num[i][j], self.den[i][j]) + + return out + # Utility function to convert a transfer function polynomial to a string # Borrowed from poly1d library def _tfpolyToString(coeffs, var='s'): @@ -448,7 +467,10 @@ return num, den def convertToTransferFunction(sys, inputs=1, outputs=1): - """Convert a system to transfer function form (if needed.)""" + """Convert a system to transfer function form (if needed). If sys is a + scalar, then the number of inputs and outputs can be specified manually.""" + + import statesp if isinstance(sys, TransferFunction): return sys This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |