3.5.2 解常微分方程组

integrate模块还提供了对常微分方程组进行积分的函数odeint( )。下面我们看看如何用它计算洛伦茨吸引子的轨迹。洛伦茨吸引子由下面的三个微分方程定义：

这三个方程定义了三维空间中各个坐标点上的速度矢量。从某个坐标开始沿着速度矢量进行积分，就可以计算出无质量点在此空间中的运动轨迹。其中σ、ρ、β为三个常数，不同的参数可以计算出不同的运动轨迹：x(t)、y(t)、z(t)。当参数为某些值时，轨迹出现混沌现象。即微小的初值差别也会显著地影响运动轨迹。下面是洛伦茨吸引子的轨迹计算和绘制程序。

❶程序中先定义一个函数lorenz( )，它的任务是计算出某个坐标点的各个方向上的微分值，它可以直接根据洛伦茨吸引子的公式得出。

❷❸使用不同的位移初始值两次调用odeint( )，对微分方程求解。odeint( )有许多参数，这里用到的4个参数分别为：

●lorenz：它是计算某个位置上的各个方向的速度的函数。

●(0.0, 1.0, 0.0)：位置初始值，它是计算常微分方程所需的各个变量的初始值。

●t：表示时间的数组，odeint( )对此数组中的每个时间点进行求解，得出所有时间点的位置。

●args：这些参数直接传递给lorenz( )，因此它们在整个积分过程中都是常量。

图3-25显示了odeint( )所得到的轨迹。由结果可知，即使初始值只相差0.01，两条运动轨迹也是完全不同的。

图3-25 洛伦茨吸引子：微小的初值差别也会显著地影响运动轨迹

    
    from scipy.integrate import odeint
    import numpy as np 
    
    def lorenz(w, t, p, r, b): ❶
        # 给出位置矢量w和三个参数p、r、b
        # 计算出dx/dt、dy/dt、dz/dt的值
        x, y, z = w.tolist( )
        # 直接与lorenz的计算公式对应
        return p*
(y-x), x*
(r-z)-y, x*
y-b*
z
    
    t = np.arange(0, 30, 0.02) # 创建时间点
    # 调用ode对lorenz进行求解, 用两个不同的初始值
    track1 = odeint(lorenz, (0.0, 1.00, 0.0), t, args=(10.0, 28.0, 3.0)) ❷
    track2 = odeint(lorenz, (0.0, 1.01, 0.0), t, args=(10.0, 28.0, 3.0)) ❸