quickjs-tart

quickjs-based runtime for wallet-core logic
Log | Files | Refs | README | LICENSE

pi_bigint.js (2830B)

      1 /*
      2  * PI computation in Javascript using the BigInt type
      3  */
      4 "use strict";
      5 
      6 /* return floor(log2(a)) for a > 0 and 0 for a = 0 */
      7 function floor_log2(a)
      8 {
      9     var k_max, a1, k, i;
     10     k_max = 0n;
     11     while ((a >> (2n ** k_max)) != 0n) {
     12         k_max++;
     13     }
     14     k = 0n;
     15     a1 = a;
     16     for(i = k_max - 1n; i >= 0n; i--) {
     17         a1 = a >> (2n ** i);
     18         if (a1 != 0n) {
     19             a = a1;
     20             k |= (1n << i);
     21         }
     22     }
     23     return k;
     24 }
     25 
     26 /* return ceil(log2(a)) for a > 0 */
     27 function ceil_log2(a)
     28 {
     29     return floor_log2(a - 1n) + 1n;
     30 }
     31 
     32 /* return floor(sqrt(a)) (not efficient but simple) */
     33 function int_sqrt(a)
     34 {
     35     var l, u, s;
     36     if (a == 0n)
     37         return a;
     38     l = ceil_log2(a);
     39     u = 1n << ((l + 1n) / 2n);
     40     /* u >= floor(sqrt(a)) */
     41     for(;;) {
     42         s = u;
     43         u = ((a / s) + s) / 2n;
     44         if (u >= s)
     45             break;
     46     }
     47     return s;
     48 }
     49 
     50 /* return pi * 2**prec */
     51 function calc_pi(prec) {
     52     const CHUD_A = 13591409n;
     53     const CHUD_B = 545140134n;
     54     const CHUD_C = 640320n;
     55     const CHUD_C3 = 10939058860032000n; /* C^3/24 */
     56     const CHUD_BITS_PER_TERM = 47.11041313821584202247; /* log2(C/12)*3 */
     57 
     58     /* return [P, Q, G] */
     59     function chud_bs(a, b, need_G) {
     60         var c, P, Q, G, P1, Q1, G1, P2, Q2, G2;
     61         if (a == (b - 1n)) {
     62             G = (2n * b - 1n) * (6n * b - 1n) * (6n * b - 5n);
     63             P = G * (CHUD_B * b + CHUD_A);
     64             if (b & 1n)
     65                 P = -P;
     66             Q = b * b * b * CHUD_C3;
     67         } else {
     68             c = (a + b) >> 1n;
     69             [P1, Q1, G1] = chud_bs(a, c, true);
     70             [P2, Q2, G2] = chud_bs(c, b, need_G);
     71             P = P1 * Q2 + P2 * G1;
     72             Q = Q1 * Q2;
     73             if (need_G)
     74                 G = G1 * G2;
     75             else
     76                 G = 0n;
     77         }
     78         return [P, Q, G];
     79     }
     80 
     81     var n, P, Q, G;
     82     /* number of serie terms */
     83     n = BigInt(Math.ceil(Number(prec) / CHUD_BITS_PER_TERM)) + 10n;
     84     [P, Q, G] = chud_bs(0n, n, false);
     85     Q = (CHUD_C / 12n) * (Q << prec) / (P + Q * CHUD_A);
     86     G = int_sqrt(CHUD_C << (2n * prec));
     87     return (Q * G) >> prec;
     88 }
     89 
     90 function main(args) {
     91     var r, n_digits, n_bits, out;
     92     if (args.length < 1) {
     93         print("usage: pi n_digits");
     94         return;
     95     }
     96     n_digits = args[0] | 0;
     97 
     98     /* we add more bits to reduce the probability of bad rounding for
     99       the last digits */
    100     n_bits = BigInt(Math.ceil(n_digits * Math.log2(10))) + 32n;
    101     r = calc_pi(n_bits);
    102     r = ((10n ** BigInt(n_digits)) * r) >> n_bits;
    103     out = r.toString();
    104     print(out[0] + "." + out.slice(1));
    105 }
    106 
    107 var args;
    108 if (typeof scriptArgs != "undefined") {
    109     args = scriptArgs;
    110     args.shift();
    111 } else if (typeof arguments != "undefined") {
    112     args = arguments;
    113 } else {
    114     /* default: 1000 digits */
    115     args=[1000];
    116 }
    117 
    118 main(args);