Skip to content

Commit 3738718

Browse files
murrayrmmurrayrm
authored
Merge pull request #245 from hungpham2511/master
Implemement lft (Linear Fractional Transform) for state-space systems
2 parents 5d680f5 + 5955fd2 commit 3738718

File tree

2 files changed

+148
-0
lines changed

2 files changed

+148
-0
lines changed

‎control/statesp.py‎

Lines changed: 101 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -602,6 +602,107 @@ def feedback(self, other=1, sign=-1):
602602

603603
return StateSpace(A, B, C, D, dt)
604604

605+
def lft(self, other, nu=-1, ny=-1):
606+
"""Return the Linear Fractional Transformation.
607+
608+
A definition of the LFT operator can be found in Appendix A.7,
609+
page 512 in the 2nd Edition, Multivariable Feedback Control by
610+
Sigurd Skogestad.
611+
612+
An alternative definition can be found here:
613+
https://www.mathworks.com/help/control/ref/lft.html
614+
615+
Parameters
616+
----------
617+
other: LTI
618+
The lower LTI system
619+
ny: int, optional
620+
Dimension of (plant) measurement output.
621+
nu: int, optional
622+
Dimension of (plant) control input.
623+
624+
"""
625+
other = _convertToStateSpace(other)
626+
# maximal values for nu, ny
627+
if ny == -1:
628+
ny = min(other.inputs, self.outputs)
629+
if nu == -1:
630+
nu = min(other.outputs, self.inputs)
631+
# dimension check
632+
# TODO
633+
634+
# Figure out the sampling time to use
635+
if (self.dt == None and other.dt != None):
636+
dt = other.dt # use dt from second argument
637+
elif (other.dt == None and self.dt != None) or \
638+
timebaseEqual(self, other):
639+
dt = self.dt # use dt from first argument
640+
else:
641+
raise ValueError("Systems have different time bases")
642+
643+
# submatrices
644+
A = self.A
645+
B1 = self.B[:, :self.inputs - nu]
646+
B2 = self.B[:, self.inputs - nu:]
647+
C1 = self.C[:self.outputs - ny, :]
648+
C2 = self.C[self.outputs - ny:, :]
649+
D11 = self.D[:self.outputs - ny, :self.inputs - nu]
650+
D12 = self.D[:self.outputs - ny, self.inputs - nu:]
651+
D21 = self.D[self.outputs - ny:, :self.inputs - nu]
652+
D22 = self.D[self.outputs - ny:, self.inputs - nu:]
653+
654+
# submatrices
655+
Abar = other.A
656+
Bbar1 = other.B[:, :ny]
657+
Bbar2 = other.B[:, ny:]
658+
Cbar1 = other.C[:nu, :]
659+
Cbar2 = other.C[nu:, :]
660+
Dbar11 = other.D[:nu, :ny]
661+
Dbar12 = other.D[:nu, ny:]
662+
Dbar21 = other.D[nu:, :ny]
663+
Dbar22 = other.D[nu:, ny:]
664+
665+
# well-posed check
666+
F = np.block([[np.eye(ny), -D22], [-Dbar11, np.eye(nu)]])
667+
if matrix_rank(F) != ny + nu:
668+
raise ValueError("lft not well-posed to working precision.")
669+
670+
# solve for the resulting ss by solving for [y, u] using [x,
671+
# xbar] and [w1, w2].
672+
TH = np.linalg.solve(F, np.block(
673+
[[C2, np.zeros((ny, other.states)), D21, np.zeros((ny, other.inputs - ny))],
674+
[np.zeros((nu, self.states)), Cbar1, np.zeros((nu, self.inputs - nu)), Dbar12]]
675+
))
676+
T11 = TH[:ny, :self.states]
677+
T12 = TH[:ny, self.states: self.states + other.states]
678+
T21 = TH[ny:, :self.states]
679+
T22 = TH[ny:, self.states: self.states + other.states]
680+
H11 = TH[:ny, self.states + other.states: self.states + other.states + self.inputs - nu]
681+
H12 = TH[:ny, self.states + other.states + self.inputs - nu:]
682+
H21 = TH[ny:, self.states + other.states: self.states + other.states + self.inputs - nu]
683+
H22 = TH[ny:, self.states + other.states + self.inputs - nu:]
684+
685+
Ares = np.block([
686+
[A + B2.dot(T21), B2.dot(T22)],
687+
[Bbar1.dot(T11), Abar + Bbar1.dot(T12)]
688+
])
689+
690+
Bres = np.block([
691+
[B1 + B2.dot(H21), B2.dot(H22)],
692+
[Bbar1.dot(H11), Bbar2 + Bbar1.dot(H12)]
693+
])
694+
695+
Cres = np.block([
696+
[C1 + D12.dot(T21), D12.dot(T22)],
697+
[Dbar21.dot(T11), Cbar2 + Dbar21.dot(T12)]
698+
])
699+
700+
Dres = np.block([
701+
[D11 + D12.dot(H21), D12.dot(H22)],
702+
[Dbar21.dot(H11), Dbar22 + Dbar21.dot(H12)]
703+
])
704+
return StateSpace(Ares, Bres, Cres, Dres, dt)
705+
605706
def minreal(self, tol=0.0):
606707
"""Calculate a minimal realization, removes unobservable and
607708
uncontrollable states"""

