使用scipy.integrate.odeint解决颂歌系统（改变常数！）？

2 回答

• 8

  是的，这是可能的 . 在 a 是常量的情况下，我猜你调用了 scipy.integrate.odeint(fun, u0, t, args) ，其中 fun 被定义为你的问题， u0 = [x0, y0, z0] 是初始条件， t 是要求ODE求解的时间点序列， args = (a, b, c) 是要传递的额外参数到 fun .

  在 a 取决于时间的情况下，您只需重新考虑 a 作为函数，例如（给定常量 a0 ）：

  def a(t):
      return a0 * t
  

  然后你必须修改 fun ，它计算每个时间步的衍生物，以考虑先前的变化：

  def fun(u, t, a, b, c):
      x = u[0]
      y = u[1]
      z = u[2]
      dx_dt = a(t) * x + y * z # A change on this line: a -> a(t)
      dy_dt = b * (y - z)
      dz_dt = - x * y + c * y - z
      return [dx_dt, dy_dt, dz_dt]
  

  最后，请注意 u0 ， t 和 args 保持不变，您可以再次调用 scipy.integrate.odeint(fun, u0, t, args) .

  关于这种方法的正确性的一句话 . 数值积分近似的性能受到影响，我不知道究竟是如何（没有理论上的保证），但这里有一个简单的例子：

  import matplotlib.pyplot as plt
  import numpy as np
  import scipy as sp
  import scipy.integrate
  
  tmax = 10.0
  
  def a(t):
      if t < tmax / 2.0:
          return ((tmax / 2.0) - t) / (tmax / 2.0)
      else:
          return 1.0
  
  def func(x, t, a):
      return - (x - a(t))
  
  x0 = 0.8
  t = np.linspace(0.0, tmax, 1000)
  args = (a,)
  y = sp.integrate.odeint(func, x0, t, args)
  
  fig = plt.figure()
  ax = fig.add_subplot(111)
  h1, = ax.plot(t, y)
  h2, = ax.plot(t, [a(s) for s in t])
  ax.legend([h1, h2], ["y", "a"])
  ax.set_xlabel("t")
  ax.grid()
  plt.show()
  

  

  我希望这能帮到您 .

  回复于 2024-05-05T14:58:20+08:00
• 1

  不，这是不可能的，因为求解器将使用您无法控制的内部时间步，并且每个时间步使用该函数的多个评估 .

  如果参数是分段常数，那么您可以从跳转点到跳转点进行积分 .

  如果是分段线性，则可以使用插值函数 .

  对于其他任何事情，您必须实现自己的逻辑，从时间 t 到这些参数的正确值 .

  注意，数值积分的近似顺序不仅受RK方法的阶数限制，而且受微分方程（部分）的微分阶数的限制 .

  回复于 2024-05-05T14:58:20+08:00