Mathematics & Programme

高精度求对数

这是一个有点难度的问题, 因为现成的公式:

    ln(1+x)=x-x²/3+x³/3-...(-1)^k-1*x^k/k+...      ( |x|<1 )   ①

不满足计算其整个定义域的要求.

1.

幸好, 有一个双曲线余弦函数与自然对数的关系式:

   ln(x)=2arctanh((x-1)/(x+1))                         ②

而 arctanh(y)= y + y^3/3 + y^5/5 + ...       (y≤1)    ③

我们先稍微研究一下②式的性质.

设 y=(x-1)/(x+1), 则②式化为ln(x)=2arctanh(y).

对于函数 y=(x-1)/(x+1) , 当 x 增大时 y 永远也大不过 1 ,完全可以用③式计算, 但是 x 较大时 y 非常接近 1 ,使得级数③,的高精度计算非常慢. 只有当 y 接近 0 (x 接近 1) 时, 速度才快.

想到对数的性质:

   ln(x)=y+ln(x/e^y)      (y 是任意实数)           ④

这样就可以通过选择适当的 y 值使 x/e^y尽量接近1,从而方便使用①或②③式来计算.并且 x 以相同程度接近 1 时(x-1)/(x+1)比(x-1)更有利. 且③比①收敛得更快.因此建议使用②③式来计算.

2.

剩下的问题是如何找这个 y .

y 应该是 ln(x) 近似值.因为这样才能使 ln(x)-y 尽可能地接近 ln(1) .

如果 x 在double或extended类型的范围内那就好办了: 把 x 转换为double (或extended)再 ln(x)就得到y的近似值.

如果没那么好采, 上述办法行不通, 怎办? 联想到常用对数的一个性质:

   lg(x)=x的位数-1+lg(x的小数点移至其第一位后的数)

事实上它是④式的特例. 当 x 以非十进制表示时(例如二进制)上式稍加修改就可以用:

   log(x)=x的位数-1+log(x的小数点移至其第一位后的数)

接着用换底公式 ln(x)=log(x)/log(e) 这个问题就解决得差不多了, 通常 log 的底取的与内部的进位制一样, 比如 2, 10, 65536, 1000000000...

(3.1)

推导方法是 令x=e^t, ln(x)=t=2ln(x^0.5) -> ln(x-1+1)=2ln(x^0.5-1+1), 便得上式.

也可以用

(同②) 和 (3.2)

推导是用tanh(x)的半角公式. 参见 双曲线三角函数 底部.

把x化小再用幂级数计算.

(3.1)式给出了两个计算式, 有前者适用于x非常大的时候(以至于可以忽略 1), 后者适于 x 较小时, 因为此时 √x+1≈1 , 用前者计算时会产生较大误差(如果用幂级数法计算开方则前式似乎更好, 但比起用牛顿迭代法开方结合后式还是慢些).

Last update 10/9/2004