‎control/tests/statesp_test.py‎

Lines changed: 47 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -448,6 +448,53 @@ def empty(shape):
448448
np.testing.assert_array_equal(empty((D.shape[0], 0)), g.C)
449449
np.testing.assert_array_equal(D, g.D)
450450

451+
def test_lft(self):
452+
""" test lft function with result obtained from matlab implementation"""
453+
# test case
454+
A = [[1, 2, 3],
455+
[1, 4, 5],
456+
[2, 3, 4]]
457+
B = [[0, 2],
458+
[5, 6],
459+
[5, 2]]
460+
C = [[1, 4, 5],
461+
[2, 3, 0]]
462+
D = [[0, 0],
463+
[3, 0]]
464+
P = StateSpace(A, B, C, D)
465+
Ak = [[0, 2, 3],
466+
[2, 3, 5],
467+
[2, 1, 9]]
468+
Bk = [[1, 1],
469+
[2, 3],
470+
[9, 4]]
471+
Ck = [[1, 4, 5],
472+
[2, 3, 6]]
473+
Dk = [[0, 2],
474+
[0, 0]]
475+
K = StateSpace(Ak, Bk, Ck, Dk)
476+
477+
# case 1
478+
pk = P.lft(K, 2, 1)
479+
Amatlab = [1, 2, 3, 4, 6, 12, 1, 4, 5, 17, 38, 61, 2, 3, 4, 9, 26, 37, 2, 3, 0, 3, 14, 18, 4, 6, 0, 8, 27, 35, 18, 27, 0, 29, 109, 144]
480+
Bmatlab = [0, 10, 10, 7, 15, 58]
481+
Cmatlab = [1, 4, 5, 0, 0, 0]
482+
Dmatlab = [0]
483+
np.testing.assert_allclose(np.array(pk.A).reshape(-1), Amatlab)
484+
np.testing.assert_allclose(np.array(pk.B).reshape(-1), Bmatlab)
485+
np.testing.assert_allclose(np.array(pk.C).reshape(-1), Cmatlab)
486+
np.testing.assert_allclose(np.array(pk.D).reshape(-1), Dmatlab)
487+
488+
# case 2
489+
pk = P.lft(K)
490+
Amatlab = [1, 2, 3, 4, 6, 12, -3, -2, 5, 11, 14, 31, -2, -3, 4, 3, 2, 7, 0.6, 3.4, 5, -0.6, -0.4, 0, 0.8, 6.2, 10, 0.2, -4.2, -4, 7.4, 33.6, 45, -0.4, -8.6, -3]
491+
Bmatlab = []
492+
Cmatlab = []
493+
Dmatlab = []
494+
np.testing.assert_allclose(np.array(pk.A).reshape(-1), Amatlab)
495+
np.testing.assert_allclose(np.array(pk.B).reshape(-1), Bmatlab)
496+
np.testing.assert_allclose(np.array(pk.C).reshape(-1), Cmatlab)
497+
np.testing.assert_allclose(np.array(pk.D).reshape(-1), Dmatlab)
451498

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

0 commit comments

Comments
 (0)