我试图解决以下两个问题:
其思想是求解q,然后绘制i(t)=dq/dt。以下是完整的代码:
import timeit
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint
start = timeit.default_timer()
#Data:
t = np.linspace(0, 10, 100)
w1 = 350 #350 rad/s
w2 = 1e3 #1000 rad/s
phi1 = 150
phi2 = 30
E = np.cos(t)
R1 = 10e3 #10 kOhm
R2 = 3.3e3 #3k3 Ohm
L = 10e-3 #10 mHy
C = 1.56e-6 #1.56 uF
#Create model -> x = [q , k]
def model(x , t):
q = x[0]
k = x[1]
dq_dt = x[1]
i = dq_dt
R = R1 * i + R2 * i**3
dk_dt = 1 / L * (E - R * k - 1/C * q)
dx_dt = np.array([dq_dt , dk_dt])
return dx_dt
#Init. cond.:
x0 = np.array([0 , 0])
#Solver ODE:
x = odeint(model, x0, t)
q = x[: , 0]
i = np.diff(q)
vr = (R1 * i + R2 * i**3) * i
vl = L * np.diff(i)
vc = 1/C * q
#Plot:
plt.style.use("bmh")
fig,(ax1,ax2) = plt.subplots(nrows = 2, ncols = 1, sharex = True)
ax1.plot(t, E, "r-", linewidth = 2)
ax2.plot(t, np.append(i,0), "b-", linewidth = 2)
ax1.grid(True)
ax1.set_title("E(t)")
ax1.set_ylabel("E(t) [V]")
ax2.grid(True)
ax2.set_title("i(t)")
ax2.set_ylabel("i(t) [A]")
ax2.set_xlabel("t [s]")
plt.tight_layout()
plt.show()
end = timeit.default_timer()
print ("\n\n\nSIM TIME = ", end - start, " s")
当我像那样运行代码时,会出现以下错误:
File "...", line 73, in <module>
x = odeint(model, x0, t)
File ".../python3.7/site-packages/scipy/integrate/odepack.py", line 244, in odeint
int(bool(tfirst)))
TypeError: Cannot cast array data from dtype('O') to dtype('float64') according to the rule 'safe'
问题在于dx_dt = np.array([dq_dt , dk_dt])
的返回。由于E
是100个元素的数组,dk_dt = dx_dt[1]
也将是100个元素的数组,dq_dt = dx_dt[0]
只是一个值。如果我这样做dx_dt = np.array([dq_dt , dk_dt] , dtype = float)
,我会得到以下错误:
TypeError: only size-1 arrays can be converted to Python scalars
The above exception was the direct cause of the following exception:
Traceback (most recent call last):
File "...", line 73, in <module>
x = odeint(model, x0, t)
File ".../python3.7/site-packages/scipy/integrate/odepack.py", line 244, in odeint
int(bool(tfirst)))
File "...", line 65, in model
dx_dt = np.array([dq_dt , dk_dt] , dtype = float)
ValueError: setting an array element with a sequence.
那么,有没有一种方法可以将dx_dt
转换为浮点,这样它就不会打扰odeint
当我通过
x0
定义复制-n-paste时:并在起点进行测试计算:
model
返回包含标量0和(100,)数组的2元素对象数组。那不是odeint
可以处理的相关问题 更多 >
编程相关推荐