我习惯于在MatLab中进行离散时间控制系统的仿真,现在我正在尝试使用Python和NumPy。
因此,我下面的代码是有效的,但我想迭代数字向量,而不是将值追加到列表中。有可能吗?
换句话说,与其使用
xl.append(xt)
ul.append(uc)
我想使用一些与MatLab等价的代码,比如x[:,k+1]=np.dot(Ad,x[:,k])+Bd*uc,但它在我的代码上不起作用。如果我这样做,我没有得到预期的两行列向量,而是得到了一个2x2矩阵和一个错误。
另一个问题:为什么必须使用plt.lot(tk,u[:,0],Label=‘u’)而不是plt.lot(tk,u,Label=‘u’)?
from control.matlab import *
import math
import numpy as np
import matplotlib.pyplot as plt
Ts = 0.1
N = 50
#x = np.zeros((2, N+1))
tk = np.zeros(N)
u = np.zeros(N)
v = np.random.randn(N)/86.6 #% measurement noise
wn = 1.12
wn2 = pow(wn, 2)
A = [[0, 1], [-1.5, -1.4]]
B = [[0], [1.5]]
C = [[1, 0]]
D = 0
# Control gains
K = np.array([2.64, 3.41071429])
# Now build a feedback with control law u = -K*x
Ad = np.eye(2) + np.multiply(A, Ts)
Bd = np.multiply(B, Ts)
Cd = C
xt = [[1.0], [0.12]] # initial states
xl = []
ul = []
for k in range(0, N):
tk[k] = k*Ts
uc = -K.dot(xt)
xt = np.dot(Ad, xt) + Bd*uc
xt[1, 0] += v[k]
xl.append(xt)
ul.append(uc)
x = np.array(xl)
u = np.array(ul)
#x = np.delete(x, N, 1) # delete the last position of x
#s = TransferFunction.s
#Gs = wn2/(s**2 + 0*s + wn2) # This is the KF solution
#yout, T = step(Gs)
plt.rcParams["figure.figsize"] = (10, 7)
plt.figure()
#plt.plot(T, yout, label='Open loop')
plt.plot(tk, x[:, 0], label='x_0')
plt.plot(tk, x[:, 1], label='x_1')
plt.plot(tk, u[:, 0], label='u')
plt.legend()
plt.title('Pendulum ex. 7.14 Franklin book')
plt.xlabel('Time')
plt.ylabel('amp.')
plt.show()
我想要的是这样的代码:
from control.matlab import *
import math
import numpy as np
import matplotlib.pyplot as plt
Ts = 0.1
N = 50
x = np.zeros((2, N+1))
tk = np.zeros(N)
u = np.zeros(N)
v = np.random.randn(N)/86.6 #% measurement noise
wn = 1.12
wn2 = pow(wn, 2)
A = [[0, 1], [-1.5, -1.4]]
B = [[0], [1.5]]
C = [[1, 0]]
D = 0
# Control gains
K = np.array([2.64, 3.41071429])
# Now build a feedback with control law u = -K*x
Ad = np.eye(2) + np.multiply(A, Ts)
Bd = np.multiply(B, Ts)
Cd = C
for k in range(0, N):
tk[k] = k*Ts
u[k] = -K.dot(x[:, k])
x[1, k] += v[k]
x[:, k+1] = np.dot(Ad, x[:, k]) + Bd*u[k]
x = np.delete(x, N, 1) # delete the last position of x
#s = TransferFunction.s
#Gs = wn2/(s**2 + 0*s + wn2) # This is the KF solution
#yout, T = step(Gs)
plt.rcParams["figure.figsize"] = (10, 7)
plt.figure()
#plt.plot(T, yout, label='Open loop')
plt.plot(tk, x[:, 0], label='x_0')
plt.plot(tk, x[:, 1], label='x_1')
plt.plot(tk, u[:, 0], label='u')
plt.legend()
plt.title('Pendulum ex. 7.14 Franklin book')
plt.xlabel('Time')
plt.ylabel('amp.')
plt.show()
但它会导致以下错误:
Traceback (most recent call last):
File "C:\Users\ ... \np_matrices_v1.py", line 46, in <module>
x[:, k+1] = np.dot(Ad, x[:, k]) + Bd*u[k]
ValueError: could not broadcast input array from shape (2,2) into shape (2,)
2条答案
按热度按时间qqrboqgw1#
我不知道为什么,但如果你试着:
您将获得:
但如果你这么做了:
你会得到正确的答案:
有谁能解释一下吗?
ymdaylpp2#
注意它们的形状。
x
是(2,1),因此y
也是。y
可以分配给(2,1)插槽,但不能分配给(2,)形状。与MatLab相比,
numpy
阵列可以是一维的,例如(2,)。此外,领先维度是最外层的,而不是拖后性维度。MatLab很容易将(2,3,1)形简化为(2,3)形,但(2,1,1)形只变成(2,1)形。broadcasting
与numpy
使用形状可能不同的数组相同。两条基本规则是因此,a(2,)可以变成a(1,2)。
如果从
x
中删除内部的[],则会得到一个一维数组:然后可以将其分配给
xa
:xa[:,2] = A.dot(x)
行reshape
和ravel
可用于删除尺寸。也索引A.dot(x)[:,0]