c中的莱布尼茨函数



我正在寻找 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;
}

我已经测试了你的功能,它似乎工作正常。对于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: %fn", leibnizPi(1000));
    return 0;
}

OP的代码很有趣地使用double数学进行计算,然后在潜在的长循环后返回float。 而不是截断为浮点数。 建议返回double。 对于float,结果预计不会比 223 中的 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("%en", f1-f2);
  printf("                     pi  %.12en", 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 %.12en",  (d-pi)/pi, d);
  }
  printf("                     pi  %.12en", 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

最新更新