RuntimeWarning с уравнениями Лоренца

Question

Я пытаюсь выполнить уравнения Лоренца в питоне (я следую упражнению 8.3 от Марка Ньюмана - Вычислительная физика (2012, независимая издательская платформа CreateSpace)) Я уже получил графику, и все выглядит "правильно". Это, вероятно, математическая проблема, а не проблема программирования, но я пишу здесь, чтобы убедиться. Во-первых, это мой код:

      from numpy import arange,array
import pylab

def f(v,t):
    s=10 
    r=28
    b=8/3
    x= v[0] 
    y= v[1]
    z= v[2]
    fx= s*(y - x)
    fy= r*x - y - x*z
    fz= x*y - b*z
    return array([fx,fy,fz],float)



def d(N):
    a=0.0
    b=50.0
    h=(b-a)/N
    r=array([0.0,1.0,0.0],float)
    tpoints=arange(a,b,h)
    xpoints= []
    ypoints= []
    zpoints= []
    for t in tpoints:
        xpoints.append(r[0])
        ypoints.append(r[1])
        zpoints.append(r[2])
        k1 = h*f(r,t)
        k2 = h*f(r+0.5*k1,t+0.5*h)    
        k3 = h*f(r+0.5*k2,t+0.5*h)
        k4 = h*f(r+k3,t+h)
        r += (k1+2*k2+2*k3+k4)*(1/6)
    return tpoints,xpoints,ypoints,zpoints



for i in range (1,6):
    N=10**i
    pylab.plot(d(N)[0],d(N)[1], label=N)
    pylab.xlabel("t")
    pylab.ylabel("x(t)")
    pylab.title("Gráficos x em função de t")
    pylab.legend()
    pylab.show()
    pylab.plot(d(N)[0],d(N)[2], label=N)
    pylab.xlabel("t")
    pylab.ylabel("y(t)")
    pylab.title("Gráficos y em função de t")
    pylab.legend()
    pylab.show()
    pylab.plot(d(N)[0],d(N)[3], label=N)
    pylab.xlabel("t")
    pylab.ylabel("z(t)")
    pylab.title("Gráficos z em função de t")
    pylab.legend()
    pylab.show()
    pylab.plot(d(N)[1],d(N)[3], label=N)
    pylab.xlabel("x")
    pylab.ylabel("z(x)")
    pylab.title("Gráficos z em função de x")
    pylab.legend()
    pylab.show()

Это дает мне графики для решения проблемы, и я думаю, что это правильно.

Когда я иду из i=1к i=3это дает мне эту ошибку:

Я думаю, что это связано с математической проблемой, но когда я ищу ошибку, я получаю что-то с массивами. Так что я проверяю это.

За iравно или больше 4, это не доставило мне проблем.

Answer 1

Для системы Лоренца с RK4 требуется размер шага около 0.05или меньше, т. N=10**4или больше для вашей конструкции. См. почти дубликат аттрактора Лоренца с питоном Рунге-Кутта .

Для больших размеров шага, то есть случаев, в которых вы получаете ошибку, метод будет возвращать хаотические результаты, которые в основном не связаны с точным решением системы и любыми связанными с ним границами. Таким образом, расхождение с переполнением с плавающей запятой возможно из-за квадратичных суперлинейных членов.