	   "Wrapper" Routines for Calling FFTW from Fortran
		   http://theory.lcs.mit.edu/~fftw

FFTW cannot be called directly from standard Fortran because of
difficulties in parameter-passing and other matters.  This directory
contains specially-written, Fortran-callable, C "wrapper" routines for
the FFTW functions.

Before reading this document, you should already be familiar with the
ordinary usage of FFTW--see the FFTW manual in the doc/ directory.

Not all of the routines in FFTW are exposed to Fortran via wrapper
functions.  The unexposed routines are left as an exercise for the
reader.

*****

The wrapper routines: 

fftw_f77.c contains the wrapper routines for complex-complex
transforms, and rfftw_f77.c contains wrappers for real-complex
transforms.  They address three main difficulties in calling FFTW from
Fortran:

	* Fortran passes all parameters by reference

	* Fortran is case-insensitive, and the linker therefore
	  expects routine names to be "munged" in a certain way
	  (see the FORTRANIZE macro, described below).

	* Some Fortran implementations seem to have difficulty in
	  calling C functions that return a value.  All return values
	  were therefore changed into parameters.

The following functions are "exposed" to Fortran via wrapper routines:

    fftw_create_plan      -> fftw_f77_create_plan
    fftw_destroy_plan     -> fftw_f77_destroy_plan
    fftw                  -> fftw_f77
    fftw_one              -> fftw_f77_one
    fftwnd_create_plan    -> fftwnd_f77_create_plan
    fftw2d_create_plan    -> fftw2d_f77_create_plan
    fftw3d_create_plan    -> fftw3d_f77_create_plan
    fftwnd_destroy_plan   -> fftwnd_f77_destroy_plan
    fftwnd                -> fftwnd_f77
    fftwnd_one            -> fftwnd_f77_one

   rfftw_create_plan      -> rfftw_f77_create_plan
   rfftw_destroy_plan     -> rfftw_f77_destroy_plan
   rfftw                  -> rfftw_f77
   rfftw_one              -> rfftw_f77_one
   rfftwnd_create_plan    -> rfftwnd_f77_create_plan
   rfftw2d_create_plan    -> rfftw2d_f77_create_plan
   rfftw3d_create_plan    -> rfftw3d_f77_create_plan
   rfftwnd_destroy_plan   -> rfftwnd_f77_destroy_plan
   rfftwnd_real_to_complex   -> rfftwnd_f77_real_to_complex
   rfftwnd_complex_to_real   -> rfftwnd_f77_complex_to_real
   rfftwnd_one_real_to_complex   -> rfftwnd_f77_one_real_to_complex
   rfftwnd_one_complex_to_real   -> rfftwnd_f77_one_complex_to_real

We use the FORTRANIZE macro (defined in fortranize.h) to convert
names into the format required by the Fortran linker.  By default,
this macro uses all lower-case names with an underscore appended.  You
may need to change this, depending upon your compiler; see how we
handled the AIX, Solaris, and Cray compilers for examples.  (If you
don't know what your linker expects, just try compiling and linking
fftw_f77.c and f77_test.F, and see if you get link errors.)  (g77 users
should compile with -DUSING_G77 in order to use the proper FORTRANIZE.)

*****

Calling the FFTW wrapper routines:

The FFTW wrapper routines may be called from Fortran in much the same
way as the corresponding routines are called in C (see the FFTW
documentation).  There are a few differences, however.

	* The routine names are changed: "fftw" and "fftwnd" become
	  "fftw_f77" and "fftwnd_f77", respectively (see above list).
	  Similarly for "rfftw" and "rfftwnd".

        * You should paste in the constants from fftw_f77.i,
	  or #include this file if you are using the C preprocessor.

	* Return values of functions become, instead, the first
	  parameter of a subroutine.  All other parameters are unchanged.

	* Plans (what would be fftw_plan, etcetera, in C) 
	  should be declared as a type that is the same size as
	  a pointer on your machine (usually "integer").

	* Fortran programs should be linked with fftw_f77.o (the
          compiled version of fftw_f77.c) in addition to -lm -lfftw.
	  If you use rfftw, you should also link with rfftw_f77.o and
	  -l rfftw.  Enterprising readers may wish to instead incorporate
	  fftw_f77.o and rfftw_f77.o into libfftw.a and librfftw.a,
	  respectively.

Remember that the default compilation of FFTW uses double-precision
floating point values, so you must do the same in your Fortran code.
(The floating-point precision is easily changed; see the FFTW
documentation.)

NOTE: Fortran multi-dimensional arrays are stored in column-major
(Fortran) order, while, ordinarily, FFTW expects them in row-major (C)
order.  However, the wrapper functions translate the array description
from column-major to row-major for you, so you can pass the array
dimensions in the natural order for Fortran.  (All the wrapper
functions really do is reverse the order of the array dimensions when
creating the plan.  See the FFTW documentation for more discussion of
array formats.)

NOTE: In C, the final dimensions of the real and complex input and
output arrays of rfftwnd are related by n -> n/2+1.  In Fortran, this
relationship applies to the leading dimensions of the arrays.  (See
below for an example.)

For example, in C you might have something like:

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

	plan = fftw_create_plan(N,FFTW_FORWARD,FFTW_ESTIMATE);
	fftw_one(plan,in,out);
	fftw_destroy_plan(plan);

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

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

	call fftw_f77_create_plan(plan,N,FFTW_FORWARD,FFTW_ESTIMATE)
	call fftw_f77_one(plan,in,out)
	call fftw_f77_destroy_plan(plan)

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

	fftw_complex arr[L][M][N];
	fftwnd_plan plan;
	int n[3] = {L,M,N};

	plan = fftwnd_create_plan(3,n,FFTW_FORWARD,
                                  FFTW_ESTIMATE | FFTW_IN_PLACE);
	fftwnd_one(plan, arr, 0);
	fftwnd_destroy_plan(plan);

In Fortran, you would use this instead:

	double complex arr
	dimension arr(L,M,N)
	integer n
	dimension n(3)
	integer plan

	n(1) = L
	n(2) = M
	n(3) = N
	call fftwnd_f77_create_plan(plan,3,n,FFTW_FORWARD,
       +                            FFTW_ESTIMATE + FFTW_IN_PLACE)
	call fftwnd_f77_one(plan, arr, 0)
	call fftwnd_f77_destroy_plan(plan)

Instead of calling fftwnd_f77_create_plan(plan,3,n,...), we could
have called fftw3d_f77_create_plan(plan,L,M,N,...).

Note that we pass the array dimensions in the "natural" order; also
note that the last argument to fftwnd_f77 is ignored since the
transform is FFTW_IN_PLACE.

* N.B. The FFTW_IN_PLACE option works differently in the 1D fftw call,
  and could cause problems there if you attempt it from Fortran since
  there's no way to pass a NULL pointer for the "out" argument.  See
  the reference section of the FFTW documentation for more detail.

An example program that calls the wrapper routines, f77_test.F, is
included in this package.

Finally, here is an example where we transform a three-dimensional
real array using rfftwnd:
	
	double precision input
	dimension input(L,M,N)
	double complex output
	dimension output(L/2+1,M,N)
	integer plan

	call rfftw3d_f77_create_plan(plan,L,M,N,FFTW_REAL_TO_COMPLEX,
       +                             FFTW_ESTIMATE)
	call rfftwnd_f77_one_real_to_complex(plan,input,output)
	call rfftwnd_f77_destroy_plan(plan)
