Sort of like halfpiarea.c, but in bc source. なんとなく bc 型ソースコードにして halfpiarea.c のように動く 。 Compute the value of pi by using numerical methods to integrate the area of a half circle. 半円の面積を数値積分を使って計算して π の値を計算する。

Format
Plain text
Post date
2016-10-04 22:01
Publication Period
Unlimited
1. /* halfpiArea_by_m.bc
2. // Computing the value of pi by integrating the area of a (half of a) circle.
3. // by Joel Matthew Rees, 4 October 2016
4. // Copyright Joel Matthew Rees
5. //
6. // Fair use acknowledged.
7. //
8. // All rights to this expression of the method reserved.
9. //
10. // (Starting from scratch isn't that hard,
11. // and you'll like the results better.)
12. //
13. // Set m to count of areas to take;
14. // set d to 1 to center the height of rectangle on the circle;
15. // set scale if you want to change accuracy in the low bits;
16. // set q to non-zero if you don't want to show the loop output:
17. //
18. // echo "m=10000;q=1;scale=41;" | cat - halfpiArea_by_m.bc | bc -l
19. //
20. // Time it with
21. //
22. // time echo "m=10000;q=1;scale=41;" | cat - halfpiArea_by_m.bc | bc -l
23. //
24. // Note that the errors with odd counts (m=11, etc.) are not relevant
25. // with accuracies obtained at counts less than 1000000.
26. */
27. if ( m < 1 ) { m = 11; }
28. s=0;
29. w = 1/m;
30. if ( d == 1 ) {
31. h = w / 2;
32. f = 1 - h;
33. n = -1 * m;
34. for (x = -1 + h; x <= f; x = x + w ) {
35. n = n + 1;
36. y = sqrt( 1 - x * x );
37. s = s + y * w;
38. if ( ( q == 0 ) && ( m <= 1000 ) ) {
39. print n, ": (", x, ", ", y, "): ", s, "\n";
40. }
41. }
42. }
43. if ( d != 1 ) { /* Some ancient versions of bc don't have else. I think. */
44. for (n = -m; n <= m; n = n + 1 ) {
45. x = n * w;
46. y = sqrt( 1 - x * x );
47. s = s + y * w;
48. if ( ( q == 0 ) && ( m <= 1000 ) ) {
49. print n, ": (", x, ", ", y, "): ", s, "\n";
50. }
51. }
52. }
53. print n, ": (", x, ", ", y, "): ", s, "\n";
54. s * 2;