Code666 (代码贴、代码片段)

[C] 求n位圆周率π 超长圆周率计算 →→→→→进入此内容的聊天室

来自 , 2019-08-29, 写在 C, 查看 116 次.
URL http://www.code666.cn/view/8929c70f
  1. #include <stdio.h>
  2. #define L 10000 //求10000位PI值
  3. #define N L/4+1
  4.  
  5. // L 为位数，N是array长度
  6.  
  7. /*圆周率后的小数位数是无止境的，如何使用电脑来计算这无止境的小数是一些数学家与程式设计师所感兴趣的，在这边介绍一个公式配合 大数运算，可以计算指定位数的圆周率。
  8.  
  9. J.Marchin的圆周率公式：
  10. PI = [16/5 - 16 / (3*53) + 16 / (5*55) - 16 / (7*57) + ......] -  [4/239 - 4/(3*2393) + 4/(5*2395) - 4/(7*2397) + ......]
  11.  
  12. */
  13. void add ( int*, int*, int* );
  14. void sub ( int*, int*, int* );
  15. void div ( int*, int, int* );
  16.  
  17. int main ( void )
  18. {
  19.         int s[N+3] = {0};
  20.         int w[N+3] = {0};
  21.         int v[N+3] = {0};
  22.         int q[N+3] = {0};
  23.         int n = ( int ) ( L/1.39793 + 1 );
  24.         int k;
  25.  
  26.         w[0] = 16*5;
  27.         v[0] = 4*239;
  28.  
  29.         for ( k = 1; k <= n; k++ )
  30.         {
  31.                 // 套用公式
  32.                 div ( w, 25, w );
  33.                 div ( v, 239, v );
  34.                 div ( v, 239, v );
  35.                 sub ( w, v, q );
  36.                 div ( q, 2*k-1, q );
  37.  
  38.                 if ( k%2 ) // 奇数项
  39.                         add ( s, q, s );
  40.                 else    // 偶数项
  41.                         sub ( s, q, s );
  42.         }
  43.  
  44.         printf ( "%d.", s[0] );
  45.         for ( k = 1; k < N; k++ )
  46.                 printf ( "%04d", s[k] );
  47.         printf ( "\n" );
  48.         return 0;
  49. }
  50.  
  51. void add ( int *a, int *b, int *c )
  52. {
  53.         int i, carry = 0;
  54.  
  55.         for ( i = N+1; i >= 0; i-- )
  56.         {
  57.                 c[i] = a[i] + b[i] + carry;
  58.                 if ( c[i] < 10000 )
  59.                         carry = 0;
  60.                 else   // 进位
  61.                 {
  62.                         c[i] = c[i] - 10000;
  63.                         carry = 1;
  64.                 }
  65.         }
  66. }
  67.  
  68. void sub ( int *a, int *b, int *c )
  69. {
  70.         int i, borrow = 0;
  71.         for ( i = N+1; i >= 0; i-- )
  72.         {
  73.                 c[i] = a[i] - b[i] - borrow;
  74.                 if ( c[i] >= 0 )
  75.                         borrow = 0;
  76.                 else   // 借位
  77.                 {
  78.                         c[i] = c[i] + 10000;
  79.                         borrow = 1;
  80.                 }
  81.         }
  82. }
  83.  
  84. void div ( int *a, int b, int *c ) // b 为除数
  85. {
  86.         int i, tmp, remain = 0;
  87.         for ( i = 0; i <= N+1; i++ )
  88.         {
  89.                 tmp = a[i] + remain;
  90.                 c[i] = tmp / b;
  91.                 remain = ( tmp % b ) * 10000;
  92.         }
  93. }
  94.  

回复 "求n位圆周率π 超长圆周率计算"

这儿你可以回复上面这条便签

captcha

贴代码片段 - Code666