Function ddAdd

src/compute-engine/numerics/special-functions.ts:2208–2213

(a: DD, b: DD)

  const signA = n % 2 === 0 ? -1 : 1; // (−1)^{n−1}
  result += (signA * factorial(n - 1)) / Math.pow(z, n);
  result += (signA * factN) / (2 * Math.pow(z, n + 1));

  // Bernoulli tail. The k-th factor fₖ = (2k+n−1)!/((2k)!·z^{2k+n}) is built
  // incrementally: f₁ = (n+1)!/(2!·z^{n+2}),
  //   f_{k+1} = fₖ·(2k+n)(2k+n+1)/((2k+1)(2k+2)·z²).
  // (A previous version used n(n+1)⋯(n+2k−1)/(2k)! — a factor (n−1)! too
  // small; only n ≤ 2 was unaffected since 1! = 0! = 1.)
  let f = (factN * (n + 1)) / (2 * Math.pow(z, n + 2));
  let prevAbs = Infinity;