NAG Library Function Document

nag_fft_multiple_real (c06fpc)


    1  Purpose
    7  Accuracy


nag_fft_multiple_real (c06fpc) computes the discrete Fourier transforms of m  sequences, each containing n  real data values.


#include <nag.h>
#include <nagc06.h>
void  nag_fft_multiple_real (Integer m, Integer n, double x[], const double trig[], NagError *fail)


Given m  sequences of n  real data values x j p , for j=0,1,,n - 1 and p=1,2,,m, this function simultaneously calculates the Fourier transforms of all the sequences defined by
z ^ k p = 1 n j=0 n-1 x j p exp -2 π ijk / n ,   for ​ k = 0 , 1 , , n - 1 ; p = 1 , 2 , , m .  
(Note the scale factor 1 / n  in this definition.)
The transformed values z ^ k p  are complex, but for each value of p  the z ^ k p  form a Hermitian sequence (i.e., z ^ n-k p  is the complex conjugate of z ^ k p ), so they are completely determined by mn  real numbers. The first call of nag_fft_multiple_real (c06fpc) must be preceded by a call to nag_fft_init_trig (c06gzc) to initialize the array trig with trigonometric coefficients according to the value of n.
The discrete Fourier transform is sometimes defined using a positive sign in the exponential term
z ^ k p = 1 n j=0 n-1 x j p exp +2 π ijk / n .  
To compute this form, this function should be followed by a call to nag_multiple_conjugate_hermitian (c06gqc) to form the complex conjugates of the z ^ k p .
The function uses a variant of the fast Fourier transform algorithm (Brigham (1974)) known as the Stockham self-sorting algorithm, which is described in Temperton (1983). Special coding is provided for the factors 2, 3, 4, 5 and 6.


Brigham E O (1974) The Fast Fourier Transform Prentice–Hall
Temperton C (1983) Fast mixed-radix real Fourier transforms J. Comput. Phys. 52 340–350


1:     m IntegerInput
On entry: the number of sequences to be transformed, m .
Constraint: m1 .
2:     n IntegerInput
On entry: the number of real values in each sequence, n .
Constraint: n1 .
3:     x[m×n] doubleInput/Output
On entry: the m  data sequences must be stored in x consecutively. If the data values of the p th sequence to be transformed are denoted by x j p , for j=0,1,,n - 1, then the mn  elements of the array x must contain the values
x 0 1 , x 1 1 , , x n-1 1 , x 0 2 , x 1 2 , , x n-1 2 , , x 0 m , x 1 m , , x n-1 m .  
On exit: the m  discrete Fourier transforms in Hermitian form, stored consecutively, overwriting the corresponding original sequences. If the n  components of the discrete Fourier transform z ^ k p  are written as a k p + ib k p , then for 0 k n / 2 , a k p  is in array element x[ p-1 × n + k ]  and for 1 k n-1 / 2 , b k p  is in array element x[ p-1 × n + n-k ] .
4:     trig[2×n] const doubleInput
On entry: trigonometric coefficients as returned by a call of nag_fft_init_trig (c06gzc). nag_fft_multiple_real (c06fpc) makes a simple check to ensure that trig has been initialized and that the initialization is compatible with the value of n
5:     fail NagError *Input/Output
The NAG error argument (see Section 3.7 in How to Use the NAG Library and its Documentation).

Error Indicators and Warnings

Dynamic memory allocation failed.
Value of n and trig array are incompatible or trig array not initialized.
On entry, m=value.
Constraint: m1.
On entry, n=value.
Constraint: n1.


Some indication of accuracy can be obtained by performing a subsequent inverse transform and comparing the results with the original sequence (in exact arithmetic they would be identical).

Parallelism and Performance

nag_fft_multiple_real (c06fpc) is not threaded in any implementation.

Further Comments

The time taken is approximately proportional to nm logn , but also depends on the factors of n . The function is fastest if the only prime factors of n  are 2, 3 and 5, and is particularly slow if n  is a large prime, or has large prime factors.


This program reads in sequences of real data values and prints their discrete Fourier transforms (as computed by nag_fft_multiple_real (c06fpc)). The Fourier transforms are expanded into full complex form using nag_multiple_hermitian_to_complex (c06gsc) and printed. Inverse transforms are then calculated by calling nag_multiple_conjugate_hermitian (c06gqc) followed by nag_fft_multiple_hermitian (c06fqc) showing that the original sequences are restored.

Program Text

Program Text (c06fpce.c)

Program Data

Program Data (c06fpce.d)

Program Results

Program Results (c06fpce.r)

© The Numerical Algorithms Group Ltd, Oxford, UK. 2017