docs: function findx0BigA(a, p, q) {
        var w, s, D, B, u, z, zbar, logpg;
  
        // snTerm -- see formula 34
        logpg = Math.log(p * gamma(a + 1));
        function snTerm(x) {
           return function(n, v) {
              return n === 0 ? 1 : v * x / (a + n);
           };
        }
        // See formula 34, where `n` is off by 1
        function f(x, n) {
           return Math.exp((logpg + x - Math.log(series(snTerm(x), n + 1))) / a);
        }
        B = q * gamma(a);
        D = Math.max(2, a * (a - 1));
        s = funS(p, q);
        w = a + s * Math.sqrt(a) + (s * s - 1) / 3 +
           (s * s * s - 7 * s) / (36 * Math.sqrt(a)) -
           (3 * s * s * s * s + 7 * s * s - 16) / (810 * a) +
           (9 * s * s * s * s * s + 256 * s * s * s - 433 * s) /
              (38880 * a * Math.sqrt(a));
        if (a >= 500 && Math.abs(1 - w / a) < 1e-6) { return w; }
        if (p > 0.5) {
           if (w < 3 * a) { return w; }
           if (B < Math.pow(10, -D)) { return findx0TinyB(a, -Math.log(B)); }
           u = -Math.log(B) + (a - 1) * Math.log(w) -
                 Math.log(1 + (1 - a) / (1 + w));
           // formula 33
           return -Math.log(B) + (a - 1) * Math.log(u) -
                 Math.log(1 + (1 - a) / (1 + u));
        }
        // handle p <= 0.5
        z = w > 0.15 * (a + 1) ? w : f(f(f(f(w, 1), 3), 3), 4);
        if (z <= 0.01 * (a + 1) || z > 0.7 * (a + 1)) { return z; }
        /* eslint-disable no-extra-parens */
        return (function() {
           // find N for formula 36
           var N, term;
  
           N = 0;
           term = 1;
           do {
              N += 1;
              term = snTerm(z)(N, term);
           } while (term > 1e-4);
           zbar = f(z, N + 1);
           return zbar * (1 - (a * Math.log(zbar) - zbar - logpg +
                    Math.log(series(snTerm(z), N + 1))) / (a - zbar));
        }());
        /* eslint-enable no-extra-parens */
     }