值错误和odepack.error使用集成.odeint()

Question

我试着写一个方程来建模，然后画出一个完整的控制系统（特别是关于巡航控制）。但是，每次运行时都会收到两个错误：

ValueError:对象太深，无法使用所需的数组 odepack.error：函数调用的结果不是正确的浮点数组。在

我读过这些问题：

• scipy curve_fit error: Result from function call is not a proper array of floats
• How to solve this differential equation using scipy odeint?
• Object Too Deep for Desired Array - scipy.integrate.odeint

这看起来应该是有帮助的，但是我不确定如何应用这些来解决我的问题。我对python相当陌生，所以如果我错过了一些显而易见的事情或者做了一些非常愚蠢的事情，请容忍我。我对绘制它没有任何问题，所以一旦我找到了如何让它工作起来，我想我已经准备好了。在

import numpy as np
import scipy.integrate as integrate

##Parameters

kp=.5 #proportional gain
ki=.1 #integral gain
vr=30 #desired velocity in m/s
Tm=190 #Max Torque in Nm
wm=420 #engine speed
B=0.4 #Beta
an=12 #at gear 4
p=1.3 #air density
Cd=0.32 #Drag coefficient
Cr=.01 #Coefficient of rolling friction
A=2.4 #frontal area

##Variables

m=18000 #weight
v=20 #starting velocity
time=np.linspace(0,10,50) #time
theta=np.radians(4) #Theta

def vderivs(state,t):
    v = state
    vel=[]
    ti=0

    while ti < t:
        v1 = an*controller(ti,vr,v)*torque(v)
        v2 = m*Cr*np.sign(v)
        v3 = 0.5*p*Cd*A*v**2
        v4 = m*np.sin(theta)

        if t < 10:
            vtot = v1+v2+v3
            vfin = np.divide(vtot,m)
        else:
            vtot = v1+v2+v3+v4
            vfin = np.divide(vtot,m)

        vel.append(vfin)
        ti+=1

    trueVel = np.array(vel, float)
    return trueVel

def uderivs(state,t):
    v = state
    deltax = vr - v
    return deltax    

def controller(time,desired,currentV):
    z = integrate.odeint(uderivs, currentV, time)
    u = kp*(vr-currentV)+ki*z
    return u.flatten()

def torque(v):    
    return Tm*(1-B*(np.divide(an*v,wm)-1)**2)   

def velocity(mass,desired,theta,t):
    v = integrate.odeint(vderivs, desired, t)
    return v.flatten()

test = velocity(m,vr,theta,time)
print(test)

如果你还需要我帮忙，请告诉我！在

Answer 1

我看到的一个错误是函数controller；您试图从一个整数（vr）中减去一个函数（velocity），这可能是某些问题的根源。在

“object too deep for desired array”错误可能来自另一个问题，其中函数controller和velocity的返回形状错误；它们必须是一维numpy数组。您可以使用flatten()方法修复此问题：

In [82]: z.shape
Out[82]: (50, 1)

In [83]: z_flat=z.flatten()

In [84]: z_flat.shape
Out[84]: (50,)

但是为了完全调试，您需要修复控制器函数中的上述错误。在

Answer 2

把这个单独贴出来，因为我有你的代码。好吧，运行并产生产出：P

实际上，一个很大的问题是一些我没有注意到的秘密广播，但我改变了很多其他的事情。在

首先，隐形广播是，如果你把一个一维函数与一个参数集成，odeint返回一个列向量，然后当你用这个结果填充行向量时，你就得到了一个2d数组（矩阵）。例如：

In [704]: a
Out[704]: array([0, 1, 2, 3, 4])

In [705]: b
Out[705]: 
array([[0],
       [1],
       [2]])

In [706]: a+b
Out[706]: 
array([[0, 1, 2, 3, 4],
       [1, 2, 3, 4, 5],
       [2, 3, 4, 5, 6]])

我们得到了速度的输出，它是一个列向量，比如b，然后把它加到其他时间函数中，得到一个矩阵。在

---

关于递归，我想我解决了这个问题。这两个导数函数应该取一个标量t，此时它们计算导数。要做到这一点，vderivs需要做积分，它应该一直做到t。所以我重新定义了它们：

^{pr2}$

当然uderivs也可以，但可以写得更简单：

def uderivs(v, t):
    return vr - v

然后，确保velocity和controller传递正确的值（使用v0而不是{}作为起始速度）：

def controller(currentV, time):
    z = integrate.odeint(uderivs, currentV, time)
    return kp*(vr-currentV) + ki*z.squeeze()

def velocity(desired, theta, time):
    return integrate.odeint(vderivs, desired, time)

谁知道物理原理是否正确，但这给出了：

注意，它还没有达到期望的速度，所以我增加了求解的时间

time = np.linspace(0,50,50) #time

{2美元^

以下是我运行的所有代码：

import matplotlib.pylab as plt
import numpy as np
import scipy.integrate as integrate

##Parameters
kp = .5 #proportional gain
ki = .1 #integral gain
vr = 30 #desired velocity in m/s
Tm = 190 #Max Torque in Nm
wm = 420 #engine speed
B = 0.4 #Beta
an = 12 #at gear 4
p = 1.3 #air density
Cd = 0.32 #Drag coefficient
Cr = .01 #Coefficient of rolling friction
A = 2.4 #frontal area

##Variables
m = 18000.0 #weight
v0 = 20. #starting velocity
t = np.linspace(0, 20, 50) #time
dt = .1
theta = np.radians(4) #Theta

def torque(v):    
    return Tm * (1 - B*(an*v/wm - 1)**2)  

def vderivs(v, t):
    ts = np.arange(0, t, dt)
    v1 = an * controller(v, ts) * torque(v)
    v2 = m*Cr*np.sign(v)
    v3 = 0.5*p*Cd*A*v**2
    v4 = m*np.sin(theta)
    vtot = v1+v2+v3+v4*(ts>=10)
    return vtot/m

def uderivs(v, t):
    return vr - v

def controller(currentV, time):
    z = integrate.odeint(uderivs, currentV, time)
    return kp*(vr-currentV) + ki*z.squeeze()

def velocity(desired, theta, time):
    return integrate.odeint(vderivs, desired, time)

plt.plot(t, velocity(v0, theta, t), 'k-', lw=2, label='velocity')
plt.plot(t, controller(v0, t), 'r', lw=2, label='controller')
plt.legend()
plt.show()

Answer 3

这并不是一个真正的答案，而是我注意到的对您代码的一些注释。在

正如@qmorgan指出的，您已经将controller函数中的一个参数命名为与另一个函数相同：velocity请尝试在变量名中保持一致，以避免发生此类情况。例如，所有常量都是简短的缩写。所以，当你命名函数中的参数时，也许用词，这样你就有了一个习惯，知道你在哪里用了什么。在

你也犯了一些类似的错误：你在函数中有一个参数，但是在函数体中却忽略了它，而是使用了一个常量。例如，您定义了velocity来获取(mass, desired, theta, t)，但是您传递了它(m, v, theta, time)，其中v是您的开始速度。小心！您没有注意到这个错误，因为事实上，velocity忽略了desired，而是使用了全局常量vr。相反，这一部分应该有

def velocity(mass, desired, theta, t):
    return integrate.odeint(vderivs, desired, t)

plt.plot(time, velocity(m, vr, theta, time), 'k-')

---

要将列表转换为numpy数组，可以使用np.array(vel, float)而不是np.array([x for x in vel], float)，因为[x for x in vel]与{}本身相同。在

---

在vderivs中，通过t的循环可以完全消除，但至少我认为它应该更像：

^{pr2}$

它不会忽略输入t。在

---

现在，uderivs取当前速度和全局定义的期望速度，并返回它们的差：

def uderivs(v, t):  
    return vr - v   

但是您总是传递它vr（请参见controller的定义），所以每次您对它进行积分时，它都会返回一个长度为t.size和值为vr的平面数组。在

---

您可能希望theta = 4而不是theta = 4

---

numpy中已经存在函数np.sign，不需要实现自己的函数。在