[风之影][Lisp大挑战第一季]圆周率

cabinsummer

　　有了电脑，计算圆周率已经不是难事。十几年前，我曾经在VB下开发了一个小小的计算圆周率的程序，当时在奔I的电脑上花了十几分钟计算出一万位，还兴奋了好久。现在电脑硬件已经远远超过了那个时代。往事重提，在明经上发出挑战贴，看各路英雄豪杰一展身手。
　　看哪位大侠的算法最快最好。以计算位数和计算时间评比两个冠军，开源有奖！！！
　　以下是我从网上搜索的20个计算圆周率的公式，大家可以任意选择公式并设计算法，程序贴出时注明公式、计算位数、计算时间、硬件配置。期待大侠风采……
　　规则限制，必须是lisp编程。

highflybir

http://bbs.mjtd.com/thread-67190-1-1.html
我以前发表过这样类似的帖子。呵呵，楼主的这个挑战还是很难。需要数学知识和lisp技巧两者的结合。
要在 lisp中实现这样的算法还真不容易。

highflybir

1. ;;;高精度计算Pi函数
   
2. (defun CalPi (digits n / b c d e f g h r s x)
   
3.   (setq c (/ (1+ n) (/ (log 2) (log 10))))    ;需要迭代的次数
   
4.   (setq c (fix c))                        ;转化为整数
   
5.   (setq e 0 r nil)                                ;存储结果的字符串赋空值
   
6.   (setq h (/ digits 5))                          ;从小数后算起
   
7.   (repeat c                                    
   
8.     (setq f (cons h f))                          ;初始余数为10000 * 2 / 10
   
9.   )
   
10.   (repeat (1+ (/ n 4))                           ;重复1+ 800/4 = 201次
    
11.     (setq d 0)                                   ;每次末位小数为0
    
12.     (setq g (+ c c))                             ;分母。因为每次循环都输出了4位，所以在后面运算时乘以了a，所以这里得 -2
    
13.     (setq b c)      ;分子
    
14.     (setq x nil)
    
15.     (while (> b 0)
    
16.       ;;根据公式，乘以分子
    
17.       (setq d (* d b))     
    
18.       (setq b (1- b))
    
19.       (setq d (+ d (* (car f) digits)))   ;因为每次外循环都输出了4位
    
20.       ;;根据公式，除以分母
    
21.       (setq f (cdr f))
    
22.       (setq g (1- g))
    
23.       (setq x (cons (rem d g) x))   ;带分数的 分子部分
    
24.       (setq d (/ d g))     ;带分数的 整数部分
    
25.       (setq g (1- g))
    
26.     )
    
27.     (setq f (reverse x))
    
28.     (repeat 13      
    
29.       (setq f (cdr f))
    
30.     )
    
31.     (setq s (+ e (/ d digits)))    ;printf("%.4d", e+d/a);
    
32.     (setq r (cons s r))     ;算出的每一项，注意表的每项如果不足4位要加零补全
    
33.     (setq e (rem d digits))    ;e = d % a;
    
34.     (setq c (- c 13))     ;因为精度固定为800位，每输出4位后，相当于精度需求降低了4位，所以每次可以少算13项
    
35.   )
    
36.   (reverse r)      ;把表项反转
    
37. )  
    

;;;高精度计算Pi函数<br /> (defun CalPi (digits n / b c d e f g h r s x)<br />   (setq c (/ (1+ n) (/ (log 2) (log 10))))    ;需要迭代的次数<br />   (setq c (fix c))                        ;转化为整数<br />   (setq e 0 r nil)                                ;存储结果的字符串赋空值<br />   (setq h (/ digits 5))                          ;从小数后算起<br />   (repeat c                                     <br />     (setq f (cons h f))                          ;初始余数为10000 * 2 / 10<br />   )<br />   (repeat (1+ (/ n 4))                           ;重复1+ 800/4 = 201次<br />     (setq d 0)                                   ;每次末位小数为0<br />     (setq g (+ c c))                             ;分母。因为每次循环都输出了4位，所以在后面运算时乘以了a，所以这里得 -2<br />     (setq b c)      ;分子<br />     (setq x nil)<br />     (while (> b 0)<br />       ;;根据公式，乘以分子<br />       (setq d (* d b))     <br />       (setq b (1- b))<br />       (setq d (+ d (* (car f) digits)))   ;因为每次外循环都输出了4位<br />       ;;根据公式，除以分母<br />       (setq f (cdr f))<br />       (setq g (1- g))<br />       (setq x (cons (rem d g) x))   ;带分数的 分子部分<br />       (setq d (/ d g))     ;带分数的 整数部分<br />       (setq g (1- g))<br />     )<br />     (setq f (reverse x))<br />     (repeat 13      <br />       (setq f (cdr f))<br />     )<br />     (setq s (+ e (/ d digits)))    ;printf("%.4d", e+d/a);<br />     (setq r (cons s r))     ;算出的每一项，注意表的每项如果不足4位要加零补全<br />     (setq e (rem d digits))    ;e = d % a;<br />     (setq c (- c 13))     ;因为精度固定为800位，每输出4位后，相当于精度需求降低了4位，所以每次可以少算13项<br />   )<br />   (reverse r)      ;把表项反转<br /> )  <br />

测试部分代码见附件：

cabinsummer

highflybir 发表于 2011-11-1 20:49
http://bbs.mjtd.com/thread-67190-1-1.html
我以前发表过这样类似的帖子。呵呵，楼主的这个挑战还是很难。 ...

看了大师程序，佩服的五体投地！不过没看懂算法，敢问是哪个公式？

highflybir

cabinsummer 发表于 2011-11-1 20:55
看了大师程序，佩服的五体投地！不过没看懂算法，敢问是哪个公式？

http://mathworld.wolfram.com/PiFormulas.html
参见上面的帖子，我的lisp算法出自欧拉改进型算法。--上面网站的23,24,25公式。这个算法比较慢。

另外可以 用Machin公式，可以把速度提高不少。
http://blog.pfan.cn/rickone/24388.html

受LISP语言限制，不太可能实现BBP算法，以及FFT之类的高级算法。

cabinsummer

highflybir 发表于 2011-11-1 21:14
http://mathworld.wolfram.com/PiFormulas.html
参见上面的帖子，我的lisp算法出自欧拉改进型算法。-- ...

我以前计算圆周率用的是麦钦公式，就是我贴出的第四个图片（马庭）。其实即使同样一个公式，算法设计的不同，计算的速度也相差很大。提出这个挑战贴就是希望大家能设计出更快更优化的算法。这个念头起源于x_s_s_1发的悬赏http://bbs.mjtd.com/thread-90030-1-2.html，引起了大家讨论，我想发起lisp算法的挑战来激发大家的思维，圆周率算是简单的了。如果这个话题得到大家的认可，我希望有更多的人参与更多的挑战。

vormittag

这里有用 Common Lisp 实现的算法，好像还用了多线程技术，可惜对 Common Lisp 知之甚少，没有看懂算法。
https://github.com/Bronsa/cl-picalc/commits/master

BTW: Lisp 的水好深啊

VBALISPER

看眼界了,希望多来点这种有挑战性的,可以看到更多的大师浮出水面!

cnks

只有顶的份了