1/*
2 * PI computation in Javascript using the QuickJS bigdecimal type
3 * (decimal floating point)
4 */
5"use strict";
6
7/* compute PI with a precision of 'prec' digits */
8function calc_pi(prec) {
9    const CHUD_A = 13591409m;
10    const CHUD_B = 545140134m;
11    const CHUD_C = 640320m;
12    const CHUD_C3 = 10939058860032000m; /* C^3/24 */
13    const CHUD_DIGITS_PER_TERM = 14.18164746272548; /* log10(C/12)*3 */
14
15    /* return [P, Q, G] */
16    function chud_bs(a, b, need_G) {
17        var c, P, Q, G, P1, Q1, G1, P2, Q2, G2, b1;
18        if (a == (b - 1n)) {
19            b1 = BigDecimal(b);
20            G = (2m * b1 - 1m) * (6m * b1 - 1m) * (6m * b1 - 5m);
21            P = G * (CHUD_B * b1 + CHUD_A);
22            if (b & 1n)
23                P = -P;
24            G = G;
25            Q = b1 * b1 * b1 * CHUD_C3;
26        } else {
27            c = (a + b) >> 1n;
28            [P1, Q1, G1] = chud_bs(a, c, true);
29            [P2, Q2, G2] = chud_bs(c, b, need_G);
30            P = P1 * Q2 + P2 * G1;
31            Q = Q1 * Q2;
32            if (need_G)
33                G = G1 * G2;
34            else
35                G = 0m;
36        }
37        return [P, Q, G];
38    }
39
40    var n, P, Q, G;
41    /* number of serie terms */
42    n = BigInt(Math.ceil(prec / CHUD_DIGITS_PER_TERM)) + 10n;
43    [P, Q, G] = chud_bs(0n, n, false);
44    Q = BigDecimal.div(Q, (P + Q * CHUD_A),
45                       { roundingMode: "half-even",
46                         maximumSignificantDigits: prec });
47    G = (CHUD_C / 12m) * BigDecimal.sqrt(CHUD_C,
48                                         { roundingMode: "half-even",
49                                           maximumSignificantDigits: prec });
50    return Q * G;
51}
52
53(function() {
54    var r, n_digits, n_bits;
55    if (typeof scriptArgs != "undefined") {
56        if (scriptArgs.length < 2) {
57            print("usage: pi n_digits");
58            return;
59        }
60        n_digits = scriptArgs[1] | 0;
61    } else {
62        n_digits = 1000;
63    }
64    /* we add more digits to reduce the probability of bad rounding for
65       the last digits */
66    r = calc_pi(n_digits + 20);
67    print(r.toFixed(n_digits, "down"));
68})();
69