c语言中的莱布尼茨函数答案

【问题标题】：Leibniz function in cc语言中的莱布尼茨函数
【发布时间】：2017-10-17 19:29:16
【问题描述】：

我正在寻找 C 中 PI 的莱布尼茨近似值。

我已经测试过这个功能，但它似乎不起作用：

float leibnizPi(float n){
    double pi=1.0;
    int i;
    int N;
    for (i=3, N=2*n+1; i<=N; i+=2)
        pi += ((i&2) ? -1.0 : 1.0) / i;
    return 4*pi;
}

【问题讨论】：

我认为 & 应该是 % 用于模运算
你的输出是什么？显示minimal reproducible example。
@bhow i 总是奇怪的，所以你的建议会产生一个恒定的。
@没错，我的错。看来我读得不够仔细。
i&2 部分是正确的。这对于 i=3, 7, 11,... 但不是 5, 9, 13.... 是正确的。为什么 n 是 float 当它显然只是一个计数？否则它看起来很好。为什么认为有错误？换句话说，你给它什么输入，你得到什么输出，你认为它们应该是什么？

标签： c math pi approximation

【解决方案1】：

我已经测试了你的功能，它似乎工作正常。对于n=1000，我得到3.142592 的结果。

我是如何测试这个功能的：

#include <stdio.h>

float leibnizPi(float n) {
    double pi=1.0;
    int i;
    int N;
    for (i=3, N=2*n+1; i<=N; i+=2)
        pi += ((i&2) ? -1.0 : 1.0) / i;
    return 4*pi;
}

int main(void)
{
    printf("Pi: %f\n", leibnizPi(1000));
    return 0;
}

【讨论】：

我认为这是对的。我用n=1000000 试过了，它似乎收敛了，只是需要很多术语。

【解决方案2】：

有趣的 OP 代码使用 double 数学进行计算，然后在可能很长的循环之后返回 float。而不是截断为浮点数。建议返回double。使用float，预计结果不会比可能 2²³ 中的 1 个部分更准确。

下图显示了在n > about 100,000 和最终大约 600,000 之后变得明显的不必要的影响。

float leibnizPif(float n){
    double pi=1.0;
    int i;
    int N;
    for (i=3, N=2*n+1; i<=N; i+=2)
        pi += ((i&2) ? -1.0 : 1.0) / i;
    return 4*pi;
}

double leibnizPid(float n){
    double pi=1.0;
    int i;
    int N;
    for (i=3, N=2*n+1; i<=N; i+=2)
        pi += ((i&2) ? -1.0 : 1.0) / i;
    return 4*pi;
}

int main() {

  double pi = acos(-1);
  float f1 = pi;
  float f2 = nextafterf(pi, 0);
  printf("%e\n", f1-f2);
  printf("                     pi  %.12e\n", pi);
  for (unsigned i=1; i<29; i++) {
    unsigned n = 1u << i;
    float f = leibnizPif(n);
    double d = leibnizPid(n);
    printf("%9u f % .5e %.12e   ", n, (f-pi)/pi, f);
    printf("d % .5e %.12e\n",  (d-pi)/pi, d);

  }
  printf("                     pi  %.12e\n", acos(-1));
  return 0;
}

输出

            error                pi
                     pi  3.141592653590e+00
        2 f  1.03474e-01 3.466666698456e+00   d  1.03474e-01 3.466666666667e+00
        4 f  6.30540e-02 3.339682579041e+00   d  6.30540e-02 3.339682539683e+00
        8 f  3.52602e-02 3.252365827560e+00   d  3.52602e-02 3.252365934719e+00
       16 f  1.87080e-02 3.200365543365e+00   d  1.87080e-02 3.200365515410e+00
       32 f  9.64357e-03 3.171888828278e+00   d  9.64354e-03 3.171888735237e+00
       64 f  4.89682e-03 3.156976461411e+00   d  4.89679e-03 3.156976358911e+00
      128 f  2.46747e-03 3.149344444275e+00   d  2.46748e-03 3.149344475125e+00
      256 f  1.23857e-03 3.145483732224e+00   d  1.23856e-03 3.145483689446e+00
      512 f  6.20513e-04 3.143542051315e+00   d  6.20487e-04 3.143541969477e+00
     1024 f  3.10574e-04 3.142568349838e+00   d  3.10546e-04 3.142568263114e+00
     2048 f  1.55377e-04 3.142080783844e+00   d  1.55349e-04 3.142080696509e+00
     4096 f  7.76643e-05 3.141836643219e+00   d  7.76934e-05 3.141836734621e+00
     8192 f  3.88840e-05 3.141714811325e+00   d  3.88514e-05 3.141714709002e+00
    16384 f  1.94559e-05 3.141653776169e+00   d  1.94269e-05 3.141653685021e+00
    32768 f  9.74187e-06 3.141623258591e+00   d  9.71375e-06 3.141623170237e+00
    65536 f  4.88485e-06 3.141607999802e+00   d  4.85695e-06 3.141607912146e+00
   131072 f  2.45634e-06 3.141600370407e+00   d  2.42849e-06 3.141600282926e+00
   262144 f  1.24208e-06 3.141596555710e+00   d  1.21425e-06 3.141596468272e+00
   524288 f  6.34955e-07 3.141594648361e+00   d  6.07127e-07 3.141594560935e+00
  1048576 f  3.31391e-07 3.141593694687e+00   d  3.03564e-07 3.141593607263e+00
  2097152 f  1.79610e-07 3.141593217850e+00   d  1.51782e-07 3.141593130427e+00
  4194304 f  1.03719e-07 3.141592979431e+00   d  7.58910e-08 3.141592892008e+00
  // No further improvement in float as precision is reached.
  8388608 f  2.78275e-08 3.141592741013e+00   d  3.79455e-08 3.141592772799e+00
 16777216 f  2.78275e-08 3.141592741013e+00   d  1.89728e-08 3.141592713194e+00
 33554432 f  2.78275e-08 3.141592741013e+00   d  9.48662e-09 3.141592683393e+00
 67108864 f  2.78275e-08 3.141592741013e+00   d  4.74328e-09 3.141592668491e+00
134217728 f  2.78275e-08 3.141592741013e+00   d  2.37136e-09 3.141592661040e+00
268435456 f  2.78275e-08 3.141592741013e+00   d  1.18567e-09 3.141592657315e+00
536870912 f  2.78275e-08 3.141592741013e+00   d  5.92798e-10 3.141592655452e+00
                     pi  3.141592653590e+00
                                              double keeps getting better

修改后的代码输出

// with changes to use 64-bit integer math for the loops, 
// the result continue to get better returning double.
// It just takes longer and longer. 
1073741824 f  2.78275e-08 3.141592741013e+00   d  2.95864e-10 3.141592654519e+00
2147483648 f  2.78275e-08 3.141592741013e+00   d  1.47599e-10 3.141592654053e+00
4294967296 f  2.78275e-08 3.141592741013e+00   d  7.34476e-11 3.141592653821e+00
8589934592 f  2.78275e-08 3.141592741013e+00   d  3.63406e-11 3.141592653704e+00

【讨论】：