python - scipy odeint : sum of conservative ode equations does not remain zero -- is this normal?

Question

假设我们有以下等式:

dy1/dt = f(y1, t)      [1]
dy2/dt = g(y2, t)      [2]

方程式是“保守的”，即应满足以下条件:

dy1/dt + dy2/dt = 0    [3]

使用 scipy.odeint，我发现我可以对像这样的保守方程进行积分，这对于简单的 ODE 系统来说非常好。

但是，对于较大的，我会遇到以下问题。

假设这是我的导数函数:

def deriv_function(y0s, t):
    ...body defines equations 1, and 2...

    print np.sum(ode)
    return ode

注意打印语句。

我在 deriv_function 上使用 scipy.odeint 如下:

odeint(deriv_fun, y0s, [0, 0.5])

由于打印语句，打印出以下内容:

-1.38555833473e-13 <--- note, close to zero
-0.00679107743937
-0.0067907211796
-0.0135814423985
-0.0135810861584
-0.416522145214
-0.416523165887
-0.818209018574
-0.818211056221
-1.21864678558
-1.21864881584
-2.86735888212
-2.8673729885
-2.46855840934
-2.46856658088
-3.70632102566
-3.70631206163
-4.93200749506
-4.93200691488
-6.14577326158
-6.14577268283
-8.53799987128
-8.53799713959
-10.8839304356
-10.8839320212
-13.1845005689
-13.1845021725
-15.4406122011
-15.4406123927
-17.6531469917
-17.653147185
-24.6238415795
-24.6238498033
-31.1628867985
-31.1628947266
-37.2974784594
-37.2974547092
-35.463527103
-35.4635192949
-39.5777426955
-39.5777472677
-43.5137135424
-43.5137108017
-47.2791485087
-47.2791483993
-50.881424906
-50.8814244751
-54.3275507164
-54.3275502654    <--- note, not close to zero

对于较小的方程组(不相同)，打印出以下内容:

-1.13686837722e-13
0.0
0.0
0.0
-1.13686837722e-13
1.13686837722e-13
0.0
0.0
0.0
0.0
0.0
0.0
0.0
-1.13686837722e-13
0.0
0.0
0.0
1.13686837722e-13
0.0
0.0
0.0
0.0
-1.13686837722e-13
0.0
0.0
-1.13686837722e-13
0.0
0.0
0.0
0.0
0.0
0.0
-1.13686837722e-13
0.0
0.0
5.68434188608e-14
-5.68434188608e-14
-5.68434188608e-14
0.0

我花了很多时间试图弄清楚这个问题是否与我定义方程式的方式有关，现在我相当确定情况并非如此。为了证实这一点，我想问:非保守方程组是否有可能首先打印出接近于零的 ode 和，然后打印出不接近于零的值？

思考这个问题的另一种方式:在更大的方程组中，ode sum initial 在增加之前打印出大约为零。求解器的内部可能会发生什么导致此问题发生？

---

致后代:have a look at the paralell discussion here at scicomp.stackexchange too .

Answer 1

您看到的是数字错误。我不知道你所说的“更大的方程”是什么意思，但如果你的意思是大的长表达式，你可能会在函数的评估中进行四舍五入(你可以很容易地检查这个)。如果是这种情况，请进行数学运算以获得更好更准确的表达式。

另一件事是您同时对两个方程进行积分。如果你能解耦它们，就独立解决它们。该算法具有自适应步长选择，该值可能不同于使用一个或多个方程式:当您有多个方程式时，步长必须足够好所有 em> 个。

关于您关于肠道的问题，请查看积分曲线。这是 an example来自 MathSE

如果您靠近 fork 点，初始条件中的一个小错误可能会给您带来截然不同的结果，并且它们可能会收敛，也可能不会收敛。

没有通用的解决方法。我认为您的第一个赌注是分析:使用方程式中存在的结构来简化它们。或者，如果您的轨迹附近有极点，请应用转换将它们推开。另一方面，这种差异可以让您估计积分误差。