使用Python的边界条件的二阶颂歌



我尝试使用 odeint在Python中以y'' + ay' + by + c = 0(二阶微分方程(的形式求解方程。据我了解,odeint仅在y(0(= y1,y'(0(= y2的形式的初始条件下工作。而我的条件是边界:y'(0(= 0,y'(pi/4(= 0。是否有任何方法可以在这种情况下使用odeint?这是我拥有的代码:

from scipy.integrate import odeint
import numpy as np
def lagrange(x, teta):
  y = x[0]
  dy = x[1]  
  xdot = [[],[]]
  xdot[0] = dy
  xdot[1] = (y + dy**2*y*(y**2 + dy**2)**(-1.5) + y*(y**2+dy**2)**(-0.5))/((y**2+dy**2)**(-0.5) - dy**2*(y**2+dy**2)**(-1.5))
  return xdot
phi = np.linspace(0,np.pi/4,100)
U = odeint(lagrange, [u1_0, u2_0], phi)

官方文档指出,它解决了一阶ode-s的僵硬或非建立系统的初始值问题。因此,它无法解决边界价值问题。

您很可能需要实现有限差异方法,有限元方法或拍摄方法来解决问题。

我会亲自推荐有限的差异方法,因为它基本上包括添加3个对角矩阵(其中2个是tridiagonal(,设置解决方案向量并使用线性求解器。

如果您需要设置一些帮助,我不久前就编写了一个扩散 - 添加求解器,您可以向其移植并扩展以包含反应术语(应该是b*Identity matrix(。

根据您的价值,您也可能需要实现人工扩散(如果Peclet号码太大(。

最新更新