Matlab的“OutputFcn”在Python中的实现

Question

我正在尝试使用Python的solve\u ivp来解决一个ODE。但是，我想根据当前解决方案和以前解决方案之间的比较动态地更改ODE的右侧。这背后的想法是，我的右手边是一个向量场，我想通过基于上一个解的方向反转右手边来确保向量场的方向性。你知道吗

实现方法如下：我想检查右手边函数定义中上一个解和向量场之间的点积。如果点积为负，则右侧乘以-1。你知道吗

因此，我需要访问ODE解算器的先前状态，并将其与当前迭代进行比较。在MATLAB中，在求解ODE时可以使用“OutputFcn”。该函数在每次积分器迭代后调用。因此，在函数中，可以简单地将状态提取为变量，并在下一次迭代中使用它。我还没有找到类似Python的东西。你知道吗

def RHS(timesnotused,x):
    out = solve_ivp(doubleGyreVar, [0,T/2, T], [x[0], x[1], 1, 0, 0, 1], rtol = 1e-10, atol=1e-10)
    output = out.y
    J = output[2:,-1].reshape(2,2)
    CG = np.matmul(J.T , J)
    lambdas, xis = np.linalg.eig(CG)
    xi_1 = xis[np.argmin(lambdas)] 
    xi_2 = xis[np.argmax(lambdas)]
    lambda_1 = np.min(lambdas)
    lambda_2 = np.max(lambdas)

    alpha = ((lambda_2-lambda_1) / (lambda_2+lambda_1))**2

    sign = 1
    if np.dot(xlast,xi_1) < 0:
        sign = -1

    return(sign*alpha*xi_1)

可以看出，我希望“xScript”是前一个解决方案，并用当前迭代的XIY1检查。不知何故，每次迭代都需要更新xlast。你知道吗

Answer 1

如果你要解的微分方程是

dx/dt = sign(g(x)) * F(x)

对于一些足够光滑的函数g, F，则有一个不连续的右侧，所有高级数值解算器在接近这个跳跃奇点时都会产生无意义的结果。你知道吗

求解这样一个多相系统的最清晰的方法是只向数值解算器呈现连续的右侧，并通过同样存在于scipy.integrate.solve_ivp中的事件机制来处理相变。你知道吗

在一个类似的问题中，我用scipy的工具探索了一种机制，在odeint returns wrong results for an ODE including descrete function中用符号函数产生不连续性