我有一个二阶微分方程,我想在python中解决它.问题是,对于其中一个变量,我没有0的初始条件,只有无穷大的值.可以告诉我应该为scipy.integrate.odeint提供哪些参数?可以解决吗?
方程:
需要在时间方面找到Theta.它的一阶导数在t = 0时等于零. θ在t = 0时是未知的,但在足够大的时间内它变为零.所有其他人都知道.作为近似值,我可以设置为零,从而删除应该使问题更容易的二阶导数.
解决方法:
这远不是一个完整的答案,而是根据OP的要求发布在这里.
我在评论中描述的方法是所谓的shooting method,它允许将边界值问题转换为初始值问题.为方便起见,我将把你的函数θ重命名为y.要以数字方式求解方程,首先要使用两个辅助函数z1 = y和z2 = y’将其转换为一阶系统,因此您的当前方程式
I y'' + g y' + k y = f(y, t)
将被重新系统化
z1' = z2
z2' = f(z1, t) - g z2 - k z1
你的边界条件是
z1(inf) = 0
z2(0) = 0
首先,我们设置函数来计算新矢量函数的导数:
def deriv(z, t) :
return np.array([z[1],
f(z[0], t) - g * z[1] - k * z[0]])
如果我们有一个条件z1 [0] = a我们可以在t = 0和t = 1000之间以数字方式解决这个问题,并在最后一次得到y的值
def y_at_inf(a) :
return scipy.integrate.odeint(deriv, np.array([a, 0]),
np.linspace(0, 1000, 10000))[0][-1, 0]
所以现在我们需要知道的是,在t = 1000时,y = 0,我们的穷人的无限,是什么价值
a = scipy.optimize.root(y_at_inf, [1])