使用 Haskell 中的系列的自然对数给出不精确的结果

Question

所以我写了两个函数来计算变量 x 的自然对数，在将增量和的上限增加到 33000 之后，函数仍然 return 在 ghci 中测试的结果不精确，与从前奏，这里是代码定义：

lnOfx :: Float -> Float
lnOfx x = netSum f 1 33000
    where f i = 2*(1/oddTerm)*((x-1)/(x+1))**oddTerm
            where oddTerm = 2*i-1
          netSum f minB maxB = sum [f i | i <- [minB .. maxB]]

lnOfx2 :: Float -> Float
lnOfx2 x = netSum f 1 33000
       where f i = (1/i)*((x-1)/x)**i
             netSum f minB maxB = sum [f i | i <- [minB .. maxB]]

测试结果：

log 3
1.0986122886681098
lnOfx 3
1.0986125
lnOfx2 3
1.0986122

log 2
0.6931471805599453
lnOfx 2
0.6931472
lnOfx2 2
0.6931473

那么为什么结果会不同，正确的方法是什么（就像 Prelude 中的 log 函数一样）正确计算自然对数？

Answer 1

浮点数学很棘手。可能导致精度损失的一件事是添加具有非常不同幅度的数字。例如，在您的算法中，从 i=25 左右开始，总和中的项足够小，以至于它们不再产生影响：

-- 25t term:    
let f x i = let oddTerm = 2*i-1 in 2*(1/oddTerm)*((x-1)/(x+1))**oddTerm
let y = f 3.0 25

-- summation up to 24 item
let s = 1.098612288668109

-- this will return True, surprisingly!
s + y == s

你可以做的一件事是将数字以相反的顺序相加，这样小的数字在被添加到大的数字之前先被加在一起。

lnOfx :: Float -> Float
lnOfx x = netSum f 1 33000
    where f i = 2*(1/oddTerm)*((x-1)/(x+1))**oddTerm
            where oddTerm = 2*i-1
          netSum f minB maxB = sum $ reverse [f i | i <- [minB .. maxB]]

在我的测试中，这足以让 print (lnOfx 3.0) 和 print (log 3.0) 显示所有相同的数字。

但总的来说，我会推荐阅读数值分析书籍以了解更多有关此类问题的信息。

使用 Haskell 中的系列的自然对数给出不精确的结果

Natural Logarithm using series in Haskell Gives Imprecise Results

math

haskell

libm