Menu
▾
▴
[python-control-discuss] SF.net SVN: python-control:[105] branches/control-0.4a/src
|
From: <kk...@us...> - 2011-02-08 22:18:29
|
Revision: 105
http://python-control.svn.sourceforge.net/python-control/?rev=105&view=rev
Author: kkchen
Date: 2011-02-08 22:18:23 +0000 (Tue, 08 Feb 2011)
Log Message:
-----------
Modified _common_den to handle complex roots more accurately.
Minor changes in TestConvert.py and TestBDAlg.py.
Kevin K. Chen <kk...@pr...>
Modified Paths:
--------------
branches/control-0.4a/src/TestBDAlg.py
branches/control-0.4a/src/TestConvert.py
branches/control-0.4a/src/xferfcn.py
Modified: branches/control-0.4a/src/TestBDAlg.py
===================================================================
--- branches/control-0.4a/src/TestBDAlg.py 2011-02-08 22:18:19 UTC (rev 104)
+++ branches/control-0.4a/src/TestBDAlg.py 2011-02-08 22:18:23 UTC (rev 105)
@@ -38,7 +38,6 @@
ans1 = feedback(self.x1, self.sys2)
ans2 = feedback(self.x1, self.sys2, 1.)
- # This one doesn't work yet. The feedback system has too many states.
np.testing.assert_array_almost_equal(ans1.A, [[-1.5, 4.], [13., 2.]])
np.testing.assert_array_almost_equal(ans1.B, [[2.5], [-10.]])
np.testing.assert_array_almost_equal(ans1.C, [[-2.5, 0.]])
@@ -65,7 +64,6 @@
ans1 = feedback(self.sys2, self.x1)
ans2 = feedback(self.sys2, self.x1, 1.)
- # This one doesn't work yet. The feedback system has too many states.
np.testing.assert_array_almost_equal(ans1.A, [[-1.5, 4.], [13., 2.]])
np.testing.assert_array_almost_equal(ans1.B, [[1.], [-4.]])
np.testing.assert_array_almost_equal(ans1.C, [[1., 0.]])
Modified: branches/control-0.4a/src/TestConvert.py
===================================================================
--- branches/control-0.4a/src/TestConvert.py 2011-02-08 22:18:19 UTC (rev 104)
+++ branches/control-0.4a/src/TestConvert.py 2011-02-08 22:18:23 UTC (rev 105)
@@ -4,11 +4,14 @@
Test state space and transfer function conversion.
-Currently, this unit test script is not complete. It converts
-several random state spaces back and forth between state space and transfer
-function representations, and asserts that the conversion outputs are correct.
-As it stands, tests pass but there is some rounding error in one of the conversions leading to imaginary numbers. See the warning message. This script may be used to diagnose errors.
+Currently, this unit test script is not complete. It converts several random
+state spaces back and forth between state space and transfer function
+representations. Ideally, it should be able to assert that the conversion
+outputs are correct. This is not yet implemented.
+Also, the conversion seems to enter an infinite loop once in a while. The cause
+of this is unknown.
+
"""
import numpy as np
@@ -22,11 +25,11 @@
"""Set up testing parameters."""
# Number of times to run each of the randomized tests.
- self.numTests = 1
+ self.numTests = 10
# Maximum number of states to test + 1
- self.maxStates = 3
+ self.maxStates = 2
# Maximum number of inputs and outputs to test + 1
- self.maxIO = 3
+ self.maxIO = 2
# Set to True to print systems to the output.
self.debug = True
Modified: branches/control-0.4a/src/xferfcn.py
===================================================================
--- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:18:19 UTC (rev 104)
+++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:18:23 UTC (rev 105)
@@ -75,8 +75,8 @@
"""
# External function declarations
-from numpy import angle, array, empty, finfo, insert, ndarray, ones, polyadd, \
- polymul, polyval, roots, sort, zeros
+from numpy import angle, any, array, empty, finfo, insert, ndarray, ones, \
+ polyadd, polymul, polyval, roots, sort, sqrt, zeros
from scipy.signal import lti
from copy import deepcopy
from slycot import tb04ad
@@ -491,7 +491,10 @@
[something] array.
"""
-
+
+ # Machine precision for floats.
+ eps = finfo(float).eps
+
# A sorted list to keep track of cumulative poles found as we scan
# self.den.
poles = []
@@ -512,8 +515,7 @@
# 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 abs(currentpoles[cp_ind] - poles[p_ind]) < (10 * eps):
# If the current element of both lists match, then we're
# good. Move to the next pair of elements.
cp_ind += 1
@@ -558,9 +560,22 @@
# Construct the common denominator.
den = 1.
- for p in poles:
- den = polymul(den, [1., -p])
+ n = 0
+ while n < len(poles):
+ if abs(poles[n].imag) > 10 * eps:
+ # To prevent buildup of imaginary part error, handle complex
+ # pole pairs together.
+ quad = polymul([1., -poles[n]], [1., -poles[n+1]])
+ assert all(quad.imag < eps), "The quadratic has a nontrivial \
+imaginary part: %g" % quad.imag.max()
+ quad = quad.real
+ den = polymul(den, quad)
+ n += 2
+ else:
+ den = polymul(den, [1., -poles[n].real])
+ n += 1
+
# Modify the numerators so that they each take the common denominator.
num = deepcopy(self.num)
for i in range(self.outputs):
@@ -585,6 +600,11 @@
zeros(largest - len(num[i][j])))
# Finally, convert the numerator to a 3-D array.
num = array(num)
+ # Remove trivial imaginary parts. Check for nontrivial imaginary parts.
+ if any(abs(num.imag) > sqrt(eps)):
+ print ("Warning: The numerator has a nontrivial nontrivial part: %g"
+ % abs(num.imag).max())
+ num = num.real
return num, den
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|
