MATLAB的ode15如何處理包含值矩陣的方程的積分？

我翻譯MATLAB代碼到Python，但擔心翻譯之前，我想了解MATLAB及其具體ODE15s求解器是如何解釋的方程式。MATLAB的ode15如何處理包含值矩陣的方程的積分？

我有一個功能腳本，它在主腳本呼籲，而這個功能腳本包含公式：

function testFun=testFunction(t,f,dmat,releasevec) 

    testFun=(dmat*f)+(releasevec.');

在testFunction，t爲時間，F到我求解值，dmat到我很好奇的常量矩陣，並釋放到一個額外的常量向量。在主腳本

的ODE15s解算器的工作原理及其魔術下面幾行：

for i=1:1461 

    [~,f]=ode15s(@(t, f) testFunction(t, f, ... 
    [dAremoval(i), dFWtoA(i), dSWtoA(i), dStoA(i), dFSedtoA(i), dSSedtoA(i); ... 
    dAtoFW(i), dFWremoval(i), dSWtoFW(i), dStoFW(i), dFSedtoFW(i), dSSedtoFW(i); ... 
    dAtoSW(i), dFWtoSW(i), dSWremoval(i), dStoSW(i), dFSedtoSW(i), dSSedtoSW(i); ... 
    dAtoS(i), dFWtoS(i), dSWtoS(i), dSremoval(i), dFSedtoS(i), dSSedtoS(i); ... 
    dAtoFSed(i), dFWtoFSed(i), dSWtoFSed(i), dStoFSed(i), dFSedremoval(i), dSSedtoFSed(i); ... 
    dAtoSSed(i), dFWtoSSed(i), dSWtoSSed(i), dStoSSed(i), dFSedtoSSed(i), dSSedremoval(i)], ... 
    [Arelease(i), FWrelease(i), SWrelease(i), Srelease(i), FSedrelease(i), SSedrelease(i)]), [i, i+1], fresults(:, i),options); 
    fresults(:, i + 1) = f(end, :).';

fresults最初，裏面的F結果爲零的表。這些選項調用odeset以獲得「非負」值。上面的d值矩陣是一個6x6矩陣。我已經計算了所有的d值和發佈值。我的問題是：ode15s如何與測試函數方程中給出的6x6矩陣進行整合？我試圖通過手工來解決這個問題，但還沒有成功。任何幫助將非常感激！！

＃

def func(y, t, params): 
    f = 5.75e-16 
    f = y 
    dmat, rvec = params 
    derivs = [(dmat*f)+rvec] 
    return derivs 

# Parameters 
dmat = np.array([[-1964977.10876756, 58831.976165, 39221.31744333, 1866923.81515922, 0.0, 0.0], 
      [58831.976165, -1.89800738e+09, 0.0, 1234.12447489, 21088.06180415, 14058.70786944], 
      [39221.31744333, 0.84352331, -7.59182852e+09, 0.0, 0.0, 0.0], 
      [1866923.81515922, 0.0, 0.0, -9.30598884e+08, 0.0, 0.0], 
      [0.0, 21088.10183616, 0.0, 0.0, -1.15427010e+09, 0.0], 
      [0.0, 0.0, 14058.73455744, 0.0, 0.0, -5.98519566e+09]], np.float) 
new_d = np.ndarray.flatten(dmat) 

rvec = np.array([[0.0], [0.0], [0.0], [0.0], [0.0], [0.0]]) 

f = 5.75e-16 

# Initial conditions for ODE 
y0 = f 

# Parameters for ODE 
params = [dmat, rvec] 

# Times 
tStop = 2.0 
tStart = 0.0 
tStep = 1.0 

t = np.arange(tStart, tStop, tStep) 

# Call the ODE Solver 
soln = odeint(func, y0, t, args=(params,)) 
#y = odeint(lambda y, t: func(y,t,params), y0, t)

來源

2015-07-10 deelsburyyy

在具有6x6矩陣和6x1向量的系統中，初始狀態y0也需要是6x1向量。您正在提供一個標量。 – LutzL

它說here是ode15s使用後向差分公式分化。你的微分方程是（據我所知）f' = testFunc(t,f)，它有一些向量矩陣計算在函數內。然後就可以通過一個向後差分公式，是更換分化：

f_next = f_prev + h*testFunc(t,f_next);

其中f_prev是向量的初始值。在這裏有在計算沒有重要的區別，只是因爲testFunc（T，F）函數包括6×6矩陣。每次它通過數字創建雅可比矩陣來解決逆問題以找到f_next。

然而，試圖碼算法，MATLAB也可能比我們想象的，因爲MATLAB更難有一些（優化相關或不相關）特殊處理的問題。你應該小心你得到的每個價值。

來源

2015-07-10 06:31:04 crbah

謝謝@cbahadir！我將把我的等式投入到scipy.integrate.ode（'vode'，method ='bdf'）中，因爲該方法允許像你所建議的那樣向後分化。我有一個關於你的帖子的問題，但是「h * testFunc（t，f_nest）」中的h是什麼？ – deelsburyyy

是步長？ – deelsburyyy

是h是步長。 – crbah

從本質上講，你需要改變很少的東西。使用numpy.ndarray作爲向量和矩陣。時間步進可以使用scipy.integrate.ode完成。您需要通過set_f_parameter重新初始化積分器，以瞭解ODE功能的每個變化或者供應矩陣和參數作爲附加功能參數。

更接近matlab界面但限於lsoda是scipy.integrate.odeint。但是，由於您使用求解器來解決棘手的問題，這可能正是您所需要的。

來源

2015-07-10 18:52:03 LutzL

謝謝@LutzL！我花了很多時間使用sp.integrate.ode求解器無濟於事。 odeint求解器似乎工作得好一點。我將把我提出的代碼添加到原始文章中。我一直遇到以下錯誤：ValueError：object too deep for desired array odepack.error：函數調用的結果不是一個適當的浮點數組。我認爲這與我的數組維度有關，但到目前爲止我無法弄清楚... – deelsburyyy

MATLAB的ode15如何處理包含值矩陣的方程的積分？

回答

相關問題