C 中的莱布尼茨函数
Leibniz function in c
我正在寻找 C 中 PI 的 Leibniz 近似值。
我已经测试了这个功能,但它似乎不起作用:
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: %f\n", leibnizPi(1000));
return 0;
}
OP 的代码有趣地使用 double
数学进行计算,然后在可能很长的循环之后 returns a float
。而不是截断为浮点数。推荐给return一个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("%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
我正在寻找 C 中 PI 的 Leibniz 近似值。
我已经测试了这个功能,但它似乎不起作用:
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: %f\n", leibnizPi(1000));
return 0;
}
OP 的代码有趣地使用 double
数学进行计算,然后在可能很长的循环之后 returns a float
。而不是截断为浮点数。推荐给return一个double
。使用 float
,预计结果不会比 223.
下面显示了在 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