微分方程- ODEINT

Question

我必须用Python中的ODEINT来求解两个微分方程，它们是：

y''(t) = (l*q)/a * (1/y(p) * [1 - z'(p)*u]
z''(t) = a * (1/y(p) * y'(p)*u 

所以我被告知：

y1=y
y2=y'
z1=z
z2=z'

和

y1' = y2
y2' = y'' = (l*q)/a * (1/y(p) * [1 - z'(p)*u]
z1' = z2
z2' = z''(t) = a * (1/y(p) * y'(p)*u

现在我要解这4个方程。已知l、q、a和u。

我尝试了这样的东西：

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

def rownanie(y, t, l, q, a, u):
    y1, y2, z1, z2 = y
    dydt = [y2, ((l*q)/a)*(1/y1)*(1-z2*u), z2, (a*y2*u)/y1]
    return dydt

l = 1
q = 1
a = 10
u = 0.25

y0 = 0
z0 = 0
t = np.linspace(0, 10, 101)
sol = odeint(rownanie, y0, z0, t, args=(l,q,a,u))
print(sol)

在这方面需要帮助

Answer 1

如果你阅读docs，你会看到odeint

解决了一阶常微分方程或非刚性系统的初值问题：

dy/dt = func(y，t，...)或者func(t，y，...)

其中y可以是一个向量

这种转换是将二阶常微分方程转换为一阶向量常微分方程的标准数学方法。

因此，您创建了一个新的向量变量(为避免混淆，我将其命名为Y)，该变量由向量Y = [y, y_prime, z, z_prime]组成:您的函数实现是正确的。

还要注意，为了进行数值求解，您需要指定所有向量的初始条件，在本例中为y0、z0、y'0和z'0。正如Thomas所指出的，在调用odeint时，您需要将这些值指定为向量的初始值。

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

def rownanie(Y, t, l, q, a, u):
    y1, y2, z1, z2 = Y
    dydt = [y2, ((l*q)/a)*(1/y1)*(1-z2*u), z2, (a*y2*u)/y1]
    return dydt

l = 1
q = 1
a = 10
u = 0.25

y0 = 0
z0 = 0
y0_prime, z0_prime = 0, 0  # you need to specify a value for these too
t = np.linspace(0, 10, 101)
sol = odeint(rownanie, [y0, y0_prime, z0, z0_prime], t, args=(l,q,a,u))
print(sol)