#include #include #include typedef long double Number; static Number circle(Number x) { return sqrtl(1.0 - x * x); } static Number monteCarloIntegration(void) { int amount = 512 * 512; int in = 0; for(int i = 0; i < amount; i++) { Number x = rand() / (Number)RAND_MAX; Number y = rand() / (Number)RAND_MAX; in += (x * x + y * y < 1); } return 4.0 * in / amount; } static Number simpson(Number a, Number b) { return ((b - a) / 6.0) * (circle(a) + 4.0 * circle((a + b) * 0.5) + circle(b)); } static Number simpsonIntegration(void) { Number sum = 0.0; for(int i = 0; i < 100; i++) { sum += simpson(i / 100.0, (i + 1) / 100.0); } return 4.0 * sum; } static Number lineSum(void) { Number sum = 0.0; const int max = 1000000; for(int i = 0; i < max; i++) { Number x1 = (Number)i / max; Number x2 = (i + 1.0) / max; Number diffX = x2 - x1; Number diffY = circle(x2) - circle(x1); sum += sqrtl(diffX * diffX + diffY * diffY); } return sum * 2.0; } static Number lineRecursion(void) { Number x1 = 0.0; Number y1 = 1.0; Number x2 = 1.0; Number y2 = 0.0; Number factor = 2.0; for(int i = 0; i < 26; i++) { Number midX = (x1 + x2) * 0.5; Number midY = (y1 + y2) * 0.5; Number length = sqrtl(midX * midX + midY * midY); midX /= length; midY /= length; x2 = midX; y2 = midY; factor *= 2; } Number diffX = x2 - x1; Number diffY = y2 - y1; return sqrtl(diffX * diffX + diffY * diffY) * factor; } static void print(Number d) { Number diff = d - M_PI; printf("%.20Lf %.20Lf\n", d, diff < 0.0 ? -diff : diff); } #define STR(x) STR2(x) #define STR2(x) #x int main(void) { print(monteCarloIntegration()); print(simpsonIntegration()); print(lineSum()); print((double)lineRecursion()); puts(STR(M_PI)); return 0; }