fftw3: Allocating aligned memory in Fortran
7.5 Allocating aligned memory in Fortran
========================================
In order to obtain maximum performance in FFTW, you should store your
data in arrays that have been specially aligned in memory (SIMD
alignment and fftw_malloc). Enforcing alignment also permits you to
Functions::) to apply a given plan to more than one pair of in/out
arrays. Unfortunately, standard Fortran arrays do _not_ provide any
alignment guarantees. The _only_ way to allocate aligned memory in
standard Fortran is to allocate it with an external C function, like the
'fftw_alloc_real' and 'fftw_alloc_complex' functions. Fortunately,
Fortran 2003 provides a simple way to associate such allocated memory
with a standard Fortran array pointer that you can then use normally.
We therefore recommend allocating all your input/output arrays using
the following technique:
1. Declare a 'pointer', 'arr', to your array of the desired type and
dimensions. For example, 'real(C_DOUBLE), pointer :: a(:,:)' for a
2d real array, or 'complex(C_DOUBLE_COMPLEX), pointer :: a(:,:,:)'
for a 3d complex array.
2. The number of elements to allocate must be an 'integer(C_SIZE_T)'.
You can either declare a variable of this type, e.g.
'integer(C_SIZE_T) :: sz', to store the number of elements to
allocate, or you can use the 'int(..., C_SIZE_T)' intrinsic
function. e.g. set 'sz = L * M * N' or use 'int(L * M * N,
C_SIZE_T)' for an L x M x N array.
3. Declare a 'type(C_PTR) :: p' to hold the return value from FFTW's
allocation routine. Set 'p = fftw_alloc_real(sz)' for a real
array, or 'p = fftw_alloc_complex(sz)' for a complex array.
4. Associate your pointer 'arr' with the allocated memory 'p' using
the standard 'c_f_pointer' subroutine: 'call c_f_pointer(p, arr,
[...dimensions...])', where '[...dimensions...])' are an array of
the dimensions of the array (in the usual Fortran order). e.g.
'call c_f_pointer(p, arr, [L,M,N])' for an L x M x N array.
(Alternatively, you can omit the dimensions argument if you
specified the shape explicitly when declaring 'arr'.) You can now
use 'arr' as a usual multidimensional array.
5. When you are done using the array, deallocate the memory by 'call
fftw_free(p)' on 'p'.
For example, here is how we would allocate an L x M 2d real array:
real(C_DOUBLE), pointer :: arr(:,:)
type(C_PTR) :: p
p = fftw_alloc_real(int(L * M, C_SIZE_T))
call c_f_pointer(p, arr, [L,M])
_...use arr and arr(i,j) as usual..._
call fftw_free(p)
and here is an L x M x N 3d complex array:
complex(C_DOUBLE_COMPLEX), pointer :: arr(:,:,:)
type(C_PTR) :: p
p = fftw_alloc_complex(int(L * M * N, C_SIZE_T))
call c_f_pointer(p, arr, [L,M,N])
_...use arr and arr(i,j,k) as usual..._
call fftw_free(p)
See Reversing array dimensions for an example allocating a
single array and associating both real and complex array pointers with
it, for in-place real-to-complex transforms.