NTT & Arbitrary-Modulus Convolution

Intuition, Vocabulary and Foundations

Convolution combines every coefficient pair: $c_k=\sum_{i+j=k}a_ib_j$. Thus $(1,2,3)*(4,5)=(4,13,22,15)$. Naive multiplication costs $O(nm)$; an NTT evaluates both polynomials at exact roots of unity in a finite field, multiplies evaluations, and interpolates without floating-point error.

For power-of-two length $N$, an NTT prime $p$ must satisfy $N\mid p-1$. If the requested modulus is not friendly, compute under several friendly primes. Their product must exceed a proven coefficient bound so CRT identifies each integer coefficient uniquely before reduction modulo the target.

An NTT turns convolution into pointwise multiplication using an exact root of unity. Multi-prime reconstruction is correct only when the selected modulus product exceeds the maximum possible unreduced coefficient.

Concept	Example	Boundary	Question answered
Primitive root	$3$ generates the nonzero residues modulo $998244353$	An arbitrary nonzero value need not generate the group	How are roots constructed?
Length constraint	$2^{20}\mid998244352$	$2^{24}\nmid998244352$	Is the transform supported?
Padding	Lengths $3$ and $2$ need $N=4$	Padding to $3$ wraps cyclically	How is linear convolution preserved?
Bound	$\lvert c_k\rvert\le\min(n,m)AB$ for magnitudes $A,B$	Average size is not a proof	How many primes are needed?
CRT range	Product $P>B$ uniquely identifies $0\le c_k\le B$	$P\le B$ permits ambiguity	Is reconstruction exact?

Step-by-Step Derivations and Algorithms

Choose the smallest power of two $N\ge n+m-1$ and define $\omega=g^{(p-1)/N}$. Bit reversal groups inputs for iterative radix-two butterflies. At stage length $\ell$, combine $u=x$ and $v=y\omega_\ell^j$ into $u+v$ and $u-v$. The inverse uses reciprocal roots and multiplies by $N^{-1}$.

For target modulus $M$, normalize inputs into $[0,M)$. With maxima $A,B$, every unreduced coefficient is at most $\min(n,m)AB$. Select compatible primes until their product exceeds this bound, convolve under each, reconstruct by CRT, and finally reduce modulo $M$.

NTT & Arbitrary-Modulus Convolution algorithm map

flowchart TD
    A[Normalize and choose transform length] --> B{Small input?}
    B -->|Yes| C[Naive convolution]
    B -->|No| D[Select compatible NTT primes]
    D --> E{Product exceeds bound?}
    E -->|No| F[Reject insufficient range]
    E -->|Yes| G[Transform and multiply]
    G --> H[Inverse transforms]
    H --> I[CRT then reduce modulo target]

Each butterfly stage computes the transform on doubled blocks, and the inverse reverses it because $N$ is invertible modulo every selected prime. With $r$ primes, time is $O(rN\log N)$, space is $O(N)$ per active transform, and CRT adds $O(rN)$ work.

Implementation, Complexity and Engineering

The reference implementation validates transform lengths, normalizes negative inputs, derives a conservative bound, and rejects an insufficient prime catalog. It uses a naive path for small inputs; production thresholds, roots, and bit-reversal permutations should be benchmarked and cached by modulus and length.

PRIMES = ((998244353, 3), (1004535809, 3), (469762049, 3))


