如何在 Python 中求解三个二次微分方程？

Question

我刚刚开始使用 Python 进行科学绘图来绘制微分方程的数值解。我知道如何使用模块来求解和绘制单微分方程，但对微分方程系统一无所知。如何绘制以下耦合系统？

我的微分方程组是：

dw/dx=y 和
dy/dx=-a-3*H*y 和
dz/dx=-H*(1+z)

a = 0.1 和 H=sqrt((1+z)**3+w+u**2/(2*a))

我的代码是：

import numpy as N 
from scipy.integrate import odeint
import matplotlib.pyplot as plt

def model(w,y,z,x,H): 
    dwdx=y
    dydx=-a-3*H*y
    dzdx=-H*(1+z)
    a=0.1
    H=sqrt((1+z)**3+w+u**2/(2*a))
    return [w,y,z,H]   

z0=1100      #initial condition
w0=-2.26e-8
y0=-.38e-4
H0=36532.63
b=0
c=10000
x=N.arange(b,c,0.01)

y=odeint(model,y0,x)    #f=Function name that returns derivative values at requested y and t values as dydt = f(y,t) 
w=odeint(model,w0,x)
z=odeint(model,z0,x)

plt.plot(w,x)
plt.plot(y,x)
plt.plot(z,x)
plt.legend(loc='best')
plt.show()

Answer 1

通用 ODE 积分器希望将动态系统简化为抽象的一阶系统。这样的系统具有状态向量空间，微分方程为该空间提供速度向量。这里的状态有 3 个标量分量，它给出了一个 3D 向量作为状态。如果要单独使用组件，ODE函数的第一步是从状态向量中提取这些组件，最后一步是从组件的导数中以正确的顺序组合返回向量。

另外，你需要按照依赖顺序排列计算步骤

def model(u,t): 
    w, y, z = u
    a=0.1
    H=sqrt((1+z)**3+w+u**2/(2*a))
    dwdx=y
    dydx=-a-3*H*y
    dzdx=-H*(1+z)
    return [dwdx, dydx, dzdx]   

然后使用组合的初始状态调用一次积分器

u0 = [ w0, y0, z0]
u = odeint(model, u0, x)

w,y,z = u.T

还请检查绘图函数的参数，一般方案是plot(x,y)。