python中具有无限初始条件的ODE答案

【问题标题】：ODEs with infinite initlal condition in pythonpython中具有无限初始条件的ODE
【发布时间】：2013-01-03 18:56:45
【问题描述】：

我有一个二阶微分方程，我想在 python 中求解它。问题是对于其中一个变量，我在0 中没有初始条件，而只有无穷大的值。谁能告诉我应该为 scipy.integrate.odeint 提供哪些参数？能解决吗？

方程：

需要根据时间找到 Theta。它的一阶导数在t=0 处为零。 t=0 不知道 theta，但它会在足够长的时间内变为零。其余的都是已知的。作为近似值I 可以设置为零，从而消除二阶导数，这应该会使问题更容易。

【问题讨论】：

你可以从这个问题stackoverflow.com/questions/1824751/…找到答案
您不能在t 的两个不同值下给出scipy.integrate.odeint 条件。如果不是无穷大的第二个条件，而是在t = t1，您可以将您的解决方案与scipy.integrate.odeint 嵌套在对scipy.optimize.root 的调用中，以在t = 0 处找到tetha 的值，从而为您提供所需的行为t = t1。也许选择一个足够大的t1 可以让你使用这个想法。您可能还想尝试scicomp.stackexchange.com 寻求帮助，找出解决问题的正确策略。
@Jaime，您能否通过使用scipy.optimize.root 并预测t1 值来提供一个粗略的答案？

【解决方案1】：

这远不是一个完整的答案，而是应 OP 的要求发布在这里。

我在评论中描述的方法是所谓的shooting method，它允许将边界值问题转换为初始值问题。为方便起见，我将把你的函数 theta 重命名为 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]

所以现在我们只需要知道a 的值是什么使y = 0 在t = 1000，我们可怜的人的无限，与

a = scipy.optimize.root(y_at_inf, [1])

【讨论】：