def ntt(values, inverse, modulus, primitive_root):
    n = len(values)
    if n == 0 or n & (n - 1) or (modulus - 1) % n:
        raise ValueError("unsupported transform length")
    j = 0
    for i in range(1, n):
        bit = n >> 1
        while j & bit:
            j ^= bit
            bit >>= 1
        j ^= bit
        if i < j:
            values[i], values[j] = values[j], values[i]
    length = 2
    while length <= n:
        root = pow(primitive_root, (modulus - 1) // length, modulus)
        if inverse:
            root = pow(root, -1, modulus)
        half = length // 2
        for start in range(0, n, length):
            factor = 1
            for offset in range(half):
                left = values[start + offset]
                right = values[start + offset + half] * factor % modulus
                values[start + offset] = (left + right) % modulus
                values[start + offset + half] = (left - right) % modulus
                factor = factor * root % modulus
        length *= 2
    if inverse:
        scale = pow(n, -1, modulus)
        for i in range(n):
            values[i] = values[i] * scale % modulus


def convolution_prime(left, right, modulus, root):
    if not left or not right:
        return []
    size = len(left) + len(right) - 1
    n = 1
    while n < size:
        n <<= 1
    first = [x % modulus for x in left] + [0] * (n - len(left))
    second = [x % modulus for x in right] + [0] * (n - len(right))
    ntt(first, False, modulus, root)
    ntt(second, False, modulus, root)
    first = [x * y % modulus for x, y in zip(first, second)]
    ntt(first, True, modulus, root)
    return first[:size]


def naive(left, right, modulus):
    if not left or not right:
        return []
    result = [0] * (len(left) + len(right) - 1)
    for i, x in enumerate(left):
        for j, y in enumerate(right):
            result[i + j] = (result[i + j] + x * y) % modulus
    return result


def reconstruct(residues, primes, target):
    value, product = 0, 1
    for residue, prime in zip(residues, primes):
        step = (residue - value) % prime
        step = step * pow(product % prime, -1, prime) % prime
        value += product * step
        product *= prime
    return value % target


def convolution(left, right, target, naive_limit=64):
    if not isinstance(target, int) or isinstance(target, bool):
        raise TypeError("target must be an integer")
    if target <= 0:
        raise ValueError("target must be positive")
    if not left or not right:
        return []
    if target == 1:
        return [0] * (len(left) + len(right) - 1)
    left = [x % target for x in left]
    right = [x % target for x in right]
    if len(left) * len(right) <= naive_limit:
        return naive(left, right, target)
    size = len(left) + len(right) - 1
    n = 1
    while n < size:
        n <<= 1
    bound = min(len(left), len(right)) * max(left) * max(right)
    selected, product = [], 1
    for prime, root in PRIMES:
        if (prime - 1) % n == 0:
            selected.append((prime, root))
            product *= prime
            if product > bound:
                break
    if product <= bound:
        raise OverflowError("insufficient CRT range")
    vectors = [convolution_prime(left, right, p, g) for p, g in selected]
    primes = [p for p, _ in selected]
    return [
        reconstruct([vector[i] for vector in vectors], primes, target)
        for i in range(size)
    ]


assert convolution_prime([1, 2, 3], [4, 5], 998244353, 3) == [4, 13, 22, 15]
assert convolution([-1, 2], [3, -4], 17) == [14, 10, 9]
assert convolution([], [1], 17) == []
assert convolution([8, 9], [7], 1) == [0, 0]

left = [i * i - 20 for i in range(40)]
right = [3 * i - 11 for i in range(35)]
for modulus in (97, 1_000_000_007, 1_000_000_000):
    expected = naive([x % modulus for x in left], [x % modulus for x in right], modulus)
    assert convolution(left, right, modulus, naive_limit=0) == expected

data = list(range(16))
transformed = data[:]
ntt(transformed, False, 998244353, 3)
ntt(transformed, True, 998244353, 3)
assert transformed == data

print("all tests passed")

Choice	Cost	Failure mode	Mitigation
Naive	$O(nm)$	Quadratic work dominates	Benchmark a crossover
Single-prime NTT	$O(N\log N)$	Length unsupported or range too small	Check divisibility and bounds
Multi-prime NTT	$O(rN\log N)$	Product does not exceed bound	Add primes or reject
CRT	$O(r)$ per coefficient	Fixed-width product overflow	Use wider or big integers
Negative inputs	Valid modulo target	Negative remainder conventions differ	Normalize into $[0,M)$
Cached roots	Reduces setup	Wrong prime-length reuse	Key by modulus, length, and direction

Applications, Case Studies and Exercises

For $(1,2,3)*(4,5)$, pad to length $4$, transform, multiply, and invert to obtain $(4,13,22,15)$. For a production modulus such as $10^9+7$, normalize inputs, prove a bound, select enough friendly primes, reconstruct exactly, and only then reduce.

Practice

Apply NTT & Arbitrary-Modulus Convolution

Case study 1. Trace every butterfly for $(1,2,3)$ and $(4,5)$ and compare with the four naive sums. Success means both methods return $(4,13,22,15)$.

Case study 2. Convolve two reproducible length-$10{,}000$ arrays modulo $10^9+7$. Record the bound, selected primes, and product; verify random coefficients by direct sums. Success requires agreement and a prime product above the bound.

Exercises. Prove the inverse formula; add root and permutation caches; test unsupported lengths, negative inputs, and insufficient CRT range; benchmark naive and NTT paths to choose the crossover empirically.

Use naive multiplication below the measured crossover, one NTT prime when its length and range suffice, and multi-prime CRT otherwise. Never skip normalization, length compatibility, or the modulus-product proof; those checks make later formal-power-series operations reliable.

Cookie Consent