如何在我的代码中出现这个溢出？答案

【问题标题】：How to come out with this Overflow in my code?如何在我的代码中出现这个溢出？
【发布时间】：2020-08-21 22:53:04
【问题描述】：

所以，我尝试使用 Bulirsch-Stoer 自适应方法解决 EDO，但我收到了溢出警告。如何解决？

    from math import sin,cos,pi
from numpy import arange, array, empty, shape,float64
import matplotlib.pyplot as plt
import time

tempo_inicial = time.time()

#Definindo as constantes do exercicio
a = 1
b = 3
ti = 0
tf = 20
prec = 10e-10
x0 = y0 = 0
H = 20
nmax = 8

#Definindo as funções

def f(r,t):
    x = r[0]
    y = r[1]
    f0 = 1-(b+1)*x+a*x**2*y
    f1 = b*x-a*x**2*y
    return(array([f0,f1],float64))

我的函数可以出现在 f(r,t) 中。我也尝试使用 numpy float64，但变化不大。

【问题讨论】：

当前形式的问题是题外话，因为您询问 ODE 系统的理论属性。没有提供可能是错误的代码。方程的多项式特征很可能会导致正反馈，从而导致解在有限时间内爆炸，这将在任何求解器中产生溢出错误。
由于使用scipy.integrate.odeint 时不会发生动态爆炸，因此错误必须出现在此处未显示的部分代码中。
Bulirsch-Stoer 的实现与 stackoverflow.com/questions/38945920/… 的实现有何不同？您是否知道这种方法不合标准，就像 Euler 不如 RK4 一样？

标签： python numerical-methods numerical-integration

【解决方案1】：

不是答案，只是扩展我的评论，要求提供必要的细节来重建错误：使用scipy.integrate.odeint with the given data and function

from numpy import linspace
from scipy.integrate import odeint

time = linspace(ti, tf, 21)
rsol = odeint(f, [x0,y0], time, atol=prec, rtol=0.01*prec)
for r,t in zip(rsol,time):
    x,y = r;
    print(f't={t:.2f}, x={x:.12f}, y={y:.12f}')

我在没有任何警告的情况下得到结果

t= 0.00, x=0.000000000000, y=0.000000000000
t= 1.00, x=0.251223552118, y=0.558443392482
t= 2.00, x=0.269384788526, y=1.278992302288
t= 3.00, x=0.285993102849, y=1.984958673784
t= 4.00, x=0.306864172872, y=2.668107280985
t= 5.00, x=0.334793283816, y=3.320108076380
t= 6.00, x=0.375712107368, y=3.925439075135
t= 7.00, x=0.446840813738, y=4.447011243525
t= 8.00, x=0.636678633858, y=4.735903265638
t= 9.00, x=3.584850495692, y=1.405051562362
t=10.00, x=1.764258113375, y=1.397239762022
t=11.00, x=0.648060470900, y=2.400877519264
t=12.00, x=0.384539021730, y=3.190476648398
t=13.00, x=0.379538742491, y=3.821208432875
t=14.00, x=0.436146645718, y=4.361001236766
t=15.00, x=0.588205110586, y=4.711540741504
t=16.00, x=2.339363929033, y=2.958090116701
t=17.00, x=2.049441930465, y=1.253998459089
t=18.00, x=0.751410490789, y=2.250138349553
t=19.00, x=0.397736992172, y=3.081666874339
t=20.00, x=0.375103612719, y=3.727915883209

这看起来很有界，您报告的错误的来源无法识别。

请注意，10e-10 是 1e-9。

【讨论】：