普通微分方程索引误差

Question

我开始学习python（这是我第一次尝试使用Python脚本）。我正在尝试求解，然后随后绘制一系列微分方程，这些方程式如何链接以下方程中的化合物浓度H2O CO2＆> o2 CH2O。

但是，当我运行脚本时：

import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
import math
y=[3.0,1.0,3.0,1.0]
g=15 
i=25 
tspan= np.linspace(0,5,100)
# h2o+co2 <-> o2+ch2o
def chemanalysis(y,t):
    x =i*math.sin(math.pi*t/12)**2
    dydt= -x*(y[0]*y[1]) + g*(y[2]*y[3]) #co2
    dydt= -x*(y[0]*y[1]) + g*(y[2]*y[3]) #h2o
    dydt= x*(y[0]*y[1])  - g*(y[2]*y[3]) #o2
    dydt= x*(y[0]*y[1])  - g*(y[2]*y[3]) #ch2o

    return dydt

y=odeint(chemanalysis,1,tspan)
print(y)

我有以下错误：

Traceback (most recent call last):
File "first.py", line 25, in <module>
        y=odeint(chemanalysis,1,tspan)
        ixpr, mxstep, mxhnil, mxordn, mxords)
      File "first.py", line 16, in chemanalysis
        dydt= -x*(y[0]*y[1]) + g*(y[2]*y[3]) #co2
    IndexError: index 1 is out of bounds for axis 0 with size 1

我知道这可能是我在Python/Programming缺乏经验所带来的非常直接的问题。我非常感谢在我出错的地方或可以改进此代码的地方的任何帮助。我在其他线程上读了其他一些类似的问题/错误消息，但是我无法自己解决该问题。我知道目前没有附上此脚本的图形组件，我认为我会等到首先解决此问题。

问候

Answer 1

您想在这里做两件事：

将y=odeint(chemanalysis, 1, tspan)更改为y2=odeint(chemanalysis, y, tspan)

，要使odeint工作，您需要返回与y相同大小的列表。

def chemanalysis(y,t):
    x =i*math.sin(math.pi*t/12)**2
    dydt = []
    dydt.append(-x*(y[0]*y[1]) + g*(y[2]*y[3])) #co2
    dydt.append(-x*(y[0]*y[1]) + g*(y[2]*y[3])) #h2o
    dydt.append(x*(y[0]*y[1])  - g*(y[2]*y[3])) #o2
    dydt.append(x*(y[0]*y[1])  - g*(y[2]*y[3])) #ch2o

    return dydt

我们可以使它看起来更加柔软，但省略了。

---

更多解释说明为什么需要在chemanalysis中使用append：

dydt= -x*(y[0]*y[1]) + g*(y[2]*y[3]) #co2
dydt= -x*(y[0]*y[1]) + g*(y[2]*y[3]) #h2o
dydt= x*(y[0]*y[1])  - g*(y[2]*y[3]) #o2
dydt= x*(y[0]*y[1])  - g*(y[2]*y[3]) #ch2o
return dydt

以上只需返回dydt作为#ch2o的最后一个分配值，您每次都会在其他数据上丢失并重写dydt。而是将其作为列表，odeint要求您无论如何要做。