bernoulli.h
#include <flint/bernoulli.h>This module computes the Bernoulli numbers, defined by the exponential generating function
as exact fractions: B_0 = 1, B_1 = -1/2, B_2 = 1/6, and B_n = 0 for odd n \ge 3. It provides isolated values (from the Riemann zeta function, or for very large n combined with Harvey's multimodular algorithm), B_n modulo a prime, ranges of consecutive values generated in reverse order, and a cache used by the special functions of arb and acb (gamma function, zeta function, Euler–Maclaurin summation) that need many Bernoulli numbers repeatedly. The arith functions arith_bernoulli_number and arith_bernoulli_number_vec are wrappers and alternatives.
All algorithms are described in Bernoulli numbers, and Isolated Bernoulli numbers explains the choice of method for one value.
AlgorithmBernoulli numbersComputes exact Bernoulli numbers from the zeta function (numerically, with the denominator from von Staudt–Clausen), from Harvey's multimodular algorithm, by a reverse zeta iterator with a table of powers for ranges, and by Ramanujan's recursion or power series modulo primes for short vectors.Isolated Bernoulli numbers#
void bernoulli_fmpq_ui(fmpq_t b, ulong n)void _bernoulli_fmpq_ui(fmpz_t num, fmpz_t den, ulong n)Set the result to B_n in lowest terms with a positive denominator. The value is copied from the cache if it is there; the cache is not extended. Otherwise odd n and n < 18000 use _bernoulli_fmpq_ui_zeta, and larger even n use _bernoulli_fmpq_ui_multi_mod with the default tuning.
void _bernoulli_fmpq_ui_zeta(fmpz_t num, fmpz_t den, ulong n)void _bernoulli_fmpq_ui_multi_mod(fmpz_t num, fmpz_t den, ulong n, double alpha)Set num and den to the numerator and denominator of B_n in lowest terms, by two methods.
Returns B_n \bmod p for a prime p < 2^{\mathtt{FLINT\_D\_BITS}} (2^{53} on 64-bit machines), computed with Harvey's algorithm [Har2010] in time O(p \log^{O(1)} p), independent of n. If p divides the denominator of B_n (that is, n \ge 2 even with (p - 1) \mid n, or n = 1 and p = 2), returns UWORD_MAX.
Ranges of Bernoulli numbers#
An iterator that generates the even-indexed Bernoulli numbers B_n, B_{n-2}, B_{n-4}, \ldots, B_0 in decreasing order, without caching. An array of length one of bernoulli_rev_struct.
void bernoulli_rev_init(bernoulli_rev_t iter, ulong n)Initialises iter so that the first call to bernoulli_rev_next returns B_n. An odd n is replaced by n - 1.
void bernoulli_rev_next(fmpz_t numer, fmpz_t denom, bernoulli_rev_t iter)Sets numer and denom to the numerator and denominator of B_k in lowest terms, where B_k is the next number of the sequence, and advances iter to B_{k-2}. Must not be called after B_0 has been returned.
void bernoulli_rev_clear(bernoulli_rev_t iter)Frees the memory used by iter.
void bernoulli_fmpq_vec_no_cache(fmpq * res, ulong a, slong num)Sets res to the num consecutive Bernoulli numbers B_a, B_{a+1}, \ldots, B_{a + \mathit{num} - 1}, without reading or writing the cache. To reuse a vector, use bernoulli_cache_compute and bernoulli_fmpq_ui instead. Throws if a > 2^{31} or \mathit{num} > 10^9.
Cache#
The cache holds B_0, \ldots, B_{N-1} as an array of fmpq. It is thread-local if FLINT is built with thread-local storage, so that each thread has its own copy, and it is freed by flint_cleanup.
The number N of cached values and the array of cached values B_0, \ldots, B_{N-1}.
void bernoulli_cache_compute(slong n)Makes sure that B_0, \ldots, B_{n-1} are cached. The cache grows by at least 32 entries (below 128) or 128 entries at a time; new entries are computed with arith_bernoulli_number_vec up to 128 and with bernoulli_fmpq_vec_no_cache beyond.
Makes sure that B_n is cached, by calling bernoulli_cache_compute(n + 1) if needed.
Bounds#
Returns an integer b with |B_n| \le 2^b, using a table for n < 512 and otherwise the inequality |B_n| < 4\, n!/(2\pi)^n < 4 (n+1)^{n+1} e^{-n}/(2\pi)^n evaluated in integer arithmetic with a table of logarithms. Returns WORD_MIN if B_n = 0. The bound exceeds the true exponent by at most one percent. Throws if n > 2^{\mathtt{FLINT\_BITS} - 6} or the computation would overflow; the function is meant for n small enough that B_n could be computed exactly.
Internal functions#
These inline helpers estimate the precision and number of terms for the zeta series used by the reverse iterator.
Returns a crude estimate of the number of bits of the denominator of B_n: by von Staudt–Clausen the denominator is at most n^{\sigma_0(n)}, and the estimate is \tfrac12 \log_2(n)\, n^{1.066/\log\log n}, using \sigma_0(n) < n^{1.066/\log\log n} for large n.
slong bernoulli_global_prec(ulong nmax)Returns the working precision of the iterator for B_{\mathit{nmax}}: arith_bernoulli_number_size(nmax) + bernoulli_denom_size(nmax). No guard bits are added, as the denominator estimate is pessimistic.
slong bernoulli_zeta_terms(ulong s, slong prec)Returns the number N of terms of the zeta series \sum_{j \le N} j^{-s} needed for an error below 2^{-\mathit{prec}}: about 2^{(\mathit{prec} + 1)/(s - 1)}, made odd.
slong bernoulli_power_prec(slong i, ulong s1, slong wp)Returns the precision for the term i^{-s_1} in a sum computed to absolute precision \mathit{wp}: \mathit{wp} - s_1 \log_2 i, at least 10.
The reverse iterator uses a table for n < BERNOULLI_REV_MIN (32).
References#
- [Har2010] D. Harvey, "A multimodular algorithm for computing Bernoulli numbers" (2010), Mathematics of Computation 79.272: 2361-2370
The FLINT team, “bernoulli.h (module overview),” FLINT: Fast Library for Number Theory, documentation, version 3.7.0-dev, 2026. https://flintlib.org/doc/bernoulli/index.html (generated 2026-10-20).
@misc{flint-doc-bernoulli-index,
author = {{The FLINT team}},
title = {{bernoulli.h (module overview)}},
howpublished = {FLINT: Fast Library for Number Theory, documentation, version 3.7.0-dev},
year = {2026},
url = {https://flintlib.org/doc/bernoulli/index.html},
note = {Generated 2026-10-20}
}