#include #include #include #include // Calculate arctan(1/x) using Taylor series // arctan(1/x) = 1/x - 1/(3*x^3) + 1/(5*x^5) - ... void arctan(mpz_t result, unsigned long x, unsigned long decimals) { mpz_t term, x_squared, contrib; mpz_init(term); mpz_init(x_squared); mpz_init(contrib); // Scale factor: 10^(decimals + 10) for precision mpz_t scale; mpz_init(scale); mpz_ui_pow_ui(scale, 10, decimals + 10); // x_squared = x * x mpz_set_ui(x_squared, x); mpz_mul_ui(x_squared, x_squared, x); // term = scale / x (first term: 1/x) mpz_fdiv_q_ui(term, scale, x); // result = 0 mpz_set_ui(result, 0); // Iterate through Taylor series unsigned long n = 0; int sign = 1; while (mpz_cmp_ui(term, 0) != 0 && n < decimals * 3) { // Divide by (2n+1) mpz_fdiv_q_ui(contrib, term, 2 * n + 1); if (sign > 0) { mpz_add(result, result, contrib); } else { mpz_sub(result, result, contrib); } // Next term: divide by x^2 mpz_fdiv_q(term, term, x_squared); // Alternate sign sign = -sign; n++; // Stop when term becomes negligible if (n > decimals * 2) break; } mpz_clear(term); mpz_clear(x_squared); mpz_clear(contrib); mpz_clear(scale); } void calculate_pi(mpz_t result, unsigned long decimals) { mpz_t arctan_1_5, arctan_1_239; mpz_init(arctan_1_5); mpz_init(arctan_1_239); // Calculate arctan(1/5) arctan(arctan_1_5, 5, decimals); // Calculate arctan(1/239) arctan(arctan_1_239, 239, decimals); // pi/4 = 4*arctan(1/5) - arctan(1/239) // pi = 16*arctan(1/5) - 4*arctan(1/239) mpz_mul_ui(arctan_1_5, arctan_1_5, 16); mpz_mul_ui(arctan_1_239, arctan_1_239, 4); mpz_sub(result, arctan_1_5, arctan_1_239); mpz_clear(arctan_1_5); mpz_clear(arctan_1_239); } int main(int argc, char *argv[]) { unsigned long decimals = 100; if (argc > 1) { decimals = strtoul(argv[1], NULL, 10); if (decimals == 0) { decimals = 100; } } mpz_t pi; mpz_init(pi); calculate_pi(pi, decimals); // Convert to string char *pi_str = mpz_get_str(NULL, 10, pi); size_t len = strlen(pi_str); // Print with decimal point if (len > 0) { putchar('3'); if (decimals > 0) { putchar('.'); // Print decimals, skipping the leading "3" size_t start = 1; size_t end = start + decimals; if (end > len) end = len; for (size_t i = start; i < end && i < len; i++) { putchar(pi_str[i]); } } putchar('\n'); } free(pi_str); mpz_clear(pi); return 0; }