首页 > 解决方案 > 使用 Lorenz 方程的运行时警告

问题描述

我正在尝试在 python 中做 Lorenz 方程(我正在关注 Mark Newman - Computational Physics(2012,CreateSpace Independent Publishing Platform)的练习 8.3)我已经得到了图形,一切看起来“正确”。这可能是一个数学问题,而不是一个真正的编程问题,但我在这里发布以确保。首先,这是我的代码:

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=1i=3我这个错误:

在此处输入图像描述

我认为这与数学问题有关,但是当我搜索错误时,它给了我一些带有数组的东西。所以我正在检查它。

对于i等于或大于 4,它没有给我带来任何问题。

标签: pythonpython-3.xruntime-errorodelorenz-system

解决方案


带有 RK4 的 Lorenz 系统需要大约或更少的步长0.05,也就是说,N=10**4对于您的构造,或更大。使用 Runge-Kutta python查看近似重复的Lorenz 吸引子

对于较大的步长,即出现错误的情况,该方法将返回混沌结果,这些结果大多与系统的精确解以及与之相关的任何界限无关。因此,由于二次超线性项,浮点溢出的发散是可能的。


推荐阅读