在 python 中使用 scipy.integrate.solve_ivp 时，有没有办法绑定变量的值？

计算科学 颂 scipy 约束
2021-12-01 02:13:19

我想用两个变量 x 和 u 在 python 中解决 IVP，但我需要 u 的值介于 0 和 1 之间。现在它给了我一个对 u 负值的解决方案。这是我的代码。

def model_case_3(t, z, Kmax, k, b, list_Kmax, sigma):
 a=1.5
 x, u= z
 m=1
 r = list_Kmax[0][0]
 dxdt =  (x)*(r*(1-a*u**2)*(1-x/(Kmax*(1-0.999*u**2)))-m/(k+b*u)-0.05)
 dudt = sigma*((-2*a*(b**2)*r*(u**3)+4*a*b*k*r*(u**2)+2*a*(k**2)*r*u-b*m)/((b*u+k)**2))
 return [dxdt, dudt]

 sol = solve_ivp(fun=model_case_3, t_span=[scaled_days[i][j][0], scaled_days[i][j][-1]], y0=[scaled_pop[i][j][0], list_u[0][0][0]], t_eval=scaled_days[i][j],  args=(list_Kmax[0][0], k0, b0, list_b, sigma0))
1个回答

如果精确解确实保持在 [0,1] 内，则求解器可能仍然过于粗略地解析动力学并“跳过”物理边界。解决此问题的一种方法是在solve_ivp调用中使用较低的绝对和相对容差： solve_ivp(..., atol=1e-9, rtol=1e-9, ...例如。

最后，如果解决方案组件超出允许的间隔，您可以通过在每一步之后简单地将解决方案组件投影回边界来解决问题。这是一个务实的解决方案，可能并不总是表现良好......

一件重要的事情是，您的解决方案限制不能“对抗” ODE，即例如，如果 f(u=0) < 0，则不应强制一个变量 u >= 0，因为这只会造成麻烦（振荡，不必要的低时间步长值）在积分。一种解决方案是修改您的 ODE 函数，以便变量的时间导数随着它们越来越接近其边界而缓慢衰减到 0。

您可以阅读本文以了解更多信息。

其它你可能感兴趣的问题
上一篇Heineken Virus C 程序（如果在我的代码中，如何简化和丢弃 else） 下一篇2D-FEM 中参考三角形上 P1 基元的对称性