首页 > 解决方案 > 如何积分耦合微分方程?

问题描述

我有一个方程组,我一直在尝试让 Python 来求解和绘图,但绘图并不正确。这是我的代码:

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


#function that returns dx/dt and dy/dt
def func(z,t):
    for r in range(-10,10):
                beta=2
                gamma=0.8

                c = z[0]
                tau = z[1]
                dcdt = r*c+c**2-c**3-beta*c*tau**2
                dtaudt = -gamma*tau+0.5*beta*c*tau
    return [dcdt,dtaudt]

#inital conditions
z0 = [2,0]

#time points
t = np.linspace(0,24,100)

#solve ODE
z = odeint(func,z0,t)

#seperating answers out

c = z[:,0]
tau = z[:,1]

print(z)


#plot results
plt.plot(t,c,'r-')
plt.plot(t,tau,'b--')
plt.legend(['c(t)','tau(t)'])
plt.show()

让我解释。我正在研究双扩散对流。我不想对 r 的值做出任何假设,但 beta 和 gamma 是正的。所以我想给它们赋值,而不是给 r。这是我得到的情节,从理解问题来看,图表不正确。tau 图绝对不应该停留在 0 上,而 c 图应该做得更多。我对 Python 比较陌生并且正在学习课程,但我真的很想了解我做错了什么,因此将不胜感激使用简单语言的帮助。

标签: pythonintegration

解决方案


我在您的功能中发现了 2 个问题,您应该检查一下。

    for r in range(-10,10):

在这里,您正在执行一个 for 循环,只是重新评估 dcdt 和 dtaudt。结果,输出值与仅计算 r=9 (循环中的最后一个值)相同

    dtaudt = -gamma*tau+0.5*beta*c*tau

在这里你有dtaudt = tau*(beta*c/2. -gamma)。您的选择tau[0]=0意味着 tau 将保持为 0。

尝试这个:

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

r = 1
beta=2
gamma=0.8

#function that returns dx/dt and dy/dt
def func(z,t):
    c = z[0]
    tau = z[1]
    dcdt = r*c+c**2-c**3-beta*c*tau**2
    dtaudt = -gamma*tau+0.5*beta*c*tau
    print(dtaudt)
    return [dcdt,dtaudt]

#inital conditions
z0 = [2,0.2] #tau[0] =!0.0

#time points
t = np.linspace(0,24,100)

#solve ODE
z = odeint(func,z0,t)

#seperating answers out

c = z[:,0]
tau = z[:,1]

#plot results
plt.plot(t,c,'r-')
plt.plot(t,tau,'b--')
plt.legend(['c(t)','tau(t)'])
plt.show()

推荐阅读