Skip to content

Add similarity transformation function #319

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
murrayrm merged 2 commits into from
Jun 30, 2019
Merged

Add similarity transformation function #319

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
69ec333
add similarity transformation function (w/ optional time rescaling)
murrayrm Jun 23, 2019
715fe55
update unit test to avoid bugs in modal form
murrayrm Jun 23, 2019
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
39 changes: 37 additions & 2 deletions control/canonical.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -6,10 +6,11 @@
from .statesp import StateSpace
from .statefbk import ctrb, obsv

from numpy import zeros, shape, poly, iscomplex, hstack
from numpy import zeros, shape, poly, iscomplex, hstack, dot, transpose
from numpy.linalg import solve, matrix_rank, eig

__all__ = ['canonical_form', 'reachable_form', 'observable_form']
__all__ = ['canonical_form', 'reachable_form', 'observable_form', 'modal_form',
'similarity_transform']

def canonical_form(xsys, form='reachable'):
"""Convert a system into canonical form
Expand Down Expand Up @@ -212,3 +213,37 @@ def modal_form(xsys):
zsys.C = xsys.C.dot(Tzx)

return zsys, Tzx


def similarity_transform(xsys, T, timescale=1):
"""Perform a similarity transformation, with option time rescaling.

Transform a linear state space system to a new state space representation
z = T x, where T is an invertible matrix.

Parameters
----------
T : 2D invertible array
The matrix `T` defines the new set of coordinates z = T x.
timescale : float
If present, also rescale the time unit to tau = timescale * t

Returns
-------
zsys : StateSpace object
System in transformed coordinates, with state 'z'

"""
# Create a new system, starting with a copy of the old one
zsys = StateSpace(xsys)

# Define a function to compute the right inverse (solve x M = y)
def rsolve(M, y):
return transpose(solve(transpose(M), transpose(y)))

# Update the system matrices
zsys.A = rsolve(T, dot(T, zsys.A)) / timescale
zsys.B = dot(T, zsys.B) / timescale
zsys.C = rsolve(T, zsys.C)

return zsys
74 changes: 72 additions & 2 deletions control/tests/canonical_test.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,9 +2,9 @@

import unittest
import numpy as np
from control import ss, tf
from control import ss, tf, tf2ss, ss2tf
from control.canonical import canonical_form, reachable_form, \
observable_form, modal_form
observable_form, modal_form, similarity_transform
from control.exception import ControlNotImplemented

class TestCanonical(unittest.TestCase):
Expand Down Expand Up @@ -218,6 +218,76 @@ def test_arguments(self):
np.testing.assert_raises(
ControlNotImplemented, canonical_form, sys, 'unknown')

def test_similarity(self):
"""Test similarty transform"""

# Single input, single output systems
siso_ini = tf2ss(tf([1, 1], [1, 1, 1]))
for form in 'reachable', 'observable':
# Convert the system to one of the canonical forms
siso_can, T_can = canonical_form(siso_ini, form)

# Use a similarity transformation to transform it back
siso_sim = similarity_transform(siso_can, np.linalg.inv(T_can))

# Make sure everything goes back to the original form
np.testing.assert_array_almost_equal(siso_sim.A, siso_ini.A)
np.testing.assert_array_almost_equal(siso_sim.B, siso_ini.B)
np.testing.assert_array_almost_equal(siso_sim.C, siso_ini.C)
np.testing.assert_array_almost_equal(siso_sim.D, siso_ini.D)

# Multi-input, multi-output systems
mimo_ini = ss(
[[-1, 1, 0, 0], [0, -2, 1, 0], [0, 0, -3, 1], [0, 0, 0, -4]],
[[1, 0], [0, 0], [0, 1], [1, 1]],
[[1, 0, 0, 0], [0, 1, 0, 0], [0, 0, 1, 0]],
np.zeros((3, 2)))

# Simple transformation: row/col flips + scaling
mimo_txf = np.array(
[[0, 1, 0, 0], [2, 0, 0, 0], [0, 0, -1, 0], [0, 0, 0, 1]])

# Transform the system and transform it back
mimo_sim = similarity_transform(mimo_ini, mimo_txf)
mimo_new = similarity_transform(mimo_sim, np.linalg.inv(mimo_txf))
np.testing.assert_array_almost_equal(mimo_new.A, mimo_ini.A)
np.testing.assert_array_almost_equal(mimo_new.B, mimo_ini.B)
np.testing.assert_array_almost_equal(mimo_new.C, mimo_ini.C)
np.testing.assert_array_almost_equal(mimo_new.D, mimo_ini.D)

# Make sure rescaling by identify does nothing
mimo_new = similarity_transform(mimo_ini, np.eye(4))
np.testing.assert_array_almost_equal(mimo_new.A, mimo_ini.A)
np.testing.assert_array_almost_equal(mimo_new.B, mimo_ini.B)
np.testing.assert_array_almost_equal(mimo_new.C, mimo_ini.C)
np.testing.assert_array_almost_equal(mimo_new.D, mimo_ini.D)

# Time rescaling
mimo_tim = similarity_transform(mimo_ini, np.eye(4), timescale=0.3)
mimo_new = similarity_transform(mimo_tim, np.eye(4), timescale=1/0.3)
np.testing.assert_array_almost_equal(mimo_new.A, mimo_ini.A)
np.testing.assert_array_almost_equal(mimo_new.B, mimo_ini.B)
np.testing.assert_array_almost_equal(mimo_new.C, mimo_ini.C)
np.testing.assert_array_almost_equal(mimo_new.D, mimo_ini.D)

# Time + transformation, in one step
mimo_sim = similarity_transform(mimo_ini, mimo_txf, timescale=0.3)
mimo_new = similarity_transform(mimo_sim, np.linalg.inv(mimo_txf),
timescale=1/0.3)
np.testing.assert_array_almost_equal(mimo_new.A, mimo_ini.A)
np.testing.assert_array_almost_equal(mimo_new.B, mimo_ini.B)
np.testing.assert_array_almost_equal(mimo_new.C, mimo_ini.C)
np.testing.assert_array_almost_equal(mimo_new.D, mimo_ini.D)

# Time + transformation, in two steps
mimo_sim = similarity_transform(mimo_ini, mimo_txf, timescale=0.3)
mimo_tim = similarity_transform(mimo_sim, np.eye(4), timescale=1/0.3)
mimo_new = similarity_transform(mimo_tim, np.linalg.inv(mimo_txf))
np.testing.assert_array_almost_equal(mimo_new.A, mimo_ini.A)
np.testing.assert_array_almost_equal(mimo_new.B, mimo_ini.B)
np.testing.assert_array_almost_equal(mimo_new.C, mimo_ini.C)
np.testing.assert_array_almost_equal(mimo_new.D, mimo_ini.D)

def suite():
return unittest.TestLoader().loadTestsFromTestCase(TestFeedback)

Expand Down