.. _acb-dft:

**acb_dft.h** -- Discrete Fourier transform
===================================================================================

*Warning: the interfaces in this module are experimental and may change
without notice.*

All functions support aliasing.

Let *G* be a finite abelian group, and `\chi` a character of *G*.
For any map `f:G\to\mathbb C`, the discrete fourier transform
`\hat f:\hat G\to \mathbb C` is defined by

.. math::

   \hat f(\chi) = \sum_{x\in G}\overline{\chi(x)}f(x)

Note that by the inversion formula

.. math::

   \widehat{\hat f}(\chi) = \#G\times f(\chi^{-1})

it is straightforward to recover `f` from its DFT `\hat f`.

Main DFT functions
-------------------------------------------------------------------------------

If `G=\mathbb Z/n\mathbb Z`, we compute the DFT according to the usual convention

.. math::

   w_x = \sum_{y\bmod n} v_y e^{-\frac{2i \pi}nxy}

.. function:: void acb_dft(acb_ptr w, acb_srcptr v, slong len, slong prec)

   Set *w* to the DFT of *v* of length *len*, using an automatic choice
   of algorithm.

.. function:: void acb_dft_inverse(acb_ptr w, acb_srcptr v, slong len, slong prec)

   Compute the inverse DFT of *v* into *w*.

If several computations are to be done on the same group, the FFT scheme
should be reused.

The precomputed schemes build their internal tables of roots of unity at the
precision passed to their *init* function. Build a scheme at (at least) the
precision at which you intend to evaluate transforms; *acb_dft* and the
precomputed schemes then return results accurate to about *prec* bits.

.. type:: acb_dft_pre_struct

.. type:: acb_dft_pre_t

    Stores a fast DFT scheme on :math:`\mathbb Z/n\mathbb Z`
    as a recursive decomposition into simpler DFT
    with some tables of roots of unity.

    An *acb_dft_pre_t* is defined as an array of *acb_dft_pre_struct*
    of length 1, permitting it to be passed by reference.

.. function:: void acb_dft_precomp_init(acb_dft_pre_t pre, slong len, slong prec)

   Initializes the fast DFT scheme of length *len*, using an automatic choice of
   algorithms depending on the factorization of *len*.

   The length *len* is stored as *pre->n*.

.. function:: void acb_dft_precomp_clear(acb_dft_pre_t pre)

   Clears *pre*.

.. function:: void acb_dft_precomp(acb_ptr w, acb_srcptr v, const acb_dft_pre_t pre, slong prec)

   Computes the DFT of the sequence *v* into *w* by applying the precomputed scheme
   *pre*. Both *v* and *w* must have length *pre->n*.

.. function:: void acb_dft_inverse_precomp(acb_ptr w, acb_srcptr v, const acb_dft_pre_t pre, slong prec)

   Compute the inverse DFT of *v* into *w*.

DFT on products
-------------------------------------------------------------------------------

A finite abelian group is isomorphic to a product of cyclic components

.. math::

   G = \bigoplus_{i=1}^r \mathbb Z/n_i\mathbb Z

Characters are product of component characters and the DFT reads

.. math::

   \hat f(x_1,\dots x_r) = \sum_{y_1\dots y_r} f(y_1,\dots y_r)
   e^{-2i \pi \sum\frac{x_i y_i}{n_i}}

We assume that `f` is given by a vector of length `\prod n_i` corresponding
to a lexicographic ordering of the values `y_1,\dots y_r`, and the computation
returns the same indexing for values of `\hat f`.

.. function:: void acb_dft_prod(acb_ptr w, acb_srcptr v, slong * cyc, slong num, slong prec)

   Computes the DFT on the group product of *num* cyclic components of sizes *cyc*. Assume the entries
   of *v* are indexed according to lexicographic ordering of the cyclic
   components.

.. type:: acb_dft_prod_struct

.. type:: acb_dft_prod_t

    Stores a fast DFT scheme on a product of cyclic groups.

    An *acb_dft_prod_t* is defined as an array of *acb_dft_prod_struct*
    of length 1, permitting it to be passed by reference.

.. function:: void acb_dft_prod_init(acb_dft_prod_t t, slong * cyc, slong num, slong prec)

   Stores in *t* a DFT scheme for the product of *num* cyclic components whose sizes are given in the array *cyc*.

.. function:: void acb_dft_prod_clear(acb_dft_prod_t t)

   Clears *t*.

.. function:: void acb_dft_prod_precomp(acb_ptr w, acb_srcptr v, const acb_dft_prod_t prod, slong prec)

   Sets *w* to the DFT of *v*. Assume the entries are lexicographically
   ordered according to the product of cyclic groups initialized in *prod*.

Convolution
-------------------------------------------------------------------------------

For functions `f` and `g` on `G` we consider the convolution

.. math::

   (f \star g)(x) = \sum_{y\in G} f(x-y)g(y)

.. function:: void acb_dft_convol_naive(acb_ptr w, acb_srcptr f, acb_srcptr g, slong len, slong prec)

.. function:: void acb_dft_convol_rad2(acb_ptr w, acb_srcptr f, acb_srcptr g, slong len, slong prec)

.. function:: void acb_dft_convol(acb_ptr w, acb_srcptr f, acb_srcptr g, slong len, slong prec)

   Sets *w* to the convolution of *f* and *g* of length *len*.

   The *naive* version simply uses the definition.

   The *rad2* version embeds the sequence into a power of 2 length and
   uses the formula

   .. math::

      \widehat{f \star g}(\chi) = \hat f(\chi)\hat g(\chi)

   to compute it using three radix 2 FFT.

   The default version uses radix 2 FFT unless *len* is a product of small
   primes where a non padded FFT is faster.

FFT algorithms
-------------------------------------------------------------------------------

Fast Fourier transform techniques allow to compute efficiently
all values `\hat f(\chi)` by reusing common computations.

Specifically, if `H\triangleleft G` is a subgroup of size `M` and index
`[G:H]=m`, then writing `f_x(h)=f(xh)` the translate of `f` by representatives
`x` of `G/H`, one has a decomposition

.. math::

   \hat f(\chi) = \sum_{x\in G/H} \overline{\chi(x)} \hat{f_x}(\chi_{H})

so that the DFT on `G` can be computed using `m` DFT  on `H` (of
appropriate translates of `f`), then `M` DFT on `G/H`, one for
each restriction `\chi_{H}`.

This decomposition can be done recursively.

Naive algorithm
...............................................................................

.. function:: void acb_dft_naive(acb_ptr w, acb_srcptr v, slong n, slong prec)

   Computes the DFT of *v* into *w*, where *v* and *w* have size *n*,
   using the naive `O(n^2)` algorithm.

.. type:: acb_dft_naive_struct

.. type:: acb_dft_naive_t

.. function:: void acb_dft_naive_init(acb_dft_naive_t t, slong len, slong prec)

.. function:: void acb_dft_naive_clear(acb_dft_naive_t t)

   Stores a table of roots of unity in *t*.
   The length *len* is stored as *t->n*.

.. function:: void acb_dft_naive_precomp(acb_ptr w, acb_srcptr v, const acb_dft_naive_t t, slong prec)

   Sets *w* to the DFT of *v* of size *t->n*, using the naive algorithm data *t*.

CRT decomposition
...............................................................................

.. function:: void acb_dft_crt(acb_ptr w, acb_srcptr v, slong n, slong prec)

   Computes the DFT of *v* into *w*, where *v* and *w* have size *len*,
   using CRT to express `\mathbb Z/n\mathbb Z` as a product of cyclic groups.

.. type:: acb_dft_crt_struct

.. type:: acb_dft_crt_t

.. function:: void acb_dft_crt_init(acb_dft_crt_t t, slong len, slong prec)

.. function:: void acb_dft_crt_clear(acb_dft_crt_t t)

   Initialize a CRT decomposition of `\mathbb Z/n\mathbb Z` as a direct product
   of cyclic groups.
   The length *len* is stored as *t->n*.

.. function:: void acb_dft_crt_precomp(acb_ptr w, acb_srcptr v, const acb_dft_crt_t t, slong prec)

   Sets *w* to the DFT of *v* of size *t->n*, using the CRT decomposition scheme *t*.

Cooley-Tukey decomposition
...............................................................................

.. function:: void acb_dft_cyc(acb_ptr w, acb_srcptr v, slong n, slong prec)

   Computes the DFT of *v* into *w*, where *v* and *w* have size *n*,
   using each prime factor of `m` of `n` to decompose with
   the subgroup `H=m\mathbb Z/n\mathbb Z`.

.. type:: acb_dft_cyc_struct

.. type:: acb_dft_cyc_t

.. function:: void acb_dft_cyc_init(acb_dft_cyc_t t, slong len, slong prec)

.. function:: void acb_dft_cyc_clear(acb_dft_cyc_t t)

   Initialize a decomposition of `\mathbb Z/n\mathbb Z` into cyclic subgroups.
   The length *len* is stored as *t->n*.

.. function:: void acb_dft_cyc_precomp(acb_ptr w, acb_srcptr v, const acb_dft_cyc_t t, slong prec)

   Sets *w* to the DFT of *v* of size *t->n*, using the cyclic decomposition scheme *t*.

Radix 2 decomposition
...............................................................................

.. function:: void acb_dft_rad2(acb_ptr w, acb_srcptr v, int e, slong prec)

   Computes the DFT of *v* into *w*, where *v* and *w* have size `2^e`,
   using a radix 2 FFT.

.. function:: void acb_dft_inverse_rad2_precomp_inplace(acb_ptr v, const acb_dft_rad2_t t, slong prec)

   Computes the inverse DFT of *v* in place, of size *t->n*, using the
   precomputed radix 2 scheme *t*.

.. type:: acb_dft_rad2_struct

.. type:: acb_dft_rad2_t

.. function:: void acb_dft_rad2_init(acb_dft_rad2_t t, int e, slong prec)

.. function:: void acb_dft_rad2_clear(acb_dft_rad2_t t)

   Initialize and clear a radix 2 FFT of size `2^e`, stored as *t->n*.

.. function:: void acb_dft_rad2_precomp(acb_ptr w, acb_srcptr v, const acb_dft_rad2_t t, slong prec)

   Sets *w* to the DFT of *v* of size *t->n*, using the precomputed radix 2 scheme *t*.

Bluestein transform
...............................................................................

.. function:: void acb_dft_bluestein(acb_ptr w, acb_srcptr v, slong n, slong prec)

   Computes the DFT of *v* into *w*, where *v* and *w* have size *n*,
   by conversion to a radix 2 one using Bluestein's convolution trick.

.. type:: acb_dft_bluestein_struct

.. type:: acb_dft_bluestein_t

   Stores a Bluestein scheme for some length *n* : that is a :type:`acb_dft_rad2_t` of size
   `2^e \geq 2n-1` and a size *n* array of convolution factors.

.. function:: void acb_dft_bluestein_init(acb_dft_bluestein_t t, slong len, slong prec)
              void acb_dft_bluestein_clear(acb_dft_bluestein_t t)

   Initialize and clear a Bluestein scheme to compute DFT of size *len*.

.. function:: void acb_dft_bluestein_precomp(acb_ptr w, acb_srcptr v, const acb_dft_bluestein_t t, slong prec)

   Sets *w* to the DFT of *v* of size *t->n*, using the precomputed Bluestein scheme *t*.

