curve_fit與未知係數的ODE

Question

我試圖求解方程Iy''+ b | y'| y'+ ky = 0並且將係數擬合到數據。

這是代碼我迄今（準備運行）：

import numpy as np 
import matplotlib.pyplot as plt 
from scipy.integrate import odeint 
import pandas as pd 
from scipy.optimize import curve_fit 

# Define derivatives of function 
def f(y, t, I, b, k): 
    theta, omega = y 
    derivs = [omega, -b/I * omega * abs(omega) - k/I * theta] 
    return derivs 

# integrate the function 
def yint(t, I, b, k, y0, y1): 
    y = odeint(f, [y0, y1], t, args=(I, b, k)) 
    return y.ravel() 

# define t and y to fit 
y0 = [.5, 0] 
t = np.arange(0, 3, .01) 
y = np.cos(t)*np.e**(-.01*t) 

# fit 
vals, cov = curve_fit(yint, t, y, p0=[0.002245, 1e-5, 0.2492, y0[0], y[1]]) 

然而，當我運行的功能，我得到的錯誤：

Traceback (most recent call last): 
    File "---", line 24, in <module> 
    vals, cov = curve_fit(yint, t, y, p0=[0.002245, 1e-5, 0.2492, y0[0], y[1]]) 
    File "---.py", line 578, in curve_fit 
    res = leastsq(func, p0, args=args, full_output=1, **kw) 
    File "---.py", line 371, in leastsq 
    shape, dtype = _check_func('leastsq', 'func', func, x0, args, n) 
    File "---.py", line 20, in _check_func 
    res = atleast_1d(thefunc(*((x0[:numinputs],) + args))) 
    File "---.py", line 447, in _general_function 
    return function(xdata, *params) - ydata 
ValueError: operands could not be broadcast together with shapes (600,) (300,) 

任何思考如何解決這個問題？

Answer 1

問題是，函數yint返回shape（300，）參數的形狀數組（600）。再考慮一下yint：它通過將它表示爲兩個一階方程組來解決二階微分方程。所以y = odeint(...)的結果有兩列，一列是解決方案本身，另一列是其衍生物。它的形狀是（300,2）。將解決方案及其衍生物與ravel混合在一起沒有任何意義。相反，你只應該採取實際的解決方案，那就是你適合的東西。

所以，更換

return y.ravel() 

與

return y[:, 0]