Generic Convolution.

Asymptotically fast convolution of lists over any commutative ring in which the multiply-by-two map is injective. (More precisely, if x \in R, and x = 2^k*y for some k \geq 0, we require that R(x/2^k) returns y.)

The main function to be exported is convolution().

EXAMPLES:

sage: convolution([1, 2, 3, 4, 5], [6, 7])
[6, 19, 32, 45, 58, 35]

The convolution function is reasonably fast, even though it is written in pure Python. For example, the following takes less than a second:

sage: v=convolution(range(1000), range(1000))

ALGORITHM: Converts the problem to multiplication in the ring S[x]/(x^M - 1), where S = R[y]/(y^K + 1) (where R is the original base ring). Performs FFT with respect to the roots of unity 1, y, y^2, \ldots, y^{2K-1} in S. The FFT/IFFT are accomplished with just additions and subtractions and rotating python lists. (I think this algorithm is essentially due to Schonhage, not completely sure.) The pointwise multiplications are handled recursively, switching to a classical algorithm at some point.

Complexity is O(n log(n) log(log(n))) additions/subtractions in R and O(n log(n)) multiplications in R.

AUTHORS:

  • David Harvey (2007-07): first implementation
  • William Stein: editing the docstrings for inclusion in Sage.
sage.rings.polynomial.convolution._combine(L, m, k)
Assumes L is a list of length 2^m, each entry a list of length 2^k. Combines together into a single list, effectively inverting _split(), but overlaying coefficients, i.e. list #i gets added in starting at position 2^{k-1} i. Note that the second half of the last list is ignored.
sage.rings.polynomial.convolution._convolution_fft(L1, L2)

Returns convolution of non-empty lists L1 and L2, using FFT algorithm. L1 and L2 may have arbitrary lengths \geq 4. Assumes all entries of L1 and L2 belong to the same ring.

EXAMPLES:

sage: from sage.rings.polynomial.convolution import _convolution_naive
sage: from sage.rings.polynomial.convolution import _convolution_fft
sage: _convolution_naive([1, 2, 3], [4, 5, 6])
[4, 13, 28, 27, 18]
sage: _convolution_fft([1, 2, 3], [4, 5, 6])
[4, 13, 28, 27, 18]
sage: for len1 in range(4, 30):
...      for len2 in range(4, 30):
...         L1 = [ZZ.random_element(100) for _ in range(len1)]
...         L2 = [ZZ.random_element(100) for _ in range(len2)]
...         assert _convolution_naive(L1, L2) == _convolution_fft(L1, L2)
sage.rings.polynomial.convolution._convolution_naive(L1, L2)

Returns convolution of non-empty lists L1 and L2, using naive algorithm. L1 and L2 may have arbitrary lengths.

EXAMPLES:

sage: from sage.rings.polynomial.convolution import _convolution_naive
sage: _convolution_naive([2], [3])
[6]
sage: _convolution_naive([2, 5], [3])
[6, 15]
sage: _convolution_naive([2], [3, 6])
[6, 12]
sage: _convolution_naive([1, 2, 3], [4, 5, 6, 7])
[4, 13, 28, 34, 32, 21]
sage: _convolution_naive([4, 5, 6, 7], [1, 2, 3])
[4, 13, 28, 34, 32, 21]
sage.rings.polynomial.convolution._fft(L, K, start, depth, root)

L is a list of length M = 2^m, each entry a list of length K = 2^k.

This function only operates on the [start, start + D) portion of L, where D = 2^\text{depth}. This portion is interpreted as a polynomial in S[x]/(x^D - y^(2*root)), where S = R[y]/(y^K + 1).

This function performs an inplace FFT, i.e. evaluates the polynomial at x = each D-th root of unity in S (namely the powers of y^{2K/D}), with results in bit-reversed order.

