haskell - 逼近haskell中的导数时的令人费解的行为

Question

我已经定义了一个类型类 Differentiable 可以由任何可以在无穷小上操作的类型来实现。这是一个例子:

class Fractional a => Differentiable a where
    dif :: (a -> a) -> (a -> a)
    difs :: (a -> a) -> [a -> a]
    difs = iterate dif

instance Differentiable Double where
    dif f x = (f (x + dx) - f(x)) / dx
        where dx = 0.000001

func :: Double -> Double
func = exp

我还定义了一个简单的 Double -> Double 函数来区分。

但是当我在 ghc 中测试时，会发生这种情况:

... $ ghci
GHCi, version 8.8.4: https://www.haskell.org/ghc/  :? for help
Prelude> :l testing
[1 of 1] Compiling Main             ( testing.hs, interpreted )
Ok, one module loaded.
*Main> :t func
func :: Double -> Double
*Main> derivatives = difs func
*Main> :t derivatives
derivatives :: [Double -> Double]
*Main> terms = map (\f -> f 0) derivatives
*Main> :t terms
terms :: [Double]
*Main> take 5 terms
[1.0,1.0000004999621837,1.000088900582341,-222.0446049250313,4.440892098500626e8]
*Main>

e^x|x=0 的 n 次导数的近似值为:

[1.0,1.0000004999621837,1.000088900582341,-222.0446049250313,4.440892098500626e8]

给定设置，一阶和二阶导数是完全合理的近似值，但突然之间，func 在 0 处的三阶导数是... -222.0446049250313! 怎么样!!？

Answer 1

您在这里使用的方法是 finite difference method一阶精度。

Layman 的翻译:它有效，但从数字上来说相当垃圾。具体来说，因为它只有一阶精度，所以即使使用精确实数运算，您也需要那些非常小的步骤才能获得良好的精度。您确实选择了小步长，这很好，但是小步长会带来另一个问题:舍入误差。您需要将 f (x+δx) - f x 与 δx 的差取值，这意味着差异很小，而单个值可能很大。这总是会导致 float 不准确——例如

Prelude> (1 + pi*1e-13) - 1
3.141931159689193e-13

这实际上可能并没有那么大的伤害，但是由于您需要除以 δx，因此会增加错误。

当您使用更高阶的导数时，这个问题会变得更糟/复杂，因为现在每个 f' x 和 f' (x+δx) 已经一个(不相同的!)boosted 错误，因此获取差异并再次提升显然是灾难的根源。

解决问题的最简单方法是切换到二阶精确方法，明显是中心差异。然后您可以将步长放大很多，从而在很大程度上避免舍入问题:

Prelude> let dif f x = (f (x + δx) - f(x - δx)) / (2*δx) where δx = 1e-3
Prelude> take 8 $ ($0) <$> iterate dif exp
[1.0,1.0000001666666813,1.0000003333454632,1.0000004990740052,0.9999917560676863,0.9957312752106873,8.673617379884035,7806.255641895632]

你看到前几个导数现在很好，但最终它也变得不稳定——当你迭代它时，任何 FD方法都会发生这种情况。但这无论如何都不是一个好方法:请注意，对第 n 个导数的每次评估都需要对第 n-1 个导数进行 2 次评估。所以，复杂度在导数上是指数。

逼近不透明函数的 n-th 导数的更好方法是将 n-th 阶多项式拟合到它并以符号/自动方式对其进行微分。或者，如果函数不是不透明的，则象征性地/自动地区分自己。