超越函数的高精度计算

winxos

现在计算pi的公式有很多不就是用的 反正切 计算的么？就是用高精度的反正切运算得出的啊？
所以我觉得 超越函数应该是可以找到类似的快速逼近公式的？

liangbch

反正切函数计算pi 的复杂为O(n^2), 仍属比较慢的计算方法。
　　对于$log((p+1)/p)$，p为单精度整数，可以使用反正切公式算。当p为$2^31$以内的自然数，可使用简单的多重精度除以（乘以）单精度整数算法。但是对于普通的多重精度浮点数。这个公式就无能为了了。
　　当n很大时，求x的初等函数，且精确到n位有效数字，4#的给出的结论是最好的。即若两个n位整数相乘的复杂度是M(n), 则计算x的初等函数并精确到n位，其复杂度为 k*Log(n)*M(n), k为常数。AGM 算法计算自然对数的复杂度服从这个规则。

　　当n相对较小时，有比这个复杂度更好的方法。我已经探索出一种比较复杂的算法，在n较小时（如小于10000位有效数字），使用基于小质数的对数表的查表算法，求x的自然对数至n位有效数字，速度可大大超过AGM算法。
　　这个算法已经成型。但算法分析和实现目前仍未完成。当n在一定范围内，相信这个算法可轻易打败apfloat,但是是否可超越mathematica,仍是一个未知数。

无心人

和MPFR比速度如何？

liangbch

不知道，我对MPFR的算法和性能不太了解。

仙剑魔

目前找到了2种AGM计算LOG的算法
我没算过具体值，只是测试得到的一些估计值
若p为精度

The logarithm constant: log 2
如同liangbch在#18所言(他的计算应该比较精确)
需要约2log(p)次迭代
麻烦的是要计算累减值

apfloat
这个的AGM和上面的很相似,迭代的次数似乎也一样
不过简单化了
从代码上看就是以下的式子
log(x)=pi/2*(1/agm(1,10^-n)-1/agm(1,x*10^-n))
也就是说前提是要得到所需精度的pi值
然后就是昨天从#13的论文里瞄到的一行字
从sqrt(an),sqrt(bn)计算sqrt(an+2),sqrt(bn+2)
这个实在是个好想法
这样一来,迭代次数就可以减半了

另有传说中的复数AGM算法
虽然我有点想法,不过还没实际测试过

liangbch

原帖由 仙剑魔 于 2009-1-7 11:49 发表 [url=http://bbs.emath.ac.cn/redirect.php?goto=findpost&pid=14382&ptid=351][img]
然后就是昨天从#13的论文里瞄到的一行字
从sqrt(an),sqrt(bn)计算sqrt(an+2),sqrt(bn+2)
这个实在是个好想法
这样一来,迭代次数就可以减半了

迭代次数减半，应该不会有这等好事！！！
不过AGM算法中，当迭代次数超过一半后$a_n$ 和 $b_n$ 很是接近，如此，
1）计算$a_{n+1}$可在$a_n$基础上进行.
   若$a_n>b_n$, 令 $a_n=b_n+d$,则
$a_{n+1}= b_n + d/2$

2） 当迭代次数超过一半后，$b_n$ 和 $b_{n+1}$ 很是接近,  计算$b_{n+1}$可在$b_n$基础上进行，可减少运算量
    因为计算$b_n$时，需要使用牛顿迭代法，即先计算 $t1= a_{n-1}*b_{n-1}$ ,然后牛顿迭代法逐步取精，得到$t2= 1/sqrt(t1)$, 最后计算 $b_n= t1*t2$
   当计算$b_{n+1}$时，是可以预先估算出$b_{n+1}$有多少为有效数字和$b_n$相同的，故在计算t2时，牛顿迭代法不用从头开始，t2直接取上一次迭代过程中的t2的前几位即可。

另外，当$a_n$和$b_n$很接近时，在计算$a_n^2 - b_n^2$ 时，也可采用这样的方法,或许可提高速度
令  $a_n=b_n+d$
   则：$a_n^2-b_n^2=  2b_n*d+d^2$

仙剑魔

用apfloat的形式的那个AGM
迭代次数确实可以减版,我试过了
但是迭代次数减半并不意味着运算量也减半呀
只是用了那个公式后可以用1次开4次方代替连续2次开平方,所以速度有一定上升

The logarithm constant: log 2的形式的AGM
由于每次要计算an^2-bn^2
所以我觉得不适合这样跳跃的计算

无心人

都别争了
一人写一个程序
不用特殊库
或者用一个库

比较下就知道了啦

gxqcn

原帖由 仙剑魔 于 2009-1-7 18:51 发表
用apfloat的形式的那个AGM
迭代次数确实可以减版,我试过了
但是迭代次数减半并不意味着运算量也减半呀
只是用了那个公式后可以用1次开4次方代替连续2次开平方,所以速度有一定上升

The logarithm constant: log ...

大数开方运算中，开平方有特殊的算法（apfloat没实现），速度很快。
所以速度是否上升可能取决于算法库的性能。

无心人

GxQ说的是牛顿二次收敛方法？