haskell - Haskell 中的 Verlet 集成

Question

我是 Haskell 的新手，作为练习，我一直在尝试实现 Joel Franklin 的《物理计算方法》一书中的一些代码(用 Mathematica 编写)。我编写了以下代码以将 lambda 表达式(加速度)作为第一个参数。一般来说，加速度的形式为 x'' = f(x',x,t)，但并不总是所有三个变量。

-- Implementation 2.1: The Verlet Method
verlet _  _  _  _  0 = [] 
verlet ac x0 v0 dt n = take n $ zip [0,dt..] $ verlet' ac x0 v0 0 dt
  where verlet' ac x0 v0 t0 dt = 
          let
            xt = x0 + dt*v0 + 0.5*(dt^2)*(ac x0 v0 t0)
            vt = v0 + dt*(ac x0 v0 t0)
          in
            xt:(verlet' ac xt vt (t0+dt) dt)

在 ghci 中，我将使用以下命令运行此代码(加速函数 a = -(2pi)2x 来自书中的一个练习):

verlet (\x v t -> -(2*pi)^2*x) 1 0 0.01 1000

我的问题是，这并不是真正的 Verlet 方法——这里 xn+1 = xn + dt*vn+1/2*a(xn,vn,n)，而维基百科中描述的 Verlet 方法是 xn+1 = 2*xn - xn-1+a(xn,vn,n)。我将如何重写此函数以更忠实地表示 Verlet 积分方法？

切线，有没有办法写得更优雅、更简洁？是否有线性代数库可以使这更容易？我很欣赏你的建议。

Answer 1

忠实的 Verlet 序列有 xn 取决于 x 的前两个值——xn-1 和 xn-2。这种序列的一个典型例子是斐波那契数列，它具有这个单行 Haskell 定义:

fibs :: [Int]
fibs = 0 : 1 : zipWith (+) fibs     (tail fibs)
                        -- f_(n-1)  f_n

这将斐波那契数列定义为一个无限(惰性)列表。自我引用 tail fibs为您提供上一个术语，以及对 fibs 的引用给你之前的术语。然后将这些术语与 (+) 组合在一起产生序列中的下一项。

您可以对您的问题采取相同的方法，如下所示:

type State = (Double, Double, Double)  -- (x, v, t) -- whatever you need here

step :: State -> State -> State
step s0 s1 = -- determine the next state based on the previous two states

verlet :: State -> State -> [State]
verlet s0 s1 = ss
  where ss = s0 : s1 : zipWith step ss (tail ss)

数据结构 State保存您需要的任何状态变量 - x, v, t, n, ...
函数 step类似于 (+)在斐波那契的情况下，并计算给定前两个状态的下一个状态。 verlet函数确定给定初始两个状态的整个状态序列。