This is fftw3.info, produced by makeinfo version 4.2 from fftw3.texi. This manual is for FFTW (version 3.0.1, 5 June 2003). Copyright (C) 2003 Matteo Frigo. Copyright (C) 2003 Massachusetts Institute of Technology. Permission is granted to make and distribute verbatim copies of this manual provided the copyright notice and this permission notice are preserved on all copies. Permission is granted to copy and distribute modified versions of this manual under the conditions for verbatim copying, provided that the entire resulting derived work is distributed under the terms of a permission notice identical to this one. Permission is granted to copy and distribute translations of this manual into another language, under the above conditions for modified versions, except that this permission notice may be stated in a translation approved by the Free Software Foundation. INFO-DIR-SECTION Texinfo documentation system START-INFO-DIR-ENTRY * fftw3: (fftw3). FFTW User's Manual. END-INFO-DIR-ENTRY  File: fftw3.info, Node: Dynamic Arrays in C, Next: Dynamic Arrays in C-The Wrong Way, Prev: Static Arrays in C, Up: Multi-dimensional Array Format Dynamic Arrays in C ------------------- We recommend allocating most arrays dynamically, with `fftw_malloc'. This isn't too hard to do, although it is not as straightforward for multi-dimensional arrays as it is for one-dimensional arrays. Creating the array is simple: using a dynamic-allocation routine like `fftw_malloc', allocate an array big enough to store N `fftw_complex' values (for a complex DFT), where N is the product of the sizes of the array dimensions (i.e. the total number of complex values in the array). For example, here is code to allocate a 5x12x27 rank 3 array: fftw_complex *an_array; an_array = fftw_malloc(5*12*27 * sizeof(fftw_complex)); Accessing the array elements, however, is more tricky--you can't simply use multiple applications of the `[]' operator like you could for static arrays. Instead, you have to explicitly compute the offset into the array using the formula given earlier for row-major arrays. For example, to reference the (i,j,k)-th element of the array allocated above, you would use the expression `an_array[k + 27 * (j + 12 * i)]'. This pain can be alleviated somewhat by defining appropriate macros, or, in C++, creating a class and overloading the `()' operator. The recent C99 standard provides a way to dynamically reinterpret the array as a static-like multi-dimensional array amenable to `[]', but this feature is not yet widely supported by compilers.  File: fftw3.info, Node: Dynamic Arrays in C-The Wrong Way, Prev: Dynamic Arrays in C, Up: Multi-dimensional Array Format Dynamic Arrays in C--The Wrong Way ---------------------------------- A different method for allocating multi-dimensional arrays in C is often suggested that is incompatible with FFTW: _using it will cause FFTW to die a painful death_. We discuss the technique here, however, because it is so commonly known and used. This method is to create arrays of pointers of arrays of pointers of ...etcetera. For example, the analogue in this method to the example above is: int i,j; fftw_complex ***a_bad_array; /* another way to make a 5x12x27 array */ a_bad_array = (fftw_complex ***) malloc(5 * sizeof(fftw_complex **)); for (i = 0; i < 5; ++i) { a_bad_array[i] = (fftw_complex **) malloc(12 * sizeof(fftw_complex *)); for (j = 0; j < 12; ++j) a_bad_array[i][j] = (fftw_complex *) malloc(27 * sizeof(fftw_complex)); } As you can see, this sort of array is inconvenient to allocate (and deallocate). On the other hand, it has the advantage that the (i,j,k)-th element can be referenced simply by `a_bad_array[i][j][k]'. If you like this technique and want to maximize convenience in accessing the array, but still want to pass the array to FFTW, you can use a hybrid method. Allocate the array as one contiguous block, but also declare an array of arrays of pointers that point to appropriate places in the block. That sort of trick is beyond the scope of this documentation; for more information on multi-dimensional arrays in C, see the `comp.lang.c' FAQ (http://www.eskimo.com/~scs/C-faq/s6.html).  File: fftw3.info, Node: Words of Wisdom-Saving Plans, Next: Caveats in Using Wisdom, Prev: Multi-dimensional Array Format, Up: Other Important Topics Words of Wisdom--Saving Plans ============================= FFTW implements a method for saving plans to disk and restoring them. In fact, what FFTW does is more general than just saving and loading plans. The mechanism is called "wisdom". Here, we describe this feature at a high level. *Note FFTW Reference::, for a less casual but more complete discussion of how to use wisdom in FFTW. Plans created with the `FFTW_MEASURE', `FFTW_PATIENT', or `FFTW_EXHAUSTIVE' options produce near-optimal FFT performance, but may require a long time to compute because FFTW must measure the runtime of many possible plans and select the best one. This setup is designed for the situations where so many transforms of the same size must be computed that the start-up time is irrelevant. For short initialization times, but slower transforms, we have provided `FFTW_ESTIMATE'. The `wisdom' mechanism is a way to get the best of both worlds: you compute a good plan once, save it to disk, and later reload it as many times as necessary. The wisdom mechanism can actually save and reload many plans at once, not just one. Whenever you create a plan, the FFTW planner accumulates wisdom, which is information sufficient to reconstruct the plan. After planning, you can save this information to disk by means of the function: void fftw_export_wisdom_to_file(FILE *output_file); The next time you run the program, you can restore the wisdom with `fftw_import_wisdom_from_file' (which returns non-zero on success), and then recreate the plan using the same flags as before. int fftw_import_wisdom_from_file(FILE *output_file); Wisdom is automatically used for any size to which it is applicable, as long as the planner flags are not more "patient" than those with which the wisdom was created. For example, wisdom created with `FFTW_MEASURE' can be used if you later plan with `FFTW_ESTIMATE' or `FFTW_MEASURE', but not with `FFTW_PATIENT'. The `wisdom' is cumulative, and is stored in a global, private data structure managed internally by FFTW. The storage space required is minimal, proportional to the logarithm of the sizes the wisdom was generated from. If memory usage is a concern, however, the wisdom can be forgotten and its associated memory freed by calling: void fftw_forget_wisdom(void); Wisdom can be exported to a file, a string, or any other medium. For details, see *Note Wisdom::.  File: fftw3.info, Node: Caveats in Using Wisdom, Prev: Words of Wisdom-Saving Plans, Up: Other Important Topics Caveats in Using Wisdom ======================= For in much wisdom is much grief, and he that increaseth knowledge increaseth sorrow. [Ecclesiastes 1:18] There are pitfalls to using wisdom, in that it can negate FFTW's ability to adapt to changing hardware and other conditions. For example, it would be perfectly possible to export wisdom from a program running on one processor and import it into a program running on another processor. Doing so, however, would mean that the second program would use plans optimized for the first processor, instead of the one it is running on. It should be safe to reuse wisdom as long as the hardware and program binaries remain unchanged. (Actually, the optimal plan may change even between runs of the same binary on identical hardware, due to differences in the virtual memory environment, etcetera. Users seriously interested in performance should worry about this problem, too.) It is likely that, if the same wisdom is used for two different program binaries, even running on the same machine, the plans may be sub-optimal because of differing code alignments. It is therefore wise to recreate wisdom every time an application is recompiled. The more the underlying hardware and software changes between the creation of wisdom and its use, the greater grows the risk of sub-optimal plans. Nevertheless, if the choice is between using `FFTW_ESTIMATE' or using possibly-suboptimal wisdom (created on the same machine, but for a different binary), the wisdom is likely to be better. For this reason, we provide a function to import wisdom from a standard system-wide location (`/etc/fftw/wisdom' on Unix): int fftw_import_system_wisdom(void); FFTW also provides a standalone program, `fftw-wisdom' (described by its own `man' page on Unix) with which users can create wisdom, e.g. for a canonical set of sizes to store in the system wisdom file. *Note Wisdom Utilities::.  File: fftw3.info, Node: FFTW Reference, Next: Parallel FFTW, Prev: Other Important Topics, Up: Top FFTW Reference ************** This chapter provides a complete reference for all sequential (i.e., one-processor) FFTW functions. For parallel transforms, *Note Parallel FFTW::. * Menu: * Data Types and Files:: * Using Plans:: * Basic Interface:: * Advanced Interface:: * Guru Interface:: * Wisdom:: * What FFTW Really Computes::  File: fftw3.info, Node: Data Types and Files, Next: Using Plans, Prev: FFTW Reference, Up: FFTW Reference Data Types and Files ==================== All programs using FFTW should include its header file: #include You must also link to the FFTW library. On Unix, this means adding `-lfftw3 -lm' at the _end_ of the link command. * Menu: * Complex numbers:: * Precision:: * Memory Allocation::  File: fftw3.info, Node: Complex numbers, Next: Precision, Prev: Data Types and Files, Up: Data Types and Files Complex numbers --------------- The default FFTW interface uses `double' precision for all floating-point numbers, and defines a `fftw_complex' type to hold complex numbers as: typedef double fftw_complex[2]; Here, the `[0]' element holds the real part and the `[1]' element holds the imaginary part. Alternatively, if you have a C compiler (such as `gcc') that supports the C99 revision of the ANSI C standard, you can use C's new native complex type (which is binary-compatible with the typedef above). In particular, if you `#include ' _before_ `', then `fftw_complex' is defined to be the native complex type and you can manipulate it with ordinary arithmetic (e.g. `x = y * (3+4*I)', where `x' and `y' are `fftw_complex' and `I' is the standard symbol for the imaginary unit); C++ has its own `complex' template class, defined in the standard `' header file. Reportedly, the C++ standards committee has recently agreed to mandate that the storage format used for this type be binary-compatible with the C99 type, i.e. an array `T[2]' with consecutive real `[0]' and imaginary `[1]' parts. (See report WG21/N1388 (http://anubis.dkuug.dk/JTC1/SC22/WG21/docs/papers/2002/1388.pdf).) Although not part of the official standard as of this writing, the proposal stated that: "This solution has been tested with all current major implementations of the standard library and shown to be working." To the extent that this is true, if you have a variable `complex *x', you can pass it directly to FFTW via `reinterpret_cast(x)'.  File: fftw3.info, Node: Precision, Next: Memory Allocation, Prev: Complex numbers, Up: Data Types and Files Precision --------- You can install single and long-double precision versions of FFTW, which replace `double' with `float' and `long double', respectively (*note Installation and Customization::). To use these interfaces, you: * Link to the single/long-double libraries; on Unix, `-lfftw3f' or `-lfftw3l' instead of (or in addition to) `-lfftw3'. (You can link to the different-precision libraries simultaneously.) * Include the _same_ `' header file. * Replace all lowercase instances of `fftw_' with `fftwf_' or `fftwl_' for single or long-double precision, respectively. (`fftw_complex' becomes `fftwf_complex', `fftw_execute' becomes `fftwf_execute', etcetera.) * Uppercase names, i.e. names beginning with `FFTW_', remain the same. * Replace `double' with `float' or `long double' for subroutine parameters. Depending upon your compiler and/or hardware, `long double' may not be any more precise than `double' (or may not be supported at all, although it is standard in C99).  File: fftw3.info, Node: Memory Allocation, Prev: Precision, Up: Data Types and Files Memory Allocation ----------------- void *fftw_malloc(size_t n); void fftw_free(void *p); These are functions that behave identically to `malloc' and `free', except that they guarantee that the returned pointer obeys any special alignment restrictions imposed by any algorithm in FFTW (e.g. for SIMD acceleration). *Note Data Alignment::. Data allocated by `fftw_malloc' _must_ be deallocated by `fftw_free' and not by the ordinary `free'. These routines simply call through to your operating system's `malloc' or, if necessary, its aligned equivalent (e.g. `memalign'), so you normally need not worry about any significant time or space overhead. You are _not required_ to use them to allocate your data, but we strongly recommend it.  File: fftw3.info, Node: Using Plans, Next: Basic Interface, Prev: Data Types and Files, Up: FFTW Reference Using Plans =========== Plans for all transform types in FFTW are stored as type `fftw_plan' (an opaque pointer type), and are created by one of the various planning routines described in the following sections. An `fftw_plan' contains all information necessary to compute the transform, including the pointers to the input and output arrays. void fftw_execute(const fftw_plan plan); This executes the `plan', to compute the corresponding transform on the arrays for which it was planned (which must still exist). The plan is not modified, and `fftw_execute' can be called as many times as desired. To apply a given plan to a different array, you can use the guru interface. *Note Guru Interface::. `fftw_execute' (and equivalents) is the only function in FFTW guaranteed to be thread-safe; see *Note Thread safety::. This function: void fftw_destroy_plan(fftw_plan plan); deallocates the `plan' and all its associated data. FFTW's planner saves some other persistent data, such as the accumulated wisdom and a list of algorithms available in the current configuration. If you want to deallocate all of that and reset FFTW to the pristine state it was in when you started your program, you can call: void fftw_cleanup(void); This does not deallocate your plans; you should still call `fftw_destroy_plan' if you want to do this. You should not execute any previously created plans after calling `fftw_cleanup', however. The following two routines are provided purely for academic purposes (that is, for entertainment). void fftw_flops(const fftw_plan plan, double *add, double *mul, double *fma); Given a `plan', set `add', `mul', and `fma' to an exact count of the number of floating-point additions, multiplications, and fused multiply-add operations involved in the plan's execution. The total number of floating-point operations (flops) is `add + mul + 2*fma', or `add + mul + fma' if the hardware supports fused multiply-add instructions (although the number of FMA operations is only approximate because of compiler voodoo). (The number of operations should be an integer, but we use `double' to avoid overflowing `int' for large transforms; the arguments are of type `double' even for single and long-double precision versions of FFTW.) void fftw_fprint_plan(const fftw_plan plan, FILE *output_file); void fftw_print_plan(const fftw_plan plan); This outputs a "nerd-readable" representation of the `plan' to the given file or to `stdout', respectively.  File: fftw3.info, Node: Basic Interface, Next: Advanced Interface, Prev: Using Plans, Up: FFTW Reference Basic Interface =============== The basic interface, which we expect to satisfy the needs of most users, provides planner routines for transforms of a single contiguous array with any of FFTW's supported transform types. * Menu: * Complex DFTs:: * Planner Flags:: * Real-data DFTs:: * Real-data DFT Array Format:: * Real-to-Real Transforms:: * Real-to-Real Transform Kinds::  File: fftw3.info, Node: Complex DFTs, Next: Planner Flags, Prev: Basic Interface, Up: Basic Interface Complex DFTs ------------ fftw_plan fftw_plan_dft_1d(int n, fftw_complex *in, fftw_complex *out, int sign, unsigned flags); fftw_plan fftw_plan_dft_2d(int nx, int ny, fftw_complex *in, fftw_complex *out, int sign, unsigned flags); fftw_plan fftw_plan_dft_3d(int nx, int ny, int nz, fftw_complex *in, fftw_complex *out, int sign, unsigned flags); fftw_plan fftw_plan_dft(int rank, const int *n, fftw_complex *in, fftw_complex *out, int sign, unsigned flags); Plan a complex input/output discrete Fourier transform (DFT) in zero or more dimensions, returning an `fftw_plan' (*note Using Plans::). Once you have created a plan for a certain transform type and parameters, then creating another plan of the same type and parameters, but for different arrays, is fast and shares constant data with the first plan (if it still exists). The planner returns `NULL' if the plan cannot be created. A non-`NULL' plan is always returned by the basic interface unless you are using a customized FFTW configuration supporting a restricted set of transforms. Arguments ......... * `rank' is the dimensionality of the transform (it should be the size of the array `*n'), and can be any non-negative integer. The `_1d', `_2d', and `_3d' planners correspond to a `rank' of `1', `2', and `3', respectively. A `rank' of zero is equivalent to a transform of size 1, i.e. a copy of one number from input to output. * `n', or `nx'/`ny'/`nz', or `n[rank]', respectively, gives the size of the transform dimensions. They can be any positive integer. - Multi-dimensional arrays are stored in row-major order with dimensions: `nx' x `ny'; or `nx' x `ny' x `nz'; or `n[0]' x `n[1]' x ... x `n[rank-1]'. *Note Multi-dimensional Array Format::. - FFTW is best at handling sizes of the form 2^a 3^b 5^c 7^d 11^e 13^f, where e+f is either 0 or 1, and the other exponents are arbitrary. Other sizes are computed by means of a slow, general-purpose algorithm (which nevertheless retains O(n log n) performance even for prime sizes). It is possible to customize FFTW for different array sizes; see *Note Installation and Customization::. Transforms whose sizes are powers of 2 are especially fast. * `in' and `out' point to the input and output arrays of the transform, which may be the same (yielding an in-place transform). These arrays are overwritten during planning, unless `FFTW_ESTIMATE' is used in the flags. (The arrays need not be initialized, but they must be allocated.) If `in == out', the transform is "in-place" and the input array is overwritten. If `in != out', the two arrays must not overlap (but FFTW does not check for this condition). * `sign' is the sign of the exponent in the formula that defines the Fourier transform. It can be -1 (= `FFTW_FORWARD') or +1 (= `FFTW_BACKWARD'). * `flags' is a bitwise OR (`|') of zero or more planner flags, as defined in *Note Planner Flags::. FFTW computes an unnormalized transform: computing a forward followed by a backward transform (or vice versa) will result in the original data multiplied by the size of the transform (the product of the dimensions). For more information, see *Note What FFTW Really Computes::.  File: fftw3.info, Node: Planner Flags, Next: Real-data DFTs, Prev: Complex DFTs, Up: Basic Interface Planner Flags ------------- All of the planner routines in FFTW accept an integer `flags' argument, which is a bitwise OR (`|') of zero or more of the flag constants defined below. These flags control the rigor (and time) of the planning process, and can also impose (or lift) restrictions on the type of transform algorithm that is employed. Planning-rigor flags .................... * `FFTW_ESTIMATE' specifies that, instead of actual measurements of different algorithms, a simple heuristic is used to pick a (probably sub-optimal) plan quickly. With this flag, the input/output arrays are not overwritten during planning. * `FFTW_MEASURE' tells FFTW to find an optimized plan by actually _computing_ several FFTs and measuring their execution time. Depending on your machine, this can take some time (often a few seconds). `FFTW_MEASURE' is the default planning option. * `FFTW_PATIENT' is like `FFTW_MEASURE', but considers a wider range of algorithms and often produces a "more optimal" plan (especially for large transforms), but at the expense of several times longer planning time (especially for large transforms). * `FFTW_EXHAUSTIVE' is like `FFTW_PATIENT', but considers an even wider range of algorithms, including many that we think are unlikely to be fast, to produce the most optimal plan but with a substantially increased planning time. Algorithm-restriction flags ........................... * `FFTW_DESTROY_INPUT' specifies that an out-of-place transform is allowed to _overwrite its input_ array with arbitrary data; this can sometimes allow more efficient algorithms to be employed. * `FFTW_PRESERVE_INPUT' specifies that an out-of-place transform must _not change its input_ array. This is ordinarily the _default_, except for c2r and hc2r (i.e. complex-to-real) transforms for which `FFTW_DESTROY_INPUT' is the default. In the latter cases, passing `FFTW_PRESERVE_INPUT' will attempt to use algorithms that do not destroy the input, at the expense of worse performance; for multi-dimensional c2r transforms, however, no input-preserving algorithms are implemented and the planner will return `NULL' if one is requested. * `FFTW_UNALIGNED' specifies that the algorithm may not impose any unusual alignment requirements on the input/output arrays (i.e. no SIMD may be used). This flag is normally _not necessary_, since the planner automatically detects misaligned arrays. The only use for this flag is if you want to use the guru interface to execute a given plan on a different array that may not be aligned like the original. (Using `fftw_malloc' makes this flag unnecessary even then.)  File: fftw3.info, Node: Real-data DFTs, Next: Real-data DFT Array Format, Prev: Planner Flags, Up: Basic Interface Real-data DFTs -------------- fftw_plan fftw_plan_dft_r2c_1d(int n, double *in, fftw_complex *out, unsigned flags); fftw_plan fftw_plan_dft_r2c_2d(int nx, int ny, double *in, fftw_complex *out, unsigned flags); fftw_plan fftw_plan_dft_r2c_3d(int nx, int ny, int nz, double *in, fftw_complex *out, unsigned flags); fftw_plan fftw_plan_dft_r2c(int rank, const int *n, double *in, fftw_complex *out, unsigned flags); Plan a real-input/complex-output discrete Fourier transform (DFT) in zero or more dimensions, returning an `fftw_plan' (*note Using Plans::). Once you have created a plan for a certain transform type and parameters, then creating another plan of the same type and parameters, but for different arrays, is fast and shares constant data with the first plan (if it still exists). The planner returns `NULL' if the plan cannot be created. A non-`NULL' plan is always returned by the basic interface unless you are using a customized FFTW configuration supporting a restricted set of transforms, or if you use the `FFTW_PRESERVE_INPUT' flag with a multi-dimensional out-of-place c2r transform (see below). Arguments ......... * `rank' is the dimensionality of the transform (it should be the size of the array `*n'), and can be any non-negative integer. The `_1d', `_2d', and `_3d' planners correspond to a `rank' of `1', `2', and `3', respectively. A `rank' of zero is equivalent to a transform of size 1, i.e. a copy of one number (with zero imaginary part) from input to output. * `n', or `nx'/`ny'/`nz', or `n[rank]', respectively, gives the size of the _logical_ transform dimensions. They can be any positive integer. This is different in general from the _physical_ array dimensions, which are described in *Note Real-data DFT Array Format::. - FFTW is best at handling sizes of the form 2^a 3^b 5^c 7^d 11^e 13^f, where e+f is either 0 or 1, and the other exponents are arbitrary. Other sizes are computed by means of a slow, general-purpose algorithm (which nevertheless retains O(n log n) performance even for prime sizes). (It is possible to customize FFTW for different array sizes; see *Note Installation and Customization::.) Transforms whose sizes are powers of 2 are especially fast, and it is generally beneficial for the _last_ dimension of an r2c/c2r transform to be _even_. * `in' and `out' point to the input and output arrays of the transform, which may be the same (yielding an in-place transform). These arrays are overwritten during planning, unless `FFTW_ESTIMATE' is used in the flags. (The arrays need not be initialized, but they must be allocated.) For an in-place transform, it is important to remember that the real array will require padding, described in *Note Real-data DFT Array Format::. * `flags' is a bitwise OR (`|') of zero or more planner flags, as defined in *Note Planner Flags::. The inverse transforms, taking complex input (storing the non-redundant half of a logically Hermitian array) to real output, are given by: fftw_plan fftw_plan_dft_c2r_1d(int n, fftw_complex *in, double *out, unsigned flags); fftw_plan fftw_plan_dft_c2r_2d(int nx, int ny, fftw_complex *in, double *out, unsigned flags); fftw_plan fftw_plan_dft_c2r_3d(int nx, int ny, int nz, fftw_complex *in, double *out, unsigned flags); fftw_plan fftw_plan_dft_c2r(int rank, const int *n, fftw_complex *in, double *out, unsigned flags); The arguments are the same as for the r2c transforms, except that the input and output data formats are reversed. FFTW computes an unnormalized transform: computing an r2c followed by a c2r transform (or vice versa) will result in the original data multiplied by the size of the transform (the product of the logical dimensions). An r2c transform produces the same output as a `FFTW_FORWARD' complex DFT of the same input, and a c2r transform is correspondingly equivalent to `FFTW_BACKWARD'. For more information, see *Note What FFTW Really Computes::.  File: fftw3.info, Node: Real-data DFT Array Format, Next: Real-to-Real Transforms, Prev: Real-data DFTs, Up: Basic Interface Real-data DFT Array Format -------------------------- The output of a DFT of real data (r2c) contains symmetries that, in principle, make half of the outputs redundant (*note What FFTW Really Computes::). (Similarly for the input of an inverse c2r transform.) In practice, it is not possible to entirely realize these savings in an efficient and understandable format that generalizes to multi-dimensional transforms. Instead, the output of the r2c transforms is _slightly_ over half of the output of the corresponding complex transform. We do not "pack" the data in any way, but store it as an ordinary array of `fftw_complex' values. In fact, this data is simply a subsection of what would be the array in the corresponding complex transform. Specifically, for a real transform of d (= `rank') dimensions n1 x n2 x n3 x ... x nd , the complex data is an n1 x n2 x n3 x ... x (nd/2 + 1) array of `fftw_complex' values in row-major order (with the division rounded down). That is, we only store the _lower_ half (non-negative frequencies), plus one element, of the last dimension of the data from the ordinary complex transform. (We could have instead taken half of any other dimension, but implementation turns out to be simpler if the last, contiguous, dimension is used.) For an out-of-place transform, the real data is simply an array with physical dimensions n1 x n2 x n3 x ... x nd in row-major order. For an in-place transform, some complications arise since the complex data is slightly larger than the real data. In this case, the final dimension of the real data must be _padded_ with extra values to accommodate the size of the complex data--two extra if the last dimension is even and one if it is odd. That is, the last dimension of the real data must physically contain 2 * (nd/2+1) `double' values (exactly enough to hold the complex data). This physical array size does not, however, change the _logical_ array size--only nd values are actually stored in the last dimension, and nd is the last dimension passed to the planner.  File: fftw3.info, Node: Real-to-Real Transforms, Next: Real-to-Real Transform Kinds, Prev: Real-data DFT Array Format, Up: Basic Interface Real-to-Real Transforms ----------------------- fftw_plan fftw_plan_r2r_1d(int n, double *in, double *out, fftw_r2r_kind kind, unsigned flags); fftw_plan fftw_plan_r2r_2d(int nx, int ny, double *in, double *out, fftw_r2r_kind kindx, fftw_r2r_kind kindy, unsigned flags); fftw_plan fftw_plan_r2r_3d(int nx, int ny, int nz, double *in, double *out, fftw_r2r_kind kindx, fftw_r2r_kind kindy, fftw_r2r_kind kindz, unsigned flags); fftw_plan fftw_plan_r2r(int rank, const int *n, double *in, double *out, const fftw_r2r_kind *kind, unsigned flags); Plan a real input/output (r2r) transform of various kinds in zero or more dimensions, returning an `fftw_plan' (*note Using Plans::). Once you have created a plan for a certain transform type and parameters, then creating another plan of the same type and parameters, but for different arrays, is fast and shares constant data with the first plan (if it still exists). The planner returns `NULL' if the plan cannot be created. A non-`NULL' plan is always returned by the basic interface unless you are using a customized FFTW configuration supporting a restricted set of transforms, or for size-1 `FFTW_REDFT00' kinds (which are not defined). Arguments ......... * `rank' is the dimensionality of the transform (it should be the size of the arrays `*n' and `*kind'), and can be any non-negative integer. The `_1d', `_2d', and `_3d' planners correspond to a `rank' of `1', `2', and `3', respectively. A `rank' of zero is equivalent to a copy of one number from input to output. * `n', or `nx'/`ny'/`nz', or `n[rank]', respectively, gives the (physical) size of the transform dimensions. They can be any positive integer. - Multi-dimensional arrays are stored in row-major order with dimensions: `nx' x `ny'; or `nx' x `ny' x `nz'; or `n[0]' x `n[1]' x ... x `n[rank-1]'. *Note Multi-dimensional Array Format::. - FFTW is generally best at handling sizes of the form 2^a 3^b 5^c 7^d 11^e 13^f, where e+f is either 0 or 1, and the other exponents are arbitrary. Other sizes are computed by means of a slow, general-purpose algorithm (which nevertheless retains O(n log n) performance even for prime sizes). (It is possible to customize FFTW for different array sizes; see *Note Installation and Customization::.) Transforms whose sizes are powers of 2 are especially fast. - For a `REDFT00' or `RODFT00' transform kind in a dimension of size n, it is n-1 or n+1, respectively, that should be factorizable in the above form. * `in' and `out' point to the input and output arrays of the transform, which may be the same (yielding an in-place transform). These arrays are overwritten during planning, unless `FFTW_ESTIMATE' is used in the flags. (The arrays need not be initialized, but they must be allocated.) * `kind', or `kindx'/`kindy'/`kindz', or `kind[rank]', is the kind of r2r transform used for the corresponding dimension. The valid kind constants are described in *Note Real-to-Real Transform Kinds::. In a multi-dimensional transform, what is computed is the separable product formed by taking each transform kind along the corresponding dimension, one dimension after another. * `flags' is a bitwise OR (`|') of zero or more planner flags, as defined in *Note Planner Flags::.  File: fftw3.info, Node: Real-to-Real Transform Kinds, Prev: Real-to-Real Transforms, Up: Basic Interface Real-to-Real Transform Kinds ---------------------------- FFTW currently supports 11 different r2r transform kinds, specified by one of the constants below. For the precise definitions of these transforms, see *Note What FFTW Really Computes::. For a more colloquial introduction to these transform kinds, see *Note More DFTs of Real Data::. For dimension of size `n', there is a corresponding "logical" dimension `N' that determines the normalization (and the optimal factorization); the formula for `N' is given for each kind below. Also, with each transform kind is listed its corrsponding inverse transform. FFTW computes unnormalized transforms: a transform followed by its inverse will result in the original data multiplied by `N' (or the product of the `N''s for each dimension, in multi-dimensions). * `FFTW_R2HC' computes a real-input DFT with output in "halfcomplex" format, i.e. real and imaginary parts for a transform of size `n' stored as: r0, r1, r2, r(n/2), i((n+1)/2-1), ..., i2, i1 (Logical `N=n', inverse is `FFTW_HC2R'.) * `FFTW_HC2R' computes the reverse of `FFTW_R2HC', above. (Logical `N=n', inverse is `FFTW_R2HC'.) * `FFTW_DHT' computes a discrete Hartley transform. (Logical `N=n', inverse is `FFTW_DHT'.) * `FFTW_REDFT00' computes an REDFT00 transform, i.e. a DCT-I. (Logical `N=2*(n-1)', inverse is `FFTW_REDFT00'.) * `FFTW_REDFT10' computes an REDFT10 transform, i.e. a DCT-II. (Logical `N=2*n', inverse is `FFTW_REDFT01'.) * `FFTW_REDFT01' computes an REDFT01 transform, i.e. a DCT-III. (Logical `N=2*n', inverse is `FFTW_REDFT=10'.) * `FFTW_REDFT11' computes an REDFT11 transform, i.e. a DCT-IV. (Logical `N=2*n', inverse is `FFTW_REDFT11'.) * `FFTW_RODFT00' computes an RODFT00 transform, i.e. a DST-I. (Logical `N=2*(n+1)', inverse is `FFTW_RODFT00'.) * `FFTW_RODFT10' computes an RODFT10 transform, i.e. a DST-II. (Logical `N=2*n', inverse is `FFTW_RODFT01'.) * `FFTW_RODFT01' computes an RODFT01 transform, i.e. a DST-III. (Logical `N=2*n', inverse is `FFTW_RODFT=10'.) * `FFTW_RODFT11' computes an RODFT11 transform, i.e. a DST-IV. (Logical `N=2*n', inverse is `FFTW_RODFT11'.)  File: fftw3.info, Node: Advanced Interface, Next: Guru Interface, Prev: Basic Interface, Up: FFTW Reference Advanced Interface ================== FFTW's "advanced" interface supplements the basic interface with four new planner routines, providing a new level of flexibility: you can plan a transform of multiple arrays simultaneously, operate on non-contiguous (strided) data, and transform a subset of a larger multi-dimensional array. Other than these additional features, the planner operates in the same fashion as in the basic interface, and the resulting `fftw_plan' is used in the same way (*note Using Plans::). * Menu: * Advanced Complex DFTs:: * Advanced Real-data DFTs:: * Advanced Real-to-real Transforms::  File: fftw3.info, Node: Advanced Complex DFTs, Next: Advanced Real-data DFTs, Prev: Advanced Interface, Up: Advanced Interface Advanced Complex DFTs --------------------- fftw_plan fftw_plan_many_dft(int rank, const int *n, int howmany, fftw_complex *in, const int *inembed, int istride, int idist, fftw_complex *out, const int *onembed, int ostride, int odist, int sign, unsigned flags); This plans multidimensional complex DFTs, and is exactly the same as `fftw_plan_dft' except for the new parameters `howmany', {`i',`o'}`nembed', {`i',`o'}`stride', and {`i',`o'}`dist'. `howmany' is the number of transforms to compute, where the `k'-th transform is of the arrays starting at `in+k*idist' and `out+k*odist'. The resulting plans can often be faster than calling FFTW multiple times for the individual transforms. The basic `fftw_plan_dft' interface corresponds to `howmany=1' (in which case the `dist' parameters are ignored). The two `nembed' parameters (which should be arrays of length `rank') indicate the sizes of the input and output array dimensions, respectively, where the transform is of a subarray of size `n'. (Each dimension of `n' should be `<=' the corresponding dimension of the `nembed' arrays.) That is, the input and output arrays are stored in row-major order with size given by `nembed' (not counting the strides and howmany multiplicities). Passing `NULL' for an `nembed' parameter is equivalent to passing `n' (i.e. same physical and logical dimensions, as in the basic interface.) The `stride' parameters indicate that the `j'-th element of the input or output arrays is located at `j*istride' or `j*ostride', respectively. (For a multi-dimensional array, `j' is the ordinary row-major index.) When combined with the `k'-th transform in a `howmany' loop, from above, this means that the (`j',`k')-th element is at `j*stride+k*dist'. (The basic `fftw_plan_dft' interface corresponds to a stride of 1.) For in-place transforms, the input and output `stride' and `dist' parameters should be the same; otherwise, the planner may return `NULL'. So, for example, to transform a sequence of contiguous arrays, stored one after another, one would use a `stride' of 1 and a `dist' of N, where N is the product of the dimensions. In another example, to transform an array of contiguous "vectors" of length M, one would use a `howmany' of M, a `stride' of M, and a `dist' of 1.  File: fftw3.info, Node: Advanced Real-data DFTs, Next: Advanced Real-to-real Transforms, Prev: Advanced Complex DFTs, Up: Advanced Interface Advanced Real-data DFTs ----------------------- fftw_plan fftw_plan_many_dft_r2c(int rank, const int *n, int howmany, double *in, const int *inembed, int istride, int idist, fftw_complex *out, const int *onembed, int ostride, int odist, unsigned flags); fftw_plan fftw_plan_many_dft_c2r(int rank, const int *n, int howmany, fftw_complex *in, const int *inembed, int istride, int idist, double *out, const int *onembed, int ostride, int odist, unsigned flags); Like `fftw_plan_many_dft', these two functions add `howmany', `nembed', `stride', and `dist' parameters to the `fftw_plan_dft_r2c' and `fftw_plan_dft_c2r' functions, but otherwise behave the same as the basic interface. The interpretation of `howmany', `stride', and `dist' are the same as for `fftw_plan_many_dft', above. Note that the `stride' and `dist' for the real array are in units of `double', and for the complex array are in units of `fftw_complex'. If an `nembed' parameter is `NULL', it is interpreted as what it would be in the basic interface, as described in *Note Real-data DFT Array Format::. That is, for the complex array the size is assumed to be the same as `n', but with the last dimension cut roughly in half. For the real array, the size is assumed to be `n' if the transform is out-of-place, or `n' with the last dimension "padded" if the transform is in-place. If an `nembed' parameter is non-`NULL', it is interpreted as the physical size of the corresponding array, in row-major order, just as for `fftw_plan_many_dft'. In this case, each dimension of `nembed' should be `>=' what it would be in the basic interface (e.g. the halved or padded `n').  File: fftw3.info, Node: Advanced Real-to-real Transforms, Prev: Advanced Real-data DFTs, Up: Advanced Interface Advanced Real-to-real Transforms -------------------------------- fftw_plan fftw_plan_many_r2r(int rank, const int *n, int howmany, double *in, const int *inembed, int istride, int idist, double *out, const int *onembed, int ostride, int odist, const fftw_r2r_kind *kind, unsigned flags); Like `fftw_plan_many_dft', this functions adds `howmany', `nembed', `stride', and `dist' parameters to the `fftw_plan_r2r' function, but otherwise behave the same as the basic interface. The interpretation of those additional parameters are the same as for `fftw_plan_many_dft'. (Of course, the `stride' and `dist' parameters are now in units of `double', not `fftw_complex'.)  File: fftw3.info, Node: Guru Interface, Next: Wisdom, Prev: Advanced Interface, Up: FFTW Reference Guru Interface ============== The "guru" interface to FFTW is intended to expose as much as possible of the flexibility in the underlying FFTW architecture. It allows one to compute multi-dimensional "vectors" (loops) of multi-dimensional transforms, where each vector/transform dimension has an independent size and stride. One can also use more general complex-number formats, e.g. separate real and imaginary arrays. In addition to the more flexible planner interface, there are also guru generalizations of `fftw_execute' that allow a given plan to be executed on a different array (sufficiently similar to the original). Note that in order to use the guru execute functions, you do not need to use the guru planner interface: you can create plans with the basic or advanced interface and execute them with the guru interface. For those users who require the flexibility of the guru interface, it is important that they pay special attention to the documentation lest they shoot themselves in the foot. * Menu: * Interleaved and split arrays:: * Guru vector and transform sizes:: * Guru Complex DFTs:: * Guru Real-data DFTs:: * Guru Real-to-real Transforms:: * Guru Execution of Plans::  File: fftw3.info, Node: Interleaved and split arrays, Next: Guru vector and transform sizes, Prev: Guru Interface, Up: Guru Interface Interleaved and split arrays ---------------------------- The guru interface supports two representations of complex numbers, which we call the interleaved and the split format. The "interleaved" format is the same one used by the basic and advanced interfaces, and it is documented in *Note Complex numbers::. In the interleaved format, you provide pointers to the real part of a complex number, and the imaginary part understood to be stored in the next memory location. The "split" format allows separate pointers to the real and imaginary parts of a complex array. Technically, the interleaved format is redundant, because you can always express an interleaved array in terms of a split array with appropriate pointers and strides. On the other hand, the interleaved format is simpler to use, and it is common in practice. Hence, FFTW supports it as a special case.  File: fftw3.info, Node: Guru vector and transform sizes, Next: Guru Complex DFTs, Prev: Interleaved and split arrays, Up: Guru Interface Guru vector and transform sizes ------------------------------- The guru interface introduces one basic new data structure, `fftw_iodim', that is used to specify sizes and strides for multi-dimensional transforms and vectors: typedef struct { int n; int is; int os; } fftw_iodim; Here, `n' is the size of the dimension, and `is' and `os' are the strides of that dimension for the input and output arrays. The meaning of the stride parameter depends on the type of the array that the stride refers to. _If the array is interleaved complex, strides are expressed in units of complex numbers (`fftw_complex'). If the array is split complex or real, strides are expressed in units of real numbers (`double')._ This convention is consistent with the usual pointer arithmetic in the C language. An interleaved array is denoted by a pointer `p' to `fftw_complex', so that `p+1' points to the next complex number. Split arrays are denoted by pointers to `double', in which case pointer arithmetic operates in units of `sizeof(double)'. The guru planner interfaces all take a (`rank', `dims[rank]') pair describing the transform size, and a (`howmany_rank', `howmany_dims[rank]') pair describing the "vector" size (a multi-dimensional loop of transforms to perform), where `dims' and `howmany_dims' are arrays of `fftw_iodim'. For example, the `howmany' parameter in the advanced complex-DFT interface corresponds to `howmany_rank' = 1, `howmany_dims[0].n' = `howmany', `howmany_dims[0].is' = `idist', and `howmany_dims[0].os' = `odist'. A row-major multidimensional array with dimensions `n[rank]' (*note Row-major Format::) corresponds to `dims[i].n' = `n[i]' and the recurrence `dims[i].is' = `n[i+1] * dims[i+1].is' (similarly for `os'). The stride of the last (`i=rank-1') dimension is the overall stride of the array. e.g. to be equivalent to the advanced complex-DFT interface, you would have `dims[rank-1].is' = `istride' and `dims[rank-1].os' = `ostride'. In general, we only guarantee FFTW to return a non-`NULL' plan if the vector and transform dimensions correspond to a set of distinct indices, and for in-place transforms the input/output strides should be the same.  File: fftw3.info, Node: Guru Complex DFTs, Next: Guru Real-data DFTs, Prev: Guru vector and transform sizes, Up: Guru Interface Guru Complex DFTs ----------------- fftw_plan fftw_plan_guru_dft( int rank, const fftw_iodim *dims, int howmany_rank, const fftw_iodim *howmany_dims, fftw_complex *in, fftw_complex *out, int sign, unsigned flags); fftw_plan fftw_plan_guru_split_dft( int rank, const fftw_iodim *dims, int howmany_rank, const fftw_iodim *howmany_dims, double *ri, double *ii, double *ro, double *io, unsigned flags); These two functions plan a complex-data, multi-dimensional DFT for the interleaved and split format, respectively. Transform dimensions are given by (`rank', `dims') over a multi-dimensional vector (loop) of dimensions (`howmany_rank', `howmany_dims'). `dims' and `howmany_dims' should point to `fftw_iodim' arrays of length `rank' and `howmany_rank', respectively. `flags' is a bitwise OR (`|') of zero or more planner flags, as defined in *Note Planner Flags::. In the `fftw_plan_guru_dft' function, the pointers `in' and `out' point to the interleaved input and output arrays, respectively. The sign can be either -1 (= `FFTW_FORWARD') or +1 (= `FFTW_BACKWARD'). If the pointers are equal, the transform is in-place. In the `fftw_plan_guru_split_dft' function, `ri' and `ii' point to the real and imaginary input arrays, and `ro' and `io' point to the real and imaginary output arrays. The input and output pointers may be the same, indicating an in-place transform. For example, for `fftw_complex' pointers `in' and `out', the corresponding parameters are: ri = (double *) in; ii = (double *) in + 1; ro = (double *) out; io = (double *) out + 1; Because `fftw_plan_guru_split_dft' accepts split arrays, strides are expressed in units of `double'. For a contiguous `fftw_complex' array, the overall stride of the transform should be 2, the distance between consecutive real parts or between consecutive imaginary parts; see *Note Guru vector and transform sizes::. Note that the dimension strides are applied equally to the real and imaginary parts; real and imaginary arrays with different strides are not supported. There is no `sign' parameter in `fftw_plan_guru_split_dft'. This function always plans for an `FFTW_FORWARD' transform. To plan for an `FFTW_BACKWARD' transform, you can exploit the identity that the backwards DFT is equal to the forwards DFT with the real and imaginary parts swapped. For example, in the case of the `fftw_complex' arrays above, the `FFTW_BACKWARD' transform is computed by the parameters: ri = (double *) in + 1; ii = (double *) in; ro = (double *) out + 1; io = (double *) out;