In [ ]:
import cvxopt as cvx
import control as con
import numpy as np
In [ ]:
Parameters of some Examples from different sources
In [ ]:
##############################
# Boris Paper UBoot #
##############################
A = np.matrix([[0., 1., 0.],
[0., 0., 1.],
[0., 0., -0.005]])
############### TODO ###################
############### CHECK ##################
#last line of Regelungsnormalform (characteristical polynom coefficients) in a colunm vector
a = -np.matrix(A[np.array(A).shape[1]-1])[np.newaxis].T
b = np.matrix([[0], [0], [1.]])
d = 0
c = np.matrix([[1], [0], [0]]).T
print np.matrix(con.matlab.pole(con.matlab.ss(A, b, c, d))).T
u_max = 2.5e-5
X0 = [np.matrix([-10, -0.05, -0.0046]).T,
np.matrix([-10, -0.05, 0.0046]).T,
np.matrix([-10, 0.05, -0.0046]).T,
np.matrix([-10, 0.05, 0.0046]).T,
np.matrix([ 10, -0.05, -0.0046]).T,
np.matrix([ 10, -0.05, 0.0046]).T,
np.matrix([ 10, 0.05, -0.0046]).T,
np.matrix([ 10, 0.05, 0.0046]).T]
x_max = [10, 0.05, 0.0046]
#simulation time
T = 1.0/200.0
#Introduced for fun
u_max_sys = u_max
a_hat = np.matrix([[4.4469e-8], [2.3073e-5], [4.9148e-3]])
a_hat_star = np.matrix([[1.073e-7], [4.919e-5], [10.4078e-3]])
R1 = np.matrix([[1.6021e-5, 3.26098e-3, 0.4031],
[3.2698e-3, 1.5666, 163.46],
[0.4031, 163.46, 40.713]])
##### Entwurf parameter #####
beta = 2 # beta >=1 !
In [ ]:
##############################
# Boris Diss Fusionsreactor #
##############################
A = np.matrix([[1, 0],
[0, -0.5]])
b = np.matrix([[1],[-0.5]])
c = np.matrix([1, 0])
u_max = 1
n = b.size
print "Poles:\n", np.matrix(con.matlab.pole(con.matlab.ss(A, b, c, 0))).T
X0 = [np.matrix([-0.9, -2.8]).T,
np.matrix([-0.9, 2.8]).T,
np.matrix([0.9, -2.8]).T,
np.matrix([0.9, 2.8]).T]
x_max = [0.9, 2.8]
#print "A:\n", A
#print "a:\n", a
#print "b:\n", b
#print "c:\n", c
##### Entwurf parameter #####
beta = 2 # beta >=1 !
In [ ]:
###################################################
# Boris Diss Hydraulisches Positioniersystem S.61 #
###################################################
A = np.matrix([[0, 1, 0],
[0, -0.126, 83.5],
[0, -6.89, -0.00707]])
b = np.matrix([[0],[0],[-55.51]])
c = np.matrix([1, 0, 0])
print "Poles:\n", np.matrix(con.matlab.pole(con.matlab.ss(A, b, c, 0))).T
u_max = 1
n = len(b)
X0 = [np.matrix([-13.38, -10.7, -2.58]).T,
np.matrix([-13.38, -10.7, 2.58]).T,
np.matrix([-13.38, 10.7, -2.58]).T,
np.matrix([-13.38, 10.7, 2.58]).T,
np.matrix([ 13.38, -10.7, -2.58]).T,
np.matrix([ 13.38, -10.7, 2.58]).T,
np.matrix([ 13.38, 10.7, -2.58]).T,
np.matrix([ 13.38, 10.7, 2.58]).T]
x_max = [13.38, 10.7, 2.58]
#print "A:\n", A
#print "a:\n", a
#print "b:\n", b
#print "c:\n", c
p_min = 0.2
mu_star = 1.5
zeta = 2.5
zeta_star = 1.5
k0 = 1e-3*np.matrix([-24.23, -0.2122, -10.5]).T
k1 = 1e-3*np.matrix([-40.13, -0.1396, -15.10]).T
k0_star = 1e-3*np.matrix([-11.96, -0.1111, -5.18]).T
k1_star = 1e-3*np.matrix([-83.35, -0.0477, -33.86]).T
##### Entwurf parameter #####
beta=2 # beta >=1 !
In [ ]:
###########################
# Hydraulischer Aktor #
###########################
A0 = np.matrix([[0, 1, 0],
[-10, -1.167, 25],
[0, 0, -0.8]])
print "Eigenvalues: {}".format(LA.eigvals(A0))
#a = -A[-1,:].T ### !!!!!
#print a
b0 = np.matrix([[0],[0],[2.4]])
c0 = np.matrix([1, 0, 0])
d0 = np.matrix([0])
u_max = 10.5
n = 3
X00 = [np.matrix([-20.0, -10.0, -10.0]).T,
np.matrix([-20.0, -10.0, 10.0]).T,
np.matrix([-20.0, 10.0, -10.0]).T,
np.matrix([20.0, -10.0, 10.0]).T,
np.matrix([20.0, -10.0, -10.0]).T,
np.matrix([20.0, 10.0, 10.0]).T]
#print "A:\n", A
#print "a:\n", a
#print "b:\n", b
#print "c:\n", c
# Convert to Normalform
(A1, b1, c1, d1), T1, Q1 = optim_tools.get_Steuerungsnormalform(A0, b0, c0.T, d0)
a1 = -A1[-1][:].T #!!!!
print "T1:\n", T1
# Convert to Normalform
ss, T2 = con.canonical_form(con.ss(A0, b0, c0, d0), form='reachable')
assert np.allclose(T1*X00[1],
optim_tools.reverse_x_order(T2*X00[1])),\
"own Steuerungsnormalform Transformation not equal python control version"
#print "x_r1:\n", T1*X00[1]
#print "x_r2(backwards):\n", optim_tools.reverse_x_order(T2*X00[1])
A = optim_tools.reverse_x_order(np.matrix(ss.A))
a = -A[-1][:].T #!!!!
b = optim_tools.reverse_x_order(np.matrix(ss.B))
c = optim_tools.reverse_x_order(np.matrix(ss.C))
d = optim_tools.reverse_x_order(np.matrix(ss.D)) # == 0!
print "A:\n", A
assert np.allclose(A, A1)
print "a:\n", a
assert np.allclose(a, a1)
print "b:\n", b
assert np.allclose(b, b1)
print "c:\n", c
assert np.allclose(c, c1.T)
X0 = [T1.dot(x0) for x0 in X00]
print "X0:\n", X0
#print "A1:\n", A1
#print "a1:\n", a1
#print "b1:\n", b1
#print "c1:\n", c1