如何在Haskell中优化数值库的速度

Posted 2021-03-31

我已经发布了一个用于求解延迟微分方程的小数值库：http://github.com/masterdezign/dde

主要技术限制：

1. 使动态变量（x（t），y（t），...）的数量灵活，即给定时刻的动力系统的State。
2. 轻松与利用Data.Vector.Storable的库（例如hmatrix）集成。因此，作为输入/输出，我广泛使用Data.Vector.Storable。

因此，与此解决方案不同：使用How do I optimize numerical integration performance in Haskell (with example)，newtype State = State { _state :: V.Vector Double }而不是data State = State {-# UNPACK #-} !Double {-# UNPACK #-} !Double。但是，库现在运行速度慢了两倍。

问题是：有没有办法将data State = State {-# UNPACK #-}...的速度和newtype State = State { _state :: V.Vector Double }的灵活性同时用于未指定数量的变量？我应该考虑模板Haskell在编译时创建类似data UNPACK的结构吗？

答案

我不会使用任何特定的矢量实现。像Data.Vector这样的可变长度类型是一个糟糕的选择，不仅因为当空间的维度较低时额外的长度信息是一个相当大的开销，也因为你失去了尺寸匹配的任何类型系统保证。

相反，你应该在矢量空间的选择上使一切参数化。即，您有效地使维度成为编译时变量，并且允许向量类型具有一些有意义的子变量。

import Data.VectorSpace
import Data.AdditiveGroup

newtype Stepper1 state = Stepper1 {
  _stepper
    ::  Double
    -> RHS state   -- parameterised in a similar way
    -> state
    -> (Double, Double)
    -> (Double, Double)
    -> state
  }

rk₄ :: VectorSpace v => Stepper1 v
rk₄ = Stepper1 _rk4
 where _rk4 hStep rhs' y₀ ... = y₀ ^+^ (h/6)*^(k₁ ^+^ 2*^k₂ ^+^ 2*^k₃ ^+^ k₄)
         where k₁ = rhs' (y₀, ...)
               k₂ = rhs' (y₀ ^+^ (h/2)*^k₁, ...)
               k₃ = rhs' (y₀ ^+^ (h/2)*^k₂, ...)
               k₄ = rhs' (y₀ ^+^ h*^k₃, ...)

然后，用户可以选择它将是什么特定实现。对于二维向量，标准选择是来自V2包的linear;它与VectorSpace的free-vector-spaces实例一起重新导出。但是对于测试，你也可以使用普通的旧元组，它们在VectorSpace中有一个vector-space实例。当然，从hmatrix环绕类型也是可能的，但这对性能来说真的不好 - 如果有必要，最好转换最终结果。

为了获得最佳性能，您可能需要使用一些{-# INLINE #-}编译指示。 Bang模式OTOH通常不会带来太多的性能优势 - 最重要的是类型是严格的，并且没有装箱。毫无疑问，在每个变量定义之前预先设置爆炸是没有用的 - 无论如何都不会产生任何影响，因为CAF仅在使用时进行评估。

我很高兴听到你最终获得的表现！如果它明显比你原来的State {-# UNPACK #-} !Double {-# UNPACK #-} !Double更糟糕，那就是我们应该调查的东西。