Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
def weighting(wb, m, a):
"""weighting(wb,m,a) -> wf
wb - design frequency (where |wf| is approximately 1)
m - high frequency gain of 1/wf; should be > 1
a - low frequency gain of 1/wf; should be < 1
wf - SISO LTI object
"""
s = tf([1, 0], [1])
return (s/m + wb) / (s + wb*a)
:param G: transfer function
:param K_guess: gain matrix guess
:param d_tc: time constant for derivative
:param verbose: show debug output
:param use_P: use p gain in design
:param use_I: use i gain in design
:param use_D: use d gain in design
:return: (K, G_comp, Gc_comp)
K: gain matrix
G_comp: open loop compensated plant
Gc_comp: closed loop compensated plant
"""
# compensator transfer function
H = []
if use_P:
H += [control.tf(1, 1)]
if use_I:
H += [control.tf((1), (1, 0))]
if use_D:
H += [control.tf((1, 0), (d_tc, 1))]
H = np.array([H]).T
H_num = [[H[i][j].num[0][0] for i in range(H.shape[0])] for j in range(H.shape[1])]
H_den = [[H[i][j].den[0][0] for i in range(H.shape[0])] for j in range(H.shape[1])]
H = control.tf(H_num, H_den)
# print('G', G)
# print('H', H)
ss_open = control.tf2ss(G * H)
if verbose:
print('optimizing controller')
:param h:
:param theta:
:param kp:
:param kd:
:param kdd:
:param tau:
:return:
"""
Hinv = tf([1], [h, 1])
K = tf([kdd, kd, kp], [1])
G = tf([1], [tau, 1, 0, 0])
KG = K * G
num_delay, den_delay = pade(theta, n=order)
D = tf(num_delay, den_delay)
one = tf([1], [1])
num_KG, den_KG = tfdata(KG + one)
invKG1 = tf(den_KG, num_KG)
sys = invKG1 * (KG + D) * Hinv
return sys
def plant():
"""plant() -> g
g - LTI object; 2x2 plant with a RHP zero, at s=0.5.
"""
den = [0.2, 1.2, 1]
gtf = tf([[[1], [1]],
[[2, 1], [2]]],
[[den, den],
[den, den]])
return ss(gtf)
g_tf = control.tf([12.8], [16.7,1])
td = 1
g_delay = control.tf(control.pade(td,1)[0],control.pade(td,1)[1])
g11 = g_tf*g_delay
g_tf = control.tf([-18.9], [21.0,1])
td = 3
g_delay = control.tf(control.pade(td,1)[0],control.pade(td,1)[1])
g12 = g_tf*g_delay
g_tf = control.tf([6.6], [10.9,1])
td = 7
g_delay = control.tf(control.pade(td,1)[0],control.pade(td,1)[1])
g21 = g_tf*g_delay
g_tf = control.tf([-19.4], [14.4,1])
td = 3
g_delay = control.tf(control.pade(td,1)[0],control.pade(td,1)[1])
g22 = g_tf*g_delay
g_tf = control.tf([3.8], [14.9,1])
td = 8
g_delay = control.tf(control.pade(td,1)[0],control.pade(td,1)[1])
g1f = g_tf*g_delay
g_tf = control.tf([4.9], [13.2,1])
td = 3
g_delay = control.tf(control.pade(td,1)[0],control.pade(td,1)[1])
g2f = g_tf*g_delay
row_1_num = [x[0][0] for x in (g11.num, g12.num, g1f.num)]
row_2_num = [x[0][0] for x in (g21.num, g22.num, g2f.num)]
def construct_system_block(h=0.0, theta=0.0, kp=0.0, kd=0.0, kdd=0.0, tau=0.0):
"""
See Ploeg2014 for the notation and variable naming used here. "inv"
indicates inverse.
:param h:
:param theta:
:param kp:
:param kd:
:param kdd:
:param tau:
:return:
"""
Hinv = tf([1], [h, 1])
K = tf([kdd, kd, kp], [1])
G = tf([1], [tau, 1, 0, 0])
KG = K * G
num_delay, den_delay = pade(theta, n=order)
D = tf(num_delay, den_delay)
one = tf([1], [1])
num_KG, den_KG = tfdata(KG + one)
invKG1 = tf(den_KG, num_KG)
sys = invKG1 * (KG + D) * Hinv
return sys
:param verbose: show debug output
:return: (G_ol, delay, k)
G_ol: open loop plant
delay: time delay (sec)
k: gain
"""
if verbose:
print('solving for plant model ', end='')
k, delay = delay_and_gain_sysid(y_acc, u_mix, verbose)
if verbose:
print('done')
# open loop, rate output, mix input plant
tf_acc = k * control.tf(*control.pade(delay, 1)) # order 1 approx
# can neglect motor, pole too fast to matter compared to delay, not used in sysid either
tf_integrator = control.tf((1), (1, 0))
G_ol = tf_acc * tf_integrator
return G_ol, delay, k
self.max_iterations)
if max_reached is True:
warnings.warn("[ARMAX ID] Max reached for:na: {} | nb: {} | nc: {} | "
"delay: {}".format(na, nb, nc, delay))
IC = information_criterion(na + nb + delay,
y.size - max(na, nb + delay, nc),
Vn, self.method)
if IC < IC_old:
self.na, self.nb, self.nc, self.delay, IC_old = na, nb, nc, delay, IC
G_num_opt, G_den_opt, H_num_opt, H_den_opt = G_num, G_den, H_num, H_den
self.Vn, self.Yid, self.max_reached = Vn, np.atleast_2d(y_id) * y_std, max_reached
G_num_opt[self.delay:self.nb + self.delay] = \
G_num_opt[self.delay:self.nb + self.delay] * y_std / u_std
self.G = cnt.tf(G_num_opt, G_den_opt, self.dt)
self.H = cnt.tf(H_num_opt, H_den_opt, self.dt)
print(self)
td = 8
g_delay = control.tf(control.pade(td,1)[0],control.pade(td,1)[1])
g1f = g_tf*g_delay
g_tf = control.tf([4.9], [13.2,1])
td = 3
g_delay = control.tf(control.pade(td,1)[0],control.pade(td,1)[1])
g2f = g_tf*g_delay
row_1_num = [x[0][0] for x in (g11.num, g12.num, g1f.num)]
row_2_num = [x[0][0] for x in (g21.num, g22.num, g2f.num)]
row_1_den = [x[0][0] for x in (g11.den, g12.den, g1f.den)]
row_2_den = [x[0][0] for x in (g21.den, g22.den, g2f.den)]
sys = control.tf([row_1_num,row_2_num],[row_1_den,row_2_den])
R = Tag('Reflux', 'Reflux flow rate', IOtype='INPUT')
S = Tag('Steam', 'Steam flow rate', IOtype='INPUT')
F = Tag('Feed', 'Feed flow rate', IOtype='INPUT')
xD = Tag('x_D', 'Distillate purity', IOtype='OUTPUT')
xB = Tag('x_B', 'Bottoms purity', IOtype='OUTPUT')
inputs = {'R':R,'S':S,'F':F}
outputs = {'xD':xD,'xB':xB}
wood_berry = Model(sys, "Wood Berry Distillation Model", inputs, outputs)
return wood_berry