Function bernoulliRational

src/compute-engine/numerics/bernoulli.ts:53–80

(n: number)

 * exact path caps |argument| at 100, i.e. at most B₁₀₁).
 */
export function bernoulliRational(n: number): [bigint, bigint] {
  if (!Number.isInteger(n) || n < 0)
    throw new RangeError(`bernoulliRational: invalid index ${n}`);

  for (let m = BERNOULLI.length; m <= n; m++) {
    // Odd m > 1: B_m = 0
    if (m % 2 === 1) {
      BERNOULLI.push([0n, 1n]);
      continue;
    }

    // B_m = -1/(m+1) * sum_{k=0}^{m-1} C(m+1, k) * B_k
    const mp1 = BigInt(m + 1);
    let sumNum = 0n;
    let sumDen = 1n;
    let binom = 1n; // C(m+1, 0) = 1
    for (let k = 0; k < m; k++) {
      if (k > 0) binom = (binom * (mp1 - BigInt(k) + 1n)) / BigInt(k);
      const [bkNum, bkDen] = BERNOULLI[k];
      if (bkNum === 0n) continue;
      sumNum = sumNum * bkDen + binom * bkNum * sumDen;
      sumDen = sumDen * bkDen;
    }
    BERNOULLI.push(reduce(-sumNum, mp1 * sumDen));
  }

  return BERNOULLI[n];
}

/**
 * The exact rational c such that ζ(2k) = c·π^{2k}, for integer k ≥ 1.