Как выйти с этим переполнением в моем коде? - PullRequest
Новая
       80

Как выйти с этим переполнением в моем коде?

0 голосов
/ 06 мая 2020

Итак, я пытаюсь решить EDO, используя адаптивный метод Bulirsch-Stoer, но получаю предупреждение о переполнении. Как это решить? 

    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, но это не сильно изменилось.

1 Ответ

0 голосов
/ 09 мая 2020

Не ответ, просто расширяю свой комментарий с запросом существенных деталей для восстановления ошибки: Использование scipy.integrate.odeint с заданными данными и функцией 

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.

...