Next: , Previous: , Up: Calling FFTW from Legacy Fortran   [Contents][Index]


8.4 Fortran Examples

In C, you might have something like the following to transform a one-dimensional complex array:

        fftw_complex in[N], out[N];
        fftw_plan plan;

        plan = fftw_plan_dft_1d(N,in,out,FFTW_FORWARD,FFTW_ESTIMATE);
        fftw_execute(plan);
        fftw_destroy_plan(plan);

In Fortran, you would use the following to accomplish the same thing:

        double complex in, out
        dimension in(N), out(N)
        integer*8 plan

        call dfftw_plan_dft_1d(plan,N,in,out,FFTW_FORWARD,FFTW_ESTIMATE)
        call dfftw_execute_dft(plan, in, out)
        call dfftw_destroy_plan(plan)

Notice how all routines are called as Fortran subroutines, and the plan is returned via the first argument to dfftw_plan_dft_1d. Notice also that we changed fftw_execute to dfftw_execute_dft (see FFTW Execution in Fortran). To do the same thing, but using 8 threads in parallel (see Multi-threaded FFTW), you would simply prefix these calls with:

        integer iret
        call dfftw_init_threads(iret)
        call dfftw_plan_with_nthreads(8)

(You might want to check the value of iret: if it is zero, it indicates an unlikely error during thread initialization.)

To check the number of threads currently being used by the planner, you can do the following:

        integer iret
        call dfftw_planner_nthreads(iret)

To transform a three-dimensional array in-place with C, you might do:

        fftw_complex arr[L][M][N];
        fftw_plan plan;

        plan = fftw_plan_dft_3d(L,M,N, arr,arr,
                                FFTW_FORWARD, FFTW_ESTIMATE);
        fftw_execute(plan);
        fftw_destroy_plan(plan);

In Fortran, you would use this instead:

        double complex arr
        dimension arr(L,M,N)
        integer*8 plan

        call dfftw_plan_dft_3d(plan, L,M,N, arr,arr,
       &                       FFTW_FORWARD, FFTW_ESTIMATE)
        call dfftw_execute_dft(plan, arr, arr)
        call dfftw_destroy_plan(plan)

Note that we pass the array dimensions in the “natural” order in both C and Fortran.

To transform a one-dimensional real array in Fortran, you might do:

        double precision in
        dimension in(N)
        double complex out
        dimension out(N/2 + 1)
        integer*8 plan

        call dfftw_plan_dft_r2c_1d(plan,N,in,out,FFTW_ESTIMATE)
        call dfftw_execute_dft_r2c(plan, in, out)
        call dfftw_destroy_plan(plan)

To transform a two-dimensional real array, out of place, you might use the following:

        double precision in
        dimension in(M,N)
        double complex out
        dimension out(M/2 + 1, N)
        integer*8 plan

        call dfftw_plan_dft_r2c_2d(plan,M,N,in,out,FFTW_ESTIMATE)
        call dfftw_execute_dft_r2c(plan, in, out)
        call dfftw_destroy_plan(plan)

Important: Notice that it is the first dimension of the complex output array that is cut in half in Fortran, rather than the last dimension as in C. This is a consequence of the interface routines reversing the order of the array dimensions passed to FFTW so that the Fortran program can use its ordinary column-major order.


Next: Wisdom of Fortran?, Previous: FFTW Execution in Fortran, Up: Calling FFTW from Legacy Fortran   [Contents][Index]