sage.rings.polynomial.convolution._forward_butterfly(L1, L2, r)
L1 and L2 are both lists of length K, and 0 \leq r \leq K. They represent polynomials in S = R[y]/(y^K + 1). This function returns (L_1 + y^r L_2, L_1 - y^r L_2), as a list.
sage.rings.polynomial.convolution._ifft(L, K, start, depth, root)
Inverse operation of _fft_trunc() (except that result is a factor of 2^depth too big)
sage.rings.polynomial.convolution._inverse_butterfly(L1, L2, r)
L1 and L2 are both lists of length K, and 0 \leq r \leq K. They represent polynomials in S = R[y]/(y^K + 1). This function returns (L_1 + L_2, y^{-r}*(L_1 - L_2)), as a list.
sage.rings.polynomial.convolution._nega_combine(L, m, k)
Same as _combine(), but doesn’t ignore the second half of the last list; instead it makes that piece wrap around negacyclically.
sage.rings.polynomial.convolution._negaconvolution(L1, L2, n)
Negacyclic convolution of L1 and L2. L1 and L2 must both be length 2^n.
sage.rings.polynomial.convolution._negaconvolution_fft(L1, L2, n)

Returns negacyclic convolution of lists L1 and L2, using FFT algorithm. L1 and L2 must both be length 2^n, where n \geq 3. Assumes all entries of L1 and L2 belong to the same ring.

EXAMPLES:

sage: from sage.rings.polynomial.convolution import _negaconvolution_naive
sage: from sage.rings.polynomial.convolution import _negaconvolution_fft
sage: _negaconvolution_naive(range(8), range(5, 13))
[-224, -234, -224, -192, -136, -54, 56, 196]
sage: _negaconvolution_fft(range(8), range(5, 13), 3)
[-224, -234, -224, -192, -136, -54, 56, 196]
sage: for n in range(3, 10):
...      L1 = [ZZ.random_element(100) for _ in range(1 << n)]
...      L2 = [ZZ.random_element(100) for _ in range(1 << n)]
...      assert _negaconvolution_naive(L1, L2) == _negaconvolution_fft(L1, L2, n)
sage.rings.polynomial.convolution._negaconvolution_naive(L1, L2)

Negacyclic convolution of L1 and L2, using naive algorithm. L1 and L2 must be the same length.

EXAMPLES:

sage: from sage.rings.polynomial.convolution import _negaconvolution_naive
sage: from sage.rings.polynomial.convolution import _convolution_naive
sage: _negaconvolution_naive([2], [3])
[6]
sage: _convolution_naive([1, 2, 3], [3, 4, 5])
[3, 10, 22, 22, 15]
sage: _negaconvolution_naive([1, 2, 3], [3, 4, 5])
[-19, -5, 22]
sage.rings.polynomial.convolution._split(L, m, k)

Assumes L is a list of length 2^{m+k-1}. Splits it into 2^m lists of length 2^{k-1}, returned as a list of lists. Each list is zero padded up to length 2^k.

EXAMPLES:

sage: from sage.rings.polynomial.convolution import _split
sage: _split([1, 2, 3, 4, 5, 6, 7, 8], 2, 2)
[[1, 2, 0, 0], [3, 4, 0, 0], [5, 6, 0, 0], [7, 8, 0, 0]]
sage: _split([1, 2, 3, 4, 5, 6, 7, 8], 1, 3)
[[1, 2, 3, 4, 0, 0, 0, 0], [5, 6, 7, 8, 0, 0, 0, 0]]
sage: _split([1, 2, 3, 4, 5, 6, 7, 8], 3, 1)
[[1, 0], [2, 0], [3, 0], [4, 0], [5, 0], [6, 0], [7, 0], [8, 0]]
sage.rings.polynomial.convolution.convolution(L1, L2)

Returns convolution of non-empty lists L1 and L2. L1 and L2 may have arbitrary lengths.

EXAMPLES:

sage: convolution([1, 2, 3], [4, 5, 6, 7])
[4, 13, 28, 34, 32, 21]
sage: R = Integers(47)
sage: L1 = [R.random_element() for _ in range(1000)]
sage: L2 = [R.random_element() for _ in range(3756)]
sage: L3 = convolution(L1, L2)
sage: L3[2000] == sum([L1[i] * L2[2000-i] for i in range(1000)])
True
sage: len(L3) == 1000 + 3756 - 1
True

Previous topic

Educational version of the d-Groebner Basis Algorithm over PIDs.

Next topic

Fast calculation of cyclotomic polynomials

This Page