Pastebin: halfpiArea_by_m.bc

Format

Plain text

Post date

2016-10-04 22:01

Publication Period

Unlimited

 /* halfpiArea_by_m.bc
 // Computing the value of pi by integrating the area of a (half of a) circle.
 // by Joel Matthew Rees, 4 October 2016
 // Copyright Joel Matthew Rees
 //
 // Fair use acknowledged.
 //
 // All rights to this expression of the method reserved.
 //
 // (Starting from scratch isn't that hard,
 // and you'll like the results better.)
 //
 // Set m to count of areas to take;
 // set d to 1 to center the height of rectangle on the circle;
 // set scale if you want to change accuracy in the low bits;
 // set q to non-zero if you don't want to show the loop output:
 //
 // echo "m=10000;q=1;scale=41;" | cat - halfpiArea_by_m.bc | bc -l
 //
 // Time it with 
 //
 // time echo "m=10000;q=1;scale=41;" | cat - halfpiArea_by_m.bc | bc -l
 //
 // Note that the errors with odd counts (m=11, etc.) are not relevant 
 // with accuracies obtained at counts less than 1000000.
 */
 
 if ( m < 1 ) { m = 11; }
 
 s=0;
 w = 1/m;
 if ( d == 1 ) {
    h = w / 2;
    f = 1 - h;
    n = -1 * m;
    for (x = -1 + h; x <= f; x = x + w ) {
       n = n + 1;
       y = sqrt( 1 - x * x );
       s = s + y * w;
       if ( ( q == 0 ) && ( m <= 1000 ) ) {
          print n, ": (", x, ", ", y, "): ", s, "\n";
       }
    }
 }
 if ( d != 1 ) { /* Some ancient versions of bc don't have else. I think. */
    for (n = -m; n <= m; n = n + 1 ) {
       x = n * w;
       y = sqrt( 1 - x * x );
       s = s + y * w;
       if ( ( q == 0 ) && ( m <= 1000 ) ) {
          print n, ": (", x, ", ", y, "): ", s, "\n";
       }
    }
 }
 print n, ": (", x, ", ", y, "): ", s, "\n";
 s * 2;