FLINT3.7.0-dev
/

Algorithms

Function documentation says what a function computes. This part of the documentation says how. Overview pages compare the methods FLINT uses for one task and say which functions use which; algorithm pages state a single method precisely, prove it correct where the proof fits on a page, and give references to the literature. Every algorithm page lists the functions and pages that refer to it.

Modular arithmetic

Arithmetic modulo word-size and multiword integers: division with precomputed inverses, Shoup, Montgomery, Barrett and Lemire reduction, divisibility tests, simultaneous inversion, multiword modular multiplication and delayed reduction, and how to choose among them.

OverviewChoosing a method for modular arithmeticFLINT has a dozen ways to compute modulo a word-size n, for single values, vectors, SIMD lanes and matrices. This page compares their requirements and costs and says which to use.AlgorithmDivision by an invariant wordDivides a two-word number by a fixed word using one multiplication by a precomputed reciprocal, one low multiplication and two cheap corrections, instead of a hardware division.AlgorithmMultiplication by a fixed operand (Shoup)Computes ab \bmod n for fixed a and n with two low products, one high product and one comparison, using the precomputed quotient \lfloor aB/n \rfloor.AlgorithmMontgomery multiplication (REDC)Represents x modulo an odd n by xR \bmod n, so that a modular product needs two word multiplications and a subtraction instead of a division. Lazy variants keep residues in [0, 2n) and skip the final correction.AlgorithmDivisibility testing by multiplicationTests whether an odd d divides n with one multiplication and one comparison, using d^{-1} \bmod 2^{\mathtt{FLINT\_BITS}}. Also gives the exact quotient for free.AlgorithmBarrett and Lemire reductionReduce a single word modulo n with two products and at most one correction, using \lfloor 2^{\mathtt{FLINT\_BITS}}/n \rfloor (Barrett) or \lceil 2^{\mathtt{FLINT\_BITS}}/n \rceil (Lemire, Kaser and Kurz).AlgorithmSimultaneous inversionInverts k elements modulo n with one modular inversion and 3(k - 1) multiplications, using prefix products (Montgomery's trick).OverviewArithmetic modulo multiword integersFLINT represents \mathbb{Z}/n\mathbb{Z} for n \ge 2^{64} in two ways, with fmpz residues in fmpz_mod and with fixed-length limb arrays in the generic ring mpn_mod. This page compares the representations, the reductions they use, and how they relate to nmod and to the generic rings.AlgorithmModular multiplication with a multiword modulusReduces a product of two residues modulo an \ell-limb n with a precomputed \ell-limb inverse, using one high and one low product and at most a few subtractions; also the two-limb Barrett reduction of fmpz_mod, Shoup's method and the matrix method for a fixed multiplier, and the fused ab + cd.AlgorithmDelayed reduction for multiword moduliHow mpn_mod evaluates sums of products modulo an \ell-limb n exactly in 2\ell - 1 to 2\ell + 1 limbs and reduces once, in dot products, polynomial multiplication and division, Euclid's algorithm, Gaussian elimination and matrix multiplication, and how this plugs into the generic algorithms.

Integer arithmetic

Multiplication, division, GCDs, roots and radix conversion of multiprecision integers, and binary splitting.

OverviewMultiplying integersHow FLINT multiplies integers of each size — its own assembly basecases, Karatsuba, GMP's Toom–Cook algorithms, and FFT multiplication with fft_small or Schönhage–Strassen — with the algorithms, their costs and the dispatch in flint_mpn_mul.AlgorithmBasecase multiplication kernelsSchoolbook multiplication of short limb arrays by a table of straight-line routines, one for each pair of small sizes (assembly with two independent carry chains on x86-64 and ARM64, C elsewhere), and Karatsuba multiplication built on them.AlgorithmHigh, low and middle products of integersComputing only the low half, an approximation of the high half, or an arbitrary window of the limbs of an integer product, with hard-coded kernels, Mulders' recursive short products and proved error bounds.AlgorithmInteger division and Newton iterationEuclidean division of multiprecision integers by schoolbook division with a 3-by-2 inverse, Burnikel–Ziegler divide and conquer, and Newton iteration with the Karp–Markstein trick; exact division by 2-adic (Hensel) methods; division by a fixed divisor with a precomputed inverse; and square roots by the same techniques.AlgorithmInteger GCDGreatest common divisors of multiprecision integers by the binary algorithm, Lehmer's algorithm and the subquadratic half-GCD, in the matrix formulation shared with continued fractions; GMP's implementations behind fmpz_gcd and fmpz_xgcd, and FLINT's own half-GCD for continued fractions and rational reconstruction.AlgorithmRadix conversionConversion between binary and decimal (or any radix b^e) in quasi-linear time by divide and conquer with multiplications only, carrying out the arithmetic in the target radix.AlgorithmMultimodular reduction and Chinese remainderingConverting an integer to and from its residues modulo many word-size primes in quasi-linear time with a subproduct tree, using fractional cofactors so that Chinese remaindering needs no large inverses and no intermediate reductions.AlgorithmBinary splitting and product treesProducts of many integers by balanced product trees, sums of hypergeometric series by binary splitting, and simultaneous reduction and Chinese remaindering by product and remainder trees; how FLINT uses them for factorials, primorials, vector products, Chinese remaindering, radix conversion, constants and series.

Primes and number theory

Primality testing and proving, integer factoring, square roots and discrete logarithms modulo primes, Dirichlet characters, rational reconstruction and continued fractions, number fields, Gaussian integers, binary quadratic forms, and Bernoulli, partition, Stirling and Bell numbers.

OverviewPrimality testing for word-size integersHow n_is_prime proves primality of any 64-bit integer with at most two strong probable prime tests, using table lookups, vectorised trial division and a hashed table of witnesses.OverviewFactoring word-size integersHow n_factor combines trial division, primality testing, perfect power detection, Pollard–Brent rho, Hart's one-line method and SQUFOF to factor any 64-bit integer.AlgorithmPrimality proving: APRCL and ECPPHow FLINT proves the primality of large integers, by n \pm 1 tests when the factorisations allow it, by the APRCL test with Jacobi sums, and by elliptic curve primality proving with verifiable certificates.OverviewFactoring integersThe factoring methods available in FLINT (trial division, Pollard rho, p \pm 1, ECM, SQUFOF and the quadratic sieve), their running times, and how fmpz_factor combines them.AlgorithmSquare roots modulo a primeHow nmod_sqrt finds a square root modulo a word-size prime n. It uses one exponentiation when n \equiv 3 \pmod 4 or n \equiv 5 \pmod 8, and otherwise Tonelli–Shanks in two forms or Müller's variant of Cipolla–Lehmer, depending on the power of 2 dividing n - 1. Every candidate root is checked by squaring.AlgorithmThe elliptic curve methodLenstra's elliptic curve method for finding a prime factor p of N, with Montgomery curves, Suyama's parametrisation, the Montgomery ladder in stage 1 and the standard continuation in stage 2.AlgorithmThe quadratic sieveFactoring N by finding a congruence of squares from many smooth values of quadratic polynomials, with self-initialising polynomials, a block and bucket sieve, the large prime variation and block Lanczos linear algebra.AlgorithmDiscrete logarithms modulo a prime (Pohlig–Hellman)Computes x with g^x = y modulo a prime p by solving the problem in each subgroup of prime-power order of (\mathbb{Z}/p\mathbb{Z})^\times, one base-q digit at a time, and combining the results by the Chinese remainder theorem. Fast when p - 1 has only small prime factors.AlgorithmDiscrete logarithm tablesPrecomputations for many discrete logarithms in a fixed cyclic group modulo a word-size prime power, combining baby-step giant-step tables sized for the number of queries, the Pohlig–Hellman decomposition and p-adic logarithms, and a sieve for the logarithms of all small integers.AlgorithmStrong probable prime testThe Miller–Rabin test for a single base. Every odd prime passes it, at most a quarter of the bases let a composite pass, and with tabulated pseudoprimes it becomes a proof of primality for word-size integers.AlgorithmComputing the partition functionComputes p(n) exactly by the Hardy–Ramanujan–Rademacher series, truncated with Rademacher's error bound, with the exponential sums factored into products of cosines and each term evaluated in ball arithmetic at the precision it needs; small values come from Euler's pentagonal recurrence.AlgorithmGaussian integersDivision with remainder, the canonical associate, and the four GCD algorithms of fmpzi (Euclidean in exact and floating-point arithmetic, the (1+i)-ary binary algorithm, and a lattice algorithm), together with primality testing and the removal of factors 1+i.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.AlgorithmArithmetic in number fieldsHow nf_elem represents elements of Q[x]/(f) with special cases for degrees 1 and 2 and for Q(i), multiplies with reduction by precomputed powers of the generator, inverts by the extended Euclidean algorithm, and computes norms as resultants and traces from power sums of the roots.AlgorithmRational reconstructionRecovers a fraction n/d with |n| \le N and 0 < d \le D from its residue modulo m > 2ND by stopping the extended Euclidean algorithm at the first remainder below N. Used to return rational results from multimodular and p-adic algorithms.AlgorithmContinued fractions of rational numbersComputes the continued fraction of a rational number with a half-GCD algorithm that finds the terms from leading bits, evaluates continued fractions by balanced matrix products, and uses the same machinery for the simplest fraction in an interval, Farey neighbours and Dedekind sums.AlgorithmDirichlet charactersRepresents the Dirichlet characters modulo q through the decomposition of the unit group into cyclic components with Conrey's generators, so that characters, values, conductors, parities and orders are computed from discrete logarithms in each component or directly from the Conrey number.AlgorithmStirling and Bell numbersIsolated Stirling numbers by binary splitting of rising factorials, power sums or power series, chosen by k/n; rows by convolution; Bell numbers by Dobiński's formula, a multimodular sum, the Bell triangle or generating functions modulo primes.AlgorithmBinary quadratic forms and class groupsReduction of positive definite binary quadratic forms, enumeration of the reduced forms of a discriminant (the class number), composition by Shanks' NUCOMP, and the computation of orders and of the exponent of the form class group by a two-stage method in the style of Pollard and Sutherland.

Generic algorithms

How the generic rings framework is designed: contexts and method tables, the three levels of operations and specialisation, inexact elements and truth values, domain and unable errors, and randomised testing of ring axioms.

OverviewThe design of generic ringsAn overview of the generic rings framework: run-time ring objects with method tables, headerless elements, explicit status codes and three-valued predicates, and how these choices make generic code both rigorous and fast.AlgorithmContexts, method tables and specialisationHow rings and elements are represented, how a generic call is dispatched through the method table, the three levels of operations with their defaults, specialised parametric types, and why contexts are cheap to create.AlgorithmInexact elements and truth valuesThe enclosure semantics that make ball arithmetic, power series with error terms and exact real numbers rigorous in generic code: inclusion-preserving operations, three-valued predicates, and weak normalisation of polynomials.AlgorithmDomain errors and unable errorsThe meaning of the status codes GR_DOMAIN and GR_UNABLE, how statuses compose, why there are no junk values, how algorithms use the distinction for fallback strategies, and dynamic evaluation in rings that pretend to be fields.AlgorithmRandomised testing of ring axiomsHow ring implementations are tested by checking algebraic identities on random elements from skewed distributions, why such tests stay valid for inexact rings, and how generic algorithms are tested over many rings at once.AlgorithmGeneric poweringBinary exponentiation and the sliding window method as implemented for any ring, with operation counts, and when each is preferable.

Vectors and matrices

Dot products and matrix multiplication modulo word-size integers, over \mathbb{Z} and in floating point; Gaussian elimination, fraction-free and multimodular algorithms, characteristic polynomials, lattice reduction, Hermite and Smith normal forms and matrices of polynomials, over generic rings and their specialisations.

OverviewLinear algebra over generic rings and its specialisationsHow the generic matrix module gr_mat organises multiplication, LU and fraction-free LU decomposition, solving, determinants, ranks, echelon forms, characteristic polynomials and eigenvalues over any ring, how it handles rings that are not fields and inexact rings, and how nmod_mat, fmpz_mat, mpn_mod, arb_mat, acb_mat and fq_mat specialise or reuse these algorithms.OverviewMultiplying matrices over Z/nZThe methods FLINT uses to multiply matrices modulo a word-size integer, from classical multiplication with delayed reduction to Strassen–Winograd, SIMD kernels, floating-point gemm with Chinese remaindering and bit-packed kernels for tiny moduli, and how nmod_mat_mul chooses among them.AlgorithmDot products modulo a word-size integerHow FLINT computes \sum a_i b_i \bmod n with a single final reduction, accumulating the exact sum in one, two or three words, or in SIMD lanes, depending on n, the length and the processor.AlgorithmStrassen–Winograd multiplicationThe 7-multiplication scheme FLINT uses for large matrix products over \mathbb{Z}/n\mathbb{Z}, with Bodrato's operation sequence, a schedule with two temporaries, and odd dimensions handled by virtual zero padding.AlgorithmMatrix multiplication kernels modulo nThe blocked SIMD kernels nmod_mat_mul_u32, _u52, _k52 and _fp50, which multiply matrices modulo moduli up to 2^{52} in one pass with delayed reduction, and the floating-point route of nmod_mat_mul_blas with Chinese remaindering.AlgorithmMatrix multiplication modulo small integersHow nmod_mat_mul_u8 multiplies matrices modulo N \le 255 with entries in bytes, replacing multiplications by in-register table lookups for N \le 16, and with bit-packed Strassen–Winograd and the method of four Russians for N = 2 and N = 3.AlgorithmFloating-point matrix multiplication (gemm)How FLINT's own sgemm and dgemm kernels block, pack and thread a matrix product around a register-tile micro-kernel, and when floating-point matrix products of integers are exact.OverviewMultiplying matrices over the integersThe methods FLINT uses to multiply integer matrices — word-size dot products, Waksman and Strassen–Winograd on multiprecision entries, multimodular products through nmod_mat, floating-point gemm modulo small primes, and transforming each entry once with the small-prime FFT — their costs, and when fmpz_mat_mul uses which.AlgorithmGaussian elimination and LU decompositionLU decomposition with row pivoting over \mathbb{Z}/p\mathbb{Z} in its classical, delayed-reduction and recursive forms, and how FLINT derives triangular solving, linear system solving, inverses, determinants, ranks, echelon forms and nullspaces from it; plus division-free determinants and the Howell form over \mathbb{Z}/n\mathbb{Z}.AlgorithmFraction-free eliminationBareiss' fraction-free Gaussian elimination over integral domains, Sylvester's identity, the fraction-free LU decomposition and how FLINT uses it for determinants, ranks, linear systems and echelon forms over \mathbb{Z}, \mathbb{Z}[x] and (\mathbb{Z}/p\mathbb{Z})[x].AlgorithmCharacteristic and minimal polynomialsBerkowitz's division-free algorithm and Danilevsky's reduction to companion form for the characteristic polynomial over \mathbb{Z}/n\mathbb{Z}, and the Krylov method with least common multiples that FLINT uses for the minimal polynomial.AlgorithmMultimodular and p-adic linear algebra over the integersHow FLINT computes determinants, solutions of linear systems, ranks, echelon forms, characteristic and minimal polynomials and products of integer matrices by working modulo word-size primes or prime powers, with Hadamard bounds, rational reconstruction, output-sensitive termination and exact certificates for every answer.AlgorithmLattice reduction (LLL)LLL-reduced bases and the algorithms FLINT implements for them — the original rational LLL, Storjohann's integral variant, and the floating-point L² algorithm of Nguyen and Stehlé in fmpz_lll with Gram–Schmidt data in doubles, heuristic exact dot products or multiprecision floats, certified reducedness tests, removals for van Hoeij's factoring, and the truncation strategy (ULLL) used by fmpz_lll.AlgorithmHermite and Smith normal formsThe Hermite and Smith normal forms of integer matrices and the algorithms FLINT implements for them — classical and gcd-based elimination, Kannan–Bachem, Domich–Kannan–Trotter and other computations modulo a determinant, the Pernet–Stein algorithm, Iliopoulos' algorithm, iterated Hermite forms, and Lübeck's p-adic method for elementary divisors — with the choices made by fmpz_mat_hnf and fmpz_mat_snf.AlgorithmMatrices of polynomialsHow FLINT multiplies matrices over \mathbb{Z}[x] and (\mathbb{Z}/n\mathbb{Z})[x] — classical products of polynomials, Kronecker substitution into one integer matrix product, evaluation and interpolation, and transforming each entry once with the small-prime FFT — and how it computes determinants, ranks, inverses and solutions by fraction-free elimination or by evaluation and interpolation.

Univariate polynomials

Dense univariate polynomials over generic rings, \mathbb{Z}/n\mathbb{Z}, \mathbb{Z} and finite fields: multiplication, Newton iteration, GCDs and resultants, evaluation, interpolation and composition, Hensel lifting, squarefree factorisation, factoring over finite fields and over \mathbb{Z}, and roots.

OverviewPolynomial arithmetic over generic rings and its specialisationsHow the generic polynomial module gr_poly organises multiplication, division, Newton iteration, GCDs, evaluation, composition, power series and factorisation over any ring through overridable ring methods, and how nmod_poly, fmpz_poly, fmpz_mod_poly, fq_poly and arb_poly specialise or reuse these algorithms.OverviewMultiplying polynomials over Z/nZHow nmod_poly multiplies: classical products built on modular dot products, Kronecker substitution (KS, KS2, KS4), and the small-prime FFT of fft_small, with the length and modulus-size cutoffs that choose between them for full, low, middle and high products.OverviewMultiplying polynomials over ZHow fmpz_poly multiplies, depending on the length and the coefficient size - word-size accumulation, classical and packed-limb products, Karatsuba with Toom-Cook interpolation in the polynomial variable, Kronecker substitution, the multimodular small-prime FFT, Schönhage-Strassen convolutions and classical multiplication in transform space - with the exactness lemmas and the dispatcher that chooses among them.AlgorithmKronecker substitutionMultiplies polynomials over \mathbb{Z}/n\mathbb{Z} by packing their coefficients into large integers, multiplying the integers and unpacking; Harvey's KS2 and KS4 variants evaluate at ±2^b and ±2^(-b) to split the integer product into two or four smaller ones.AlgorithmNewton iteration for power series and divisionInverse, quotient, square root, reciprocal square root, exponential and logarithm of power series by Newton iteration with middle products and the Karp–Markstein trick, and polynomial division with remainder via reversed polynomials, as implemented for nmod_poly.AlgorithmMiddle products and the transposition principleThe middle product, a window of coefficients of a product, is the transpose of polynomial multiplication and costs no more than an ordinary product of the corresponding size; FLINT computes it with dot products, a cyclic convolution whose wrap-around misses the window, or transposed FFTs, and uses it in Newton iteration and geometric evaluation.OverviewPolynomial GCDs and resultantsThe Euclidean algorithm, the half-gcd algorithm and the subresultant pseudo-remainder sequence as implemented generically in gr_poly, how the generic code chooses between them from the properties of the coefficient ring, how resultants are derived from the same sequences, and the specialised algorithms of nmod_poly, fmpz_mod_poly, fq_poly and fmpz_poly.AlgorithmGCDs of polynomials over ZThe three GCD algorithms of fmpz_poly - the subresultant pseudo-remainder sequence, the heuristic GCD through one integer GCD, and the multimodular algorithm with adaptive termination and certification by multimodular division - with the proofs that make their results correct, and the multimodular extended GCD.AlgorithmResultants and discriminants over ZResultants of integer polynomials by the subresultant pseudo-remainder sequence and by the multimodular algorithm with Hadamard's bound, discriminants as resultants with the derivative, and resultants with a known divisor, as implemented in fmpz_poly.AlgorithmHalf-GCDThe half-gcd algorithm computes the quotients of the first half of the Euclidean remainder sequence from the top halves of the inputs, recursively, giving GCDs, extended GCDs and resultants in O(M(n) \log n) operations; with a correctness proof and FLINT's cutoffs for nmod_poly.AlgorithmMultipoint evaluation and interpolationEvaluation of a polynomial at many points and interpolation through many points in O(M(n) \log n) operations with the subproduct tree and the remainder tree, and in O(M(n)) for points in geometric progression with Bluestein's chirp transform and the method of Bostan and Schost.AlgorithmComposition and modular compositionComputing f(g), f(g) \bmod x^n and f(g) \bmod h by Horner's rule, divide and conquer, Taylor shifts, the baby-step giant-step method of Brent and Kung, and the algorithm of Kinoshita and Li, with the choices made by nmod_poly.AlgorithmHensel liftingLifting a factorisation modulo a prime p into coprime monic factors to a factorisation modulo p^N, for two factors with the Bézout cofactors and for many factors along a factor tree, with quadratic convergence; the correctness proof and the fused tree nodes used by fmpz_poly.AlgorithmFactoring polynomials over the integersThe Zassenhaus algorithm (factor modulo a prime, Hensel lift, recombine subsets of the modular factors) and van Hoeij's knapsack method with CLD data and gradual feeding to LLL, as implemented in fmpz_poly_factor, with the reductions that come first (content, powers of x, deflation with Capelli certificates, squarefree factorisation), the choice of the prime and the bounds that make the result provable.AlgorithmSquarefree factorisationYun's algorithm in characteristic zero and the gcd iteration with p-th roots over finite fields, which split a polynomial into squarefree, pairwise coprime factors with distinct multiplicities using only gcds and exact divisions.AlgorithmFactoring polynomials over finite fieldsSquarefree, distinct-degree and equal-degree factorisation over \mathbb{F}_q; the Cantor–Zassenhaus, Berlekamp and Kaltofen–Shoup algorithms, the baby-step giant-step distinct-degree factorisation with modular composition, the trace and minimal-polynomial splittings, and how the generic implementation in gr_poly chooses between them.AlgorithmIrreducibility testing over finite fieldsRabin's test, Ben-Or's test and the distinct-degree test for irreducibility over \mathbb{F}_q, their costs on irreducible and on random inputs, how FLINT combines them, and the construction of irreducible polynomials.AlgorithmRoots of polynomials over finite fieldsRoots in \mathbb{F}_q by squarefree factorisation, the split \gcd(f, x^{(q-1)/2} \pm 1) and random shifts, with powers of x + c computed through Taylor shifts; roots modulo prime powers by Hensel lifting of all roots, and modulo composite n by the Chinese remainder theorem.AlgorithmCounting and isolating real roots over ZRoot bounds, Descartes' rule of signs and the Vincent-Collins-Akritas bisection algorithm, Sturm sequences computed as subresultant or primitive remainder sequences, and sign-change isolation for real-rooted polynomials, as used by fmpz_poly to count and isolate the real roots of squarefree integer polynomials exactly.

Multivariate polynomials

Sparse multivariate polynomials: the packed representation, heap and dense multiplication and division, powering and square roots, modular GCD algorithms with dense and sparse interpolation, Hensel lifting, and factorisation over the integers and over finite fields.

OverviewMultivariate polynomials: representations and algorithmsAn overview of FLINT's sparse multivariate polynomial modules: which coefficient rings they cover, how polynomials are represented, and which algorithms are used for multiplication, division, powering, square roots, GCDs, resultants and factorisation, with links to the detailed pages.AlgorithmMonomial orderings and exponent packingHow FLINT stores sparse multivariate polynomials: terms sorted by a monomial ordering, exponent vectors packed into machine words so that comparing, multiplying and dividing monomials are word operations, with guard bits for overflow detection; and the deflations, repackings and changes of variables that the algorithms apply to this representation.AlgorithmSparse and dense multiplicationMultiplication of sparse multivariate polynomials with Johnson's heap algorithm and the improvements of Monagan and Pearce (packed monomials, chaining, small heaps, three-word accumulation, parallel splitting of the output), dense array accumulation, Kronecker substitution for dense inputs, and powering by the FPS recurrence.AlgorithmDivision, exact division and square rootsMultivariate division with remainder with respect to a monomial ordering, and its heap implementation by Monagan and Pearce, in which a heap of size at most the number of terms of the divisor supplies the terms of the current remainder; early termination of exact division, quasi-division over the integers, division by several polynomials, parallel exact division, and the analogous heap algorithm for square roots with its non-square certificates.AlgorithmGreatest common divisorsHow FLINT computes GCDs of multivariate polynomials over the integers: reductions that make most inputs easy (monomial content, deflation, missing variables, univariate cases), degree bounds from random images, and the four modular algorithms it chooses between (Brown's dense algorithm, Zippel's sparse algorithm, the Berlekamp–Massey based algorithm of Hu and Monagan, and EEZ-GCD Hensel lifting), with their correctness proofs and the cost model behind the choice.AlgorithmGCDs over small finite fieldsHow the multivariate GCD algorithms over \mathbb{F}_p and \mathbb{F}_q cope with fields too small for random evaluation points, by computing images over extension fields \mathbb{F}_{p^k} or modulo irreducible polynomials of increasing degree; Brown's algorithm over \mathbb{F}_p with pairs of points and early termination; and how nmod_mpoly_gcd chooses among Brown's, Zippel's and Hensel-lifting algorithms.AlgorithmSparse interpolationRecovering a sparse polynomial from its values at a few points: Zippel's method for a known support (transposed Vandermonde systems solved in quadratic time, with extra unknowns for scaled images), and the method of Ben-Or and Tiwari for an unknown support (Berlekamp–Massey, roots of the recurrence polynomial, discrete logarithms modulo a smooth prime after Kronecker substitution) with the Berlekamp–Massey algorithm computed incrementally by half-gcd steps, as used by FLINT's sparse GCD and factoring algorithms.AlgorithmMultivariate Hensel liftingLifting a factorisation of A(x_0, \alpha_1, \ldots, \alpha_m) to a factorisation of A one variable at a time: uniqueness given the leading coefficients, the multivariate partial fraction (diophantine) equations solved by Taylor expansion, the quartic and quintic error updates of FLINT's fmpz_mpoly_hlift, bivariate lifting modulo prime powers, sparse lifting, and the leading coefficient methods of Wang and Kaltofen.AlgorithmFactorisationHow FLINT factors multivariate polynomials over the integers: removal of contents and powers of variables, squarefree decomposition by Yun's algorithm, compression of the exponents by unimodular changes of variables, irreducibility certificates from Newton polytopes, direct formulas for low degree, dense bivariate factorisation with Hensel lifting and recombination, degree pruning with univariate images, and the three lifting strategies of Wang, Zippel and Zassenhaus.AlgorithmFactoring over finite fieldsThe multivariate factorisation pipeline of nmod_mpoly_factor and fq_nmod_mpoly_factor, with what is specific to characteristic p: separable factorisation with p-th roots, deflation by powers of p, factoring over extension fields with recombination of Frobenius orbits when the field is too small for good evaluation points, quadratic equations in characteristic 2, bivariate factorisation with recombination by linear algebra, and Wang's leading coefficients along random curves.AlgorithmBivariate resultants by evaluation and interpolationThe resultant \mathrm{Res}_y(A, B) \in \mathbb{F}_p[x] of two bivariate polynomials, computed from univariate resultants at many values of x (a geometric progression or roots of unity) and interpolation, with the degree bound, the correction for specialisations that lower the degree, and FLINT's cutoffs.

Fast Fourier transforms

The FFTs behind fast multiplication of integers and polynomials: Schönhage–Strassen, small-prime FFTs in floating point, and generic DFTs over any ring.

OverviewFFTs in FLINTFLINT has three FFT implementations — Schönhage–Strassen over \mathbb{Z}/(2^N+1)\mathbb{Z}, small-prime transforms in floating point, and generic DFTs over any ring. This page compares them, says which one each multiplication routine uses, and maps the techniques they share.AlgorithmSchönhage–Strassen multiplicationMultiplies large integers with FFTs over \mathbb{Z}/(2^N + 1)\mathbb{Z}, where the roots of unity are powers of 2 and every twiddle factor is a shift; FLINT adds the \sqrt 2 trick, truncation, negacyclic pointwise products and tuned parameters.AlgorithmTruncated FFTsVan der Hoeven's truncated Fourier transform evaluates a polynomial of length t at the first t points of an FFT of length L = 2^\ell \ge t in about (t/2) \log_2 L butterflies, and its inverse recovers the polynomial from these values, which removes the jumps of radix-2 convolution at powers of two.AlgorithmThe matrix Fourier algorithmBailey's four-step algorithm computes an FFT of length n_1 n_2 as n_1 transforms of length n_2, twiddle factors and n_2 transforms of length n_1, so that each partial transform fits in cache; FLINT fuses it with the pointwise products of a convolution, truncates at whole rows and runs the partial transforms in threads.AlgorithmFFTs modulo word-size primes in floating pointThe number-theoretic transforms of fft_small modulo primes p < 2^{50}, on integers stored in doubles: the order of the roots, radix-4 butterflies with lazy reduction, blocking and truncation, and proofs that every value stays in range.AlgorithmMultiplication with several primes and CRTHow fft_small multiplies integers and polynomials whose product coefficients exceed one prime: chunk sizes and the number of primes, the two-prime and multi-prime pipelines, Chinese remaindering with the normalisation folded in, windows computed with wrap-around, sums of products in transform space, precomputed and shared operands, and transposed products for middle products.AlgorithmDiscrete Fourier transforms over a generic ringThe DFT algorithms of gr_dft (naive, radix-2 and split-radix Cooley–Tukey, four-step, mixed-radix with direct prime kernels, Good–Thomas and Bluestein), the conditions each one needs on the root of unity, and their costs.AlgorithmRigorous numerical DFTs in fixed pointHow gr_dft computes complex DFTs in multi-limb fixed-point arithmetic with an a priori bound on every intermediate value and on the output error, and how this yields ball enclosures for transforms of acb vectors.

Finite fields

Arithmetic in finite fields of non-prime order: the polynomial and Zech logarithm representations, reduction modulo sparse defining polynomials, the Frobenius map, norm and trace, embeddings between fields, and the choice of defining polynomials.

OverviewArithmetic in finite fieldsFLINT represents \mathbb{F}_{p^d} by Zech logarithms, by polynomials over a word-size or a multiprecision prime field, or, for d = 1, by residues; this page compares the representations, says how fq_default chooses among them, and lists the algorithms for reduction, inversion, powering, Frobenius, norm, trace, embeddings and defining polynomials.AlgorithmReduction modulo the defining polynomialElements of \mathbb{F}_p[x]/(f) are reduced after each product either by eliminating the high coefficients with the few nonzero terms of a sparse f, or by Newton division with a precomputed inverse of the reversal of f; sums of products are reduced only once.AlgorithmZech logarithmsIn a small field \mathbb{F}_q with a primitive element g, storing each element as its discrete logarithm turns multiplication into an addition modulo q - 1, and a table of the Zech logarithms Z(k) with 1 + g^k = g^{Z(k)} turns addition into one lookup; the table is built in O(q) operations, which also tests whether g is primitive.AlgorithmFrobenius map, norm, trace and rootsHow FLINT applies the Frobenius automorphism a \mapsto a^{p^e}, computes the trace from power sums given by Newton's identities and the norm as a resultant, and computes p-th roots, square roots (Tonelli–Shanks) and multiplicative orders in \mathbb{F}_{p^d}.AlgorithmEmbeddings of finite fieldsTo embed \mathbb{F}_{p^m} into \mathbb{F}_{p^n} given by unrelated defining polynomials, FLINT finds a root of the smaller field's polynomial in the larger field and expresses the embedding, and a left inverse built from the trace and the dual basis, as matrices over \mathbb{F}_p.AlgorithmDefining polynomialsWhich irreducible polynomial defines \mathbb{F}_{p^d} in FLINT: Conway polynomials from Lübeck's tables when available, otherwise polynomials with as few terms as possible (fq_nmod), random sparse irreducible polynomials (fq) or random primitive polynomials (fq_zech); with the facts about their existence that these choices rely on.

p-adic numbers

Arithmetic with p-adic numbers and unramified extensions: precision models, Newton and Hensel lifting, exponentials and logarithms, Frobenius and Teichmüller lifts, norms and traces.

Overviewp-adic arithmetic in FLINTThe modules for p-adic numbers, their precision models and representations, and which algorithm each function uses for inversion, square roots, the exponential and logarithm, Teichmüller lifts, the Frobenius map, norms and traces.AlgorithmPrecision and representation of p-adic numbersHow FLINT represents a p-adic number as p^v u modulo p^N, the difference between absolute and relative precision, how much precision each arithmetic operation preserves, and the two precision models used by padic/qadic and by padic_radix.AlgorithmHensel lifting of inverses, square roots and Teichmüller liftsNewton iteration in Z_p and Z_q with the precision sequence N, ⌈N/2⌉, ..., 1, applied to inverses, inverse square roots (with the final step producing the square root) and Teichmüller lifts, including the modifications needed for p = 2.AlgorithmThe p-adic exponential and logarithmConvergence and truncation bounds for the p-adic exponential and logarithm series, their evaluation by rectangular splitting with deferred denominators, the balanced (bit-burst) method with binary splitting, and the Satoh–Skjernaa–Taguchi reduction of the logarithm, as implemented in padic, qadic and padic_radix.AlgorithmArithmetic in unramified extensionsHow qadic represents the unramified extension Q_q = Q_p[X]/(f) and computes products with reduction modulo a sparse defining polynomial, inverses, powers, traces by Newton's identities, norms by resultants or as exp(Tr(log x)), and square roots including the Artin–Schreier step for p = 2.AlgorithmFrobenius automorphism and Teichmüller lifts in unramified extensionsHow qadic computes the Frobenius automorphism of Q_q by Newton lifting the image of the generator and evaluating by rectangular splitting, and how Teichmüller lifts are computed, with the relations between the two.

Ball arithmetic

Rigorous real and complex arithmetic with midpoint-radius balls: representation and error propagation, dot products, polynomial and matrix multiplication, power series, roots of polynomials and real functions, linear systems, eigenvalues and integration.

OverviewBall arithmeticHow FLINT represents real and complex numbers as midpoint-radius balls, how each operation computes a midpoint and a bound for the propagated and rounding errors, why the radius has its own fast 30-bit type, and which pages describe the algorithms.AlgorithmRepresentation of radii and midpointsThe formats of mag_t (30-bit mantissa, upper-bound rounding) and arf_t (variable-length mantissa, correct rounding), with proofs that the basic mag operations give upper bounds, and how a rounding error becomes part of the radius.AlgorithmError propagationBounds for the propagated error of arithmetic operations, roots and elementary functions on balls, with proofs, together with the tail bounds for series that are used throughout the special function code.AlgorithmDot products of ballsHow arb_dot sums many products of balls in one fixed-point accumulator with a single final rounding, bounds all radius contributions with 30-bit word arithmetic, and lowers the working precision to what the radii justify.AlgorithmMultiplying polynomials with ball coefficientsHow arb_poly and acb_poly multiply polynomials with ball coefficients by scaling, splitting into blocks of similar magnitude and exact integer polynomial products, with separate products of nonnegative polynomials for the radii; with proofs that the results are enclosures.AlgorithmPower series in ball arithmeticHow truncated power series with ball coefficients give rigorous Taylor expansions, why their correctness is automatic, and how the choice between recurrences and Newton iteration affects the width of the enclosures; with the algorithms used by arb_poly, acb_poly and the generic gr_poly code for inversion, exponentials, roots, logarithms, trigonometric functions and composition.AlgorithmCertified roots of polynomialsHow acb_poly_find_roots computes all complex roots of a polynomial by the Weierstrass–Durand–Kerner iteration from Newton polygon starting values and certifies them with inclusion disks; how real roots are recognised; root bounds; and how arb_fmpz_poly combines this with exact real root isolation.AlgorithmIsolating and refining real rootsBisection with exclusion and monotonicity tests for isolating the simple real roots of an analytic function given by a ball evaluator, Newton iteration with a rigorous quadratic error bound, the interval Newton operator, and precision doubling, as used by arb_calc, arb_poly and arb_fmpz_poly.AlgorithmBall matrix multiplicationHow arb_mat_mul and acb_mat_mul multiply ball matrices by splitting the inner dimension into blocks of entries with similar exponents, multiplying the midpoints exactly as integer matrices, and bounding the radius product in double precision, with the resulting error bounds.AlgorithmVerified linear algebraHow arb_mat and acb_mat solve linear systems, invert matrices and compute determinants with rigorous enclosures, by interval Gaussian elimination, by preconditioning with an approximate inverse, and by certifying an approximate solution with a residual bound.AlgorithmEigenvalue enclosuresHow acb_mat computes rigorous enclosures of the eigenvalues and eigenvectors of complex matrices by certifying floating-point approximations from the QR algorithm, with a global perturbation bound, a Gershgorin test after an approximate diagonalisation, and Rump's fixed-point method for eigenvalue clusters.AlgorithmRigorous numerical integrationHow acb_calc_integrate computes enclosures of integrals of analytic functions along line segments in the complex plane, combining direct interval enclosures, Gauss–Legendre quadrature with error bounds from bounds of the integrand on Bernstein ellipses, and adaptive bisection; with the stopping rules and the role of the user's evaluation function.

Special functions

Rigorous evaluation of elementary functions, constants, the gamma function, hypergeometric functions and their special cases, zeta and L-functions, and elliptic, modular and theta functions in ball arithmetic.

OverviewEvaluating special functionsHow FLINT evaluates special functions with rigorous error bounds — the common principles (enclosures, error propagation, cancellation-aware working precision, argument reduction, convergent and asymptotic expansions, analytic continuation, caching) and a map of the methods to the functions and modules that use them.AlgorithmElementary functionsHow arb evaluates exp, log, sin, cos and atan with rigorous error bounds - argument reduction modulo log 2 and π/4, precomputed tables, Taylor series by rectangular splitting on fixed-point numbers, reduction by logarithms of primes, Newton iteration and the bit-burst algorithm - and how the algorithm is chosen by precision.AlgorithmMathematical constantsThe formulas and algorithms behind arb_const_pi, arb_const_e, arb_const_log2, arb_const_log10, arb_const_euler, arb_const_catalan, arb_const_apery, arb_const_zeta5, arb_const_khinchin, arb_const_glaisher and arb_const_reciprocal_fibonacci, with their error bounds, and how the values are cached.AlgorithmThe gamma functionRigorous evaluation of Γ, 1/Γ and log Γ in ball arithmetic — the Stirling series with Olver's error bound and argument reduction, rising factorials by rectangular and binary splitting, a tabulated Taylor series of 1/Γ for low and medium precision, fast formulas for rational arguments, and the principal branch of log Γ.AlgorithmHypergeometric functionsThe general machinery for hypergeometric functions in ball arithmetic — convergent series with rigorous tail bounds (also for power series parameters), asymptotic expansions with Olver's error bounds, the choice between them by estimated cancellation, connection formulas and regularization, and the summation algorithms (forward recurrence, rectangular splitting, binary splitting, fast multipoint evaluation, bit-burst continuation).AlgorithmSums of hypergeometric series with rational parametersPartial sums of hypergeometric series whose parameters are small rational numbers and whose argument is a real ball — the forward recurrence with a common denominator, rectangular splitting with exact integer block coefficients and decreasing precision, binary splitting for short arguments, and the variant that splits real and imaginary parts.AlgorithmError functions, incomplete gamma functions and exponential integralsHow the error functions, Fresnel integrals, incomplete gamma and beta functions and the exponential, logarithmic, sine and cosine integrals are reduced to 1F1, U and 2F1, how FLINT chooses between their convergent and asymptotic expansions, the real-variable algorithms with rectangular splitting, and the bit-burst algorithm for the incomplete gamma function at very high precision.AlgorithmBessel, Airy and Coulomb wave functionsHow the Bessel functions J, Y, I, K, the Airy functions and the Coulomb wave functions are reduced to the confluent hypergeometric functions 0F1 and U, how the method and the working precision are chosen, and the error bounds of the asymptotic expansions of the Airy functions, with the connection formulas that make them valid in the whole plane.AlgorithmThe Gauss hypergeometric functionEvaluation of 2F1(a,b;c;z) on the whole complex plane — the linear fractional transformations and how one is chosen, limits at integer parameter differences by power series perturbation, polynomial and special cases, analytic continuation by Taylor series with Cauchy–Kovalevskaya majorant bounds near exp(±πi/3), and numerical integration for real arguments.AlgorithmLegendre functions and orthogonal polynomialsLegendre polynomials of large degree by the method of Johansson and Mezzarobba (hypergeometric series at 0 and 1, an asymptotic expansion with explicit error bounds, and a recurrence, chosen by cost and cancellation), Gauss–Legendre nodes and weights, associated Legendre functions of both kinds through 2F1 and connection formulas, spherical harmonics, and the classical orthogonal polynomials and their extensions to non-integer degree.AlgorithmThe dilogarithmEvaluation of the dilogarithm Li₂(z) on the whole complex plane — the series at 0 after one of four transformations, a series in log z with Bernoulli numbers near exp(±πi/3), Taylor continuation with explicit coefficient bounds, and the bit-burst algorithm at very high precision.AlgorithmZeta functions, polylogarithms and L-functionsEuler–Maclaurin summation for the Hurwitz zeta function and its derivatives with a rigorous remainder bound, power sums, Taylor expansion in the parameter, the Riemann zeta function at integers, polylogarithms, and the methods used for Dirichlet L-functions.AlgorithmZeros of the Riemann zeta functionThe Hardy Z-function, Gram points and Turing's method for proving the count of zeros, the Riemann–Siegel formula with Arias de Reyna's error bounds, and Platt's method that isolates thousands of zeros at large height with one FFT-based multi-evaluation.AlgorithmElliptic integrals, elliptic functions and modular formsThe complex arithmetic-geometric mean with Dupont's error bound, complete and incomplete elliptic integrals through the AGM and Carlson's duplication algorithm with high-order Taylor series, Weierstrass functions via theta functions, and modular forms by reduction to the fundamental domain and q-series with short addition sequences.AlgorithmTheta functionsJacobi and Riemann theta functions in any genus: characteristics, the action of the symplectic group and reduction, summation over ellipsoids with rigorous tail bounds, and the genus 1 summation with addition sequences.AlgorithmQuasi-linear evaluation of theta functionsThe algorithm of Elkies and Kieffer for Riemann theta functions in any genus: duplication formulas as a generalised AGM, square roots chosen from low-precision summations, an auxiliary real shift for robustness, dimension lowering for skewed matrices, and derivatives by finite differences.AlgorithmThe Lambert W functionEvaluation of all branches of the Lambert W function in real and complex ball arithmetic — starting values from asymptotic expansions, the series at the branch point and double-precision approximations, Halley iteration, and an a posteriori error bound from a bound for W' that also certifies the branch.

Floating-point arithmetic

Fixed-precision multiprecision floating-point and decimal arithmetic, double-word and interval arithmetic with machine doubles.

OverviewFloating-point arithmetic in FLINTA map of FLINT's floating-point types - arf and arb, nfloat, mp_real, decimal, dfloat and machine doubles - with what each represents, its precision range, whether it is rigorous, and which to use.AlgorithmFixed-precision floating point with nfloatHow nfloat gets its speed - a shallow fixed-size representation selected at run time, kernels specialised by the number of words, fixed-point accumulation for dot products, and fixed-point or exact integer matrix multiplication with rigorous scaling bounds.AlgorithmBall arithmetic with limb-granular radiiThe representation of mp_real balls, with the radius counted in ulps of the last limb and the mantissa length following the accuracy; how errors are bounded and installed, and what this costs compared with arb.AlgorithmElementary functions in mp_realHow mp_real evaluates exp, log, the trigonometric functions and their inverses from one limb to millions of limbs: argument reduction of balls, table-driven and dense-table kernels, bitwise (BKM-type) reduction with signed digits, rectangular splitting, bit-burst evaluation, Newton–Taylor inversion, the AGM logarithm, and sines of rational multiples of π.AlgorithmMulti-prime argument reductionReduction of exp and sin/cos arguments by integer combinations of logarithms of primes (or arguments of Gaussian primes) found from precomputed integer relations, with the logarithms computed by Machin-type formulas whose terms are summed by Zuniga's fast series.AlgorithmBinary splitting of hypergeometric series in mp_realThe general hypergeometric binary splitting of mp_real: exact leaves, removal of common factors, truncation near the root, rigorous tail bounds, transform reuse between the products of a merge and precomputed transforms of power tables, and parallelism.AlgorithmConstants in mp_realThe series and formulas behind the mp_real constants (π, log 2, e, Euler's constant, log 10, Catalan's constant, ζ(3), ζ(5), Γ(1/3), Γ(1/4)), the dual-number binary splitting for Euler's constant, and the cache of exact floors.AlgorithmDecimal floating-point arithmeticCorrectly rounded decimal arithmetic with limb-aligned exponents in radix 10^19: rounding with sticky bits, sums of far-apart operands, truncated products with interval rounding, division and square roots, decimal radii, conversions to and from binary, and correctly rounded functions by Ziv's strategy with detection of the hard cases.AlgorithmDouble-word and multi-double arithmeticError-free transformations, floating-point expansions of two to four doubles, their renormalisation, multiplication, long division and square roots, and how dfloat turns the dropped terms into rigorous ball radii.AlgorithmElementary functions in double-word to quad-double precisionTable-driven argument reduction, polynomials evaluated in several chains with tapered precision, and static error bounds for exp, log, sin, cos and atan on floating-point expansions, with the propagation of ball radii.AlgorithmInterval arithmetic with doublesFast rigorous enclosures with double endpoints computed in round-to-nearest mode and widened by a relative perturbation, with a proof that the widening is sufficient, and the mean value form used to bound integrands.

Exact real and complex numbers

Algebraic numbers and exact real and complex numbers: representations, arithmetic and decision procedures.

OverviewChoosing a representation for exact numbersFLINT represents exact real and complex numbers as rationals, elements of a fixed number field, individual algebraic numbers (qqbar) or elements of formal fields with transcendental generators (ca). This page compares what each representation can decide, what it costs, and which algorithms the main operations use.AlgorithmAlgebraic numbers by minimal polynomial and isolating ballThe canonical representation of an algebraic number by its minimal polynomial over Z and a complex ball isolating one root, its certification by the interval Newton method, refinement, and exact decisions of equality, realness and signs, as implemented in qqbar.AlgorithmArithmetic with algebraic numbersHow qqbar computes x + y, x - y, xy and x/y exactly: composed polynomials from power sums, factorisation and numerical selection of the right factor, closed forms and certified norms for quadratic operands, half-degree polynomials for conjugates, linear algebra in a common number field or cyclotomic field found by LLL, and guessed minimal polynomials with a proof of correctness.AlgorithmFormal fields and the ideal of relationsHow ca represents numbers as rational functions in extension numbers, modulo an ideal of algebraic relations that is built when a field is created, from minimal polynomials, defining equations of roots and powers, Vieta's formulas, integer relations between logarithms and exponentials found by LLL and proved exactly, functional equations and a Gröbner basis.AlgorithmZero testingWhy deciding x = 0 is hard for exact real and complex numbers, what is known about decidability, and how ca combines numerical evaluation in ball arithmetic (which proves x ≠ 0) with exact methods (formal reduction, evaluation as an algebraic number, rewriting in a normal form, factoring) that prove x = 0.AlgorithmSimplification and closed formsThe simplifications that ca applies automatically (canonical number field representations, reduction modulo relations, rewriting rules for exponentials, logarithms, powers and roots), and the LLL-based methods with exact verification that qqbar uses to recognise algebraic numbers, express them in a given field and find closed forms with radicals and trigonometric functions.AlgorithmLinear algebra over exact numbersGaussian elimination, fraction-free elimination, division-free determinants and characteristic polynomials, eigenvalues, Jordan forms and matrix functions over exact real and complex numbers, where every pivot decision is a zero test that may fail; how ca_mat builds on the generic gr_mat algorithms and specialises to rational and number field matrices.

SIMD and floating-point arithmetic

Arithmetic on vectors of machine words and doubles, including exact modular arithmetic carried out in floating point.

AlgorithmModular arithmetic in floating pointExact arithmetic modulo n < 2^{50} on integers stored in doubles, using a precomputed 1/n and a fused multiply-add to split the product, as in FLINT's small-prime FFT.

Parallel computation

How FLINT distributes work over threads: the thread pool, parallel loops and parallel binary splitting.

AlgorithmThe thread poolA fixed set of sleeping threads, handed out on request and woken with a function to run, with per-thread worker budgets that keep nested parallel code within the limit set by the user.AlgorithmParallel loops and binary splittingHow flint_parallel_do schedules independent calls over threads, why strided scheduling balances monotone costs to within one call, and how parallel binary splitting divides threads over a recursion tree.AlgorithmLock-free shared cachesRead-mostly tables shared by all threads, read with one acquire load and grown by publishing a new version under a lock; retired versions are kept, and geometric growth bounds the memory they use.
Cite this page

The FLINT team, “Algorithms,” FLINT: Fast Library for Number Theory, documentation, version 3.7.0-dev, 2026. https://flintlib.org/doc/algorithms/index.html (generated 2026-10-20).