International
Tables for
Crystallography
Volume B
Reciprocal space
Edited by U. Shmueli

International Tables for Crystallography (2006). Vol. B. ch. 1.3, p. 52   | 1 | 2 |

Section 1.3.3.2.2.3. The prime factor algorithm

G. Bricognea

a MRC Laboratory of Molecular Biology, Hills Road, Cambridge CB2 2QH, England, and LURE, Bâtiment 209D, Université Paris-Sud, 91405 Orsay, France

1.3.3.2.2.3. The prime factor algorithm

| top | pdf |

The CRT will now be used to factor the N-point DFT into a tensor product of d transforms, the jth of length [N_{j}].

Let the indices k and [k^{*}] be subjected to the following mappings:

  • (i) [k \;\longmapsto\; (k_{1}, k_{2}, \ldots, k_{d}), k_{j} \in {\bb Z}/N_{j} {\bb Z}], by [k_{j} \equiv k] mod [N_{j}] for each j, with reconstruction formula [k = {\textstyle\sum\limits_{i = 1}^{d}} \;k_{i} q_{i} Q_{i} \hbox{ mod } N\hbox{;}]

  • (ii) [k^{*} \;\longmapsto\; (k_{1}^{*}, k_{2}^{*}, \ldots, k_{d}^{*}), k_{j}^{*} \in {\bb Z}/N_{j} {\bb Z}], by [k_{j}^{*} \equiv q_{j} k^{*}] mod [N_{j}] for each j, with reconstruction formula [k^{*} = {\textstyle\sum\limits_{i = 1}^{d}} \;k_{i}^{*} Q_{i} \hbox{ mod } N.]

Then [\eqalign{k^{*} k &= \left({\textstyle\sum\limits_{i = 1}^{d}}\; k_{i}^{*} Q_{i}\right) \left({\textstyle\sum\limits_{j = 1}^{d}} \;k_{j} q_{j} Q_{j}\right) \hbox{ mod } N\cr &= {\textstyle\sum\limits_{i, \, j = 1}^{d}} k_{i}^{*} k_{j} Q_{i} q_{j} Q_{j} \hbox{ mod } N.}] Cross terms with [i \neq j] vanish since they contain all the factors of N, hence [\eqalign{k^{*} k &= {\textstyle\sum\limits_{j = 1}^{d}}\; q_{j} Q_{j}^{2} k_{j}^{*} k_{j} \hbox{ mod } N\cr &= {\textstyle\sum\limits_{j = 1}^{d}} (1 - n_{j} N_{j}) Q_{j} k_{j}^{*} k_{j} \hbox{ mod } N.}] Dividing by N, which may be written as [N_{j} Q_{j}] for each j, yields [\eqalign{{k^{*} k \over N} &= {\sum\limits_{j = 1}^{d}} (1 - n_{j} N_{j}) {Q_{j} \over N_{j} Q_{j}} k_{j}^{*} k_{j} \hbox{ mod } 1\cr &= {\sum\limits_{j = 1}^{d}} \left({1 \over N_{j}} - n_{j}\right) k_{j}^{*} k_{j} \hbox{ mod } 1,}] and hence [{k^{*} k \over N} \equiv {\sum\limits_{j = 1}^{d}} {k_{j}^{*} k_{j} \over N_{j}} \hbox{ mod } 1.] Therefore, by the multiplicative property of [e(.)], [e \left({k^{*} k \over N}\right) \equiv \bigotimes\limits_{j = 1}^{d} e \left({k_{j}^{*} k_{j} \over N_{j}}\right).]

Let [{\bf X} \in L ({\bb Z}/N {\bb Z})] be described by a one-dimensional array [X(k)] indexed by k. The index mapping (i)[link] turns X into an element of [L ({\bb Z}/N_{1} {\bb Z} \times \ldots \times {\bb Z}/N_{d} {\bb Z})] described by a d-dimensional array [X (k_{1}, \ldots, k_{d})]; the latter may be transformed by [\bar{F} (N_{1}) \bigotimes \ldots \bigotimes \bar{F} (N_{d})] into a new array [X^{*} (k_{1}^{*}, k_{2}^{*}, \ldots, k_{d}^{*})]. Finally, the one-dimensional array of results [X^{*} (k^{*})] will be obtained by reconstructing [k^{*}] according to (ii)[link].

The prime factor algorithm, like the Cooley–Tukey algorithm, reindexes a 1D transform to turn it into d separate transforms, but the use of coprime factors and CRT index mapping leads to the further gain that no twiddle factors need to be applied between the successive transforms (see Good, 1971[link]). This makes up for the cost of the added complexity of the CRT index mapping.

The natural factorization of N for the prime factor algorithm is thus its factorization into prime powers: [\bar{F}(N)] is then the tensor product of separate transforms (one for each prime power factor [N_{j} = p_{j}^{\nu_{j}}]) whose results can be reassembled without twiddle factors. The separate factors [p_{j}] within each [N_{j}] must then be dealt with by another algorithm (e.g. Cooley–Tukey, which does require twiddle factors). Thus, the DFT on a prime number of points remains undecomposable.

References

First citation Good, I. J. (1971). The relationship between two fast Fourier transforms. IEEE Trans. C-20, 310–317.Google Scholar








































to end of page
to top of page