时间戳

首页 > 解决方案 > 在solve_ivp中使用雅可比矩阵

python - 在solve_ivp中使用雅可比矩阵

问题描述

我正在尝试将solve_ivp 内置函数与雅可比矩阵一起使用。我收到此错误消息：

UserWarning：以下参数对所选求解器无效：jac。.format(", ".join(" {}".format(x) for x in extraneous)))

我不确定这是语法错误还是雅可比矩阵中的错误。

我的一些代码

def MyIntFun(t,y):
    return DivWatFlux(y,t,wPar, sPar, RPar,dzIN, dzN)

def jacfunc(t,y):
    jac = Richardsmatrix(y, t, wPar, sPar, RPar, dzIN, dzN)
    return jac

mt.tic()
hwODE = spi.solve_ivp(MyIntFun, [tout[0], tout[-1]], hw0, method='RK45', 
    vectorized=True,rtol=1e-4, jac=jacfunc) 
mt.toc()

这是雅可比行列式的代码：

def Richardsmatrix(hw, t, wPar, sPar, RPar, dzIN, dzN): #the jacobian (a,b,c)

    K=Kfun(hw, wPar, m)
    kRes =RPar.kRobBotR

    a = np.zeros(nN)
    b = np.zeros(nN)
    c = np.zeros(nN)


    ii=np.arange(1, nN-1)

    a[ii] = K[ii, 0]/(dzIN[ii, 0]*dzN[ii-1, 0]) #middle nodes
    a[0] = 0 #a1l is zero and not in matrix 
    a[nN-1] = K[nN-1, 0] / (dzIN[nN-1, 0] * dzN[nN-2, 0]) 

    b[ii] = -K[ii, 0] / (dzIN[ii,0] * dzN[ii-1,0])- K[ii+1,0] / (dzIN[ii,0] * dzN[ii, 0]) #middle nodes
    b[0] = -K[1, 0]/(dzIN[0 , 0]*dzN[0, 0])- kRes / dzIN[0 ,0]
    b[nN-1] = -K[nN-1, 0] / (dzIN[nN-1, 0] * dzN[nN-2, 0])

    c[ii] = K[ii+1,0] / (dzIN[ii,0] * dzN[ii, 0]) #middle nodes
    c[0] = K[1,0] / (dzIN[0, 0] * dzN[0, 0]) 
    c[nN-1] = 0

    B = np.diag(a[1:nN], -1)+ np.diag(b, 0) + np.diag(c[0:nN-1], 1)
    sB = sp.sparse.csc_matrix(B)
    return

我和我的团队在编程方面还很陌生，因此非常感谢任何帮助！

标签： pythonsolver

解决方案

有两个不同的问题。

首先，当method='RK45'传递给时solve_ivp，求解器（在本例中为 Runge-Kutta 4/5）无法使用雅可比行列式。尝试传递solver='Radau',solver='BDF'或solver='LSODA', 因为这些使用雅可比行列式，根据文档（特别是jac关键字参数已记录）。

我看到的另一个问题是Richardsmatrix你在那里的函数有一个不返回任何东西的 return 语句。您可能想尝试一下return sB（假设这sB是您希望函数返回的矩阵）。

推荐阅读