Leibniz's formula for Ο€ / 4

You are asked to print the summation of the Leibniz formula up to the n-th term of the series, correct up to 15 decimal places. In calculus, the Leibniz formula for Ο€ is given as follows: 1 - 1/3 + 1/5 -1/7 + ... = Ο€ / 4

This is my code

#include<stdio.h>
#include<math.h>
int main()
{
    int n,i;
    long double s=0;
    scanf("%d",&n);
    for(i=0;i<n;i++){
        s+=(long double)pow(-1,i)/(2*i+1);
    }
    printf("%Lf\n",s);
    return 0;
}

      

Can anyone tell me why I cannot achieve precision to fifteenth decimal place? My goal is not to print the pi / 4 value, I just need to print the summation for a given n

+3


source to share


2 answers


Q: why ... achieve precision to 15th decimal place?
A: To display 15 decimal places after using the decimal point, use "%0.15f"

. To calculate 15 decimal places, the convergence n

must be at least very large.

As @ user3386109 mentioned, "the error in the result is limited to 1 / (2n + 1)", so it will take about 5e14 iterations. (Rough estimate: 10 days on my PC.) Since typical double

has an accuracy of about 1 part in pow (2.53) or 1 in 9e15, the calculation limits are reached double

. The code below compares the order of calculation to reduce the error, but at best the error will still be at least 0.5 part in 9e15.

Since the members of the series fluctuate around the limit, stopping after n

iterations, it is possible to add a final 1/2 iteration of the next iteration. this will have about 1 bit of precision.

As others have noted, there are other methods for computing Ο€ that converge faster.




Updated for good observation by @ user3386109.

By summarizing terms, the code can summarize them in different orders. The methods below illustrate that a more moderate, more stable result is achieved earlier when the small terms are summed up first. I only expect 1 or 2 bits of the best answer at best.

This was re-done with help float

as the slight improvement doesn't show until the last few bits are stable. Even double

with this slowly converging series, it will take too long.

//Leibniz formula for pi/4
typedef float fp;

fp LeibnizForward(unsigned n) {
  volatile fp sum = 0.0;
  fp sign = 1.0;
  unsigned i = 1;
  while (n-- > 0) {
    sum += sign / i;
    sign = -sign;
    i = (i + 2);
  }
  return sum;
}

fp LeibnizReverse(unsigned n) {
  volatile fp sum = 0.0;
  fp sign = 1.0;
  unsigned i = 2 * n - 1;
  if (n % 2 == 0)
    sign = -sign;
  while (n-- > 0) {
    sum += sign / i;
    sign = -sign;
    i = (i - 2);
  }
  return sum;
}

void PiTest(unsigned n) {
  printf("%u\n", n);
  static const fp pic = 3.1415926535897932384626433832795;
  const char *format = "%s %0.9f\n";
  printf(format, "pi-", nextafterf(pic,0));
  printf(format, "pi ", pic);
  printf(format, "pi+", nextafterf(pic,4));
  fp pif = LeibnizForward(n) * 4;
  printf(format, "pif", pif);
  fflush(stdout);
  fp pir = LeibnizReverse(n) * 4;
  printf(format, "pir", pir);
  fflush(stdout);
}

int main(void) {
  PiTest(0);
  PiTest(1);
  PiTest(10);
  PiTest(100);
  PiTest(1000);
  PiTest(10000);
  PiTest(100000);
  PiTest(1000000);
  PiTest(10000000);
  PiTest(100000000);

  return 0;
}

0
pi- 3.141592503
pi  3.141592741
pi+ 3.141592979
pif 0.000000000
pir 0.000000000
1
pi- 3.141592503
pi  3.141592741
pi+ 3.141592979
pif 4.000000000
pir 4.000000000
10
pi- 3.141592503
pi  3.141592741
pi+ 3.141592979
pif 3.041839600
pir 3.041839600
25
pi- 3.141592503
pi  3.141592741
pi+ 3.141592979
pif 3.181576490
pir 3.181576729
100
pi- 3.141592503
pi  3.141592741
pi+ 3.141592979
pif 3.131592512
pir 3.131592751
1000
pi- 3.14 1592503
pi  3.14 1592741
pi+ 3.14 1592979
pif 3.14 0592575
pir 3.14 0592575
10000
pi- 3.141 592503
pi  3.141 592741
pi+ 3.141 592979
pif 3.141 498566
pir 3.141 492605
100000
pi- 3.1415 92503
pi  3.1415 92741
pi+ 3.1415 92979
pif 3.1415 85827
pir 3.1415 82489
1000000
pi- 3.14159 2503
pi  3.14159 2741
pi+ 3.14159 2979
pif 3.14159 5364
pir 3.14159 1549
10000000
pi- 3.14159 2503 previous float
pi  3.14159 2741 machine float pi
pi+ 3.14159 2979 next float
pif 3.14159 6794
pir 3.14159 2503
100000000 Additional iterations do not improve the result.
pi- 3.14159 2503
pi  3.14159 2741 
pi+ 3.14159 2979
pif 3.14159 6794
pir 3.14159 2503
1000000000
pi- 3.14159 2503
pi  3.14159 2741
pi+ 3.14159 2979
pif 3.14159 6794
pir 3.14159 2503

      

+4


source


Modify your operator printf

to print more values ​​after the decimal point.

As an example using

printf("%30.28Lf\n",s);

      



will print a number up to 28 places after the decimal point. With n=25

I got the following output.

0.7953941713587578038252220991

      

0


source







All Articles