【飞鸟集】圆周率小数后800位！（更新到2011.11.1）

highflybir

下面给出一个程序，纯lisp的，能把圆周率精确到小数后800位。稍加修改就可以到10000位，当然远没有C语言快。
不过程序运行命令 test,估计10秒中左右就能把这个数打印出来。

2. 
   
3. ;;;highflybird    2008.4.20 Haikou
   
4. (defun c:test (/ a c p t1 t2 time)
   
5.   (setq a 10000 c 2800)
   
6.   (setq t1 (getvar "TDUSRTIMER"))
   
7.   (setq p  (CalPi a c))
   
8.   (setq t2 (getvar "TDUSRTIMER"))
   
9.   (setq time (* 86400 (- t2 t1)))
   
10.   (setq p (strcat "3." (substr p 2)))
    
11.   (princ p)
    
12.   (princ "\n共耗时<秒>:")
    
13.   (princ time)
    
14.   (princ)
    
15. )
    
16. ;;;计算函数
    
17. (defun CalPi (a c / b d e f g h p s x)
    
18.   (setq b 0 e 0 P "")
    
19.   (setq h (/ a 5))
    
20.   (repeat (1+ c)
    
21.     (setq f (cons h f))
    
22.   )
    
23.   (while (> c 0)
    
24.     (setq d 0)
    
25.     (setq g (+ c c))
    
26.     (setq b c)
    
27.     (setq x nil)
    
28.     (while (> b 0)
    
29.       (setq d (* d b))
    
30.       (setq b (1- b))
    
31.       (setq d (+ d (* (car f) a)))
    
32.       (setq f (cdr f))
    
33.       (setq g (1- g))
    
34.       (setq x (cons (rem d g) x))
    
35.       (setq d (/ d g))
    
36.       (setq g (1- g))
    
37.     )
    
38.     (setq f (reverse x))
    
39.     (repeat 14
    
40.       (setq f (cdr f))
    
41.     )
    
42.     (setq s (+ e (/ d a)))
    
43.     (setq s (itoa s))
    
44.     (while (< (strlen s) 4)
    
45.       (setq s (strcat "0" s))
    
46.     )
    
47.     (setq P (strcat p s))
    
48.     (setq e (rem d a))
    
49.     (setq c (- c 14))
    
50.   )
    
51.   p
    
52. )  
    
53. ;;;以后给出完全注释版的
    

<br /> ;;;highflybird    2008.4.20 Haikou<br /> (defun c:test (/ a c p t1 t2 time)<br />   (setq a 10000 c 2800)<br />   (setq t1 (getvar "TDUSRTIMER"))<br />   (setq p  (CalPi a c))<br />   (setq t2 (getvar "TDUSRTIMER"))<br />   (setq time (* 86400 (- t2 t1)))<br />   (setq p (strcat "3." (substr p 2)))<br />   (princ p)<br />   (princ "\n共耗时<秒>:")<br />   (princ time)<br />   (princ)<br /> )<br /> ;;;计算函数<br /> (defun CalPi (a c / b d e f g h p s x)<br />   (setq b 0 e 0 P "")<br />   (setq h (/ a 5))<br />   (repeat (1+ c)<br />     (setq f (cons h f))<br />   )<br />   (while (> c 0)<br />     (setq d 0)<br />     (setq g (+ c c))<br />     (setq b c)<br />     (setq x nil)<br />     (while (> b 0)<br />       (setq d (* d b))<br />       (setq b (1- b))<br />       (setq d (+ d (* (car f) a)))<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 14<br />       (setq f (cdr f))<br />     )<br />     (setq s (+ e (/ d a)))<br />     (setq s (itoa s))<br />     (while (< (strlen s) 4)<br />       (setq s (strcat "0" s))<br />     )<br />     (setq P (strcat p s))<br />     (setq e (rem d a))<br />     (setq c (- c 14))<br />   )<br />   p<br /> )  <br /> ;;;以后给出完全注释版的<br />

今天重新优化了这个程序，使得速度提高了不少，而且编译成vlx文件，速度更是大步提升。
你若对此有趣，不妨对编译前后做一下比较。
在我的机器上，运行它不到0.5秒。
lisp文件：

vlx文件：


2011.11.1 -----更新在8楼

highflybir

原理介绍参见这篇帖子：http://blog.csdn.net/zyl910/article/details/1368387

详细注视的lisp文件在这里：
部分源代码：

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 />

ZZXXQQ

山巅一寺一壶酒……

carrot1983

死机啦。。。呵呵，机子配置不好，

dinosaurhxe

支持，

同时期待完全注释版，呵呵。

小菜123

我只用了1.28秒

zml84

很强悍，算法还没有看明白。

http://wg853785.gjjblog.com/archives/748463/

egos

还算比较快的  不用一秒时间