Rev 3 | Details | Compare with Previous | Last modification | View Log | RSS feed
Rev | Author | Line No. | Line |
---|---|---|---|
2 | pj | 1 | <!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 3.2//EN"> |
2 | <HTML> |
||
3 | <HEAD> |
||
4 | <!-- This HTML file has been created by texi2html 1.52 |
||
5 | from fftw.texi on 18 May 1999 --> |
||
6 | |||
7 | <TITLE>FFTW - Calling FFTW from Fortran</TITLE> |
||
8 | </HEAD> |
||
9 | <BODY TEXT="#000000" BGCOLOR="#FFFFFF"> |
||
10 | Go to the <A HREF="fftw_1.html">first</A>, <A HREF="fftw_4.html">previous</A>, <A HREF="fftw_6.html">next</A>, <A HREF="fftw_10.html">last</A> section, <A HREF="fftw_toc.html">table of contents</A>. |
||
11 | <P><HR><P> |
||
12 | |||
13 | |||
14 | <H1><A NAME="SEC62">Calling FFTW from Fortran</A></H1> |
||
15 | |||
16 | <P> |
||
17 | <A NAME="IDX290"></A> |
||
18 | The standard FFTW libraries include special wrapper functions that allow |
||
19 | Fortran programs to call FFTW subroutines. This chapter describes how |
||
20 | those functions may be employed to use FFTW from Fortran. We assume |
||
21 | here that the reader is already familiar with the usage of FFTW in C, as |
||
22 | described elsewhere in this manual. |
||
23 | |||
24 | |||
25 | <P> |
||
26 | In general, it is not possible to call C functions directly from |
||
27 | Fortran, due to Fortran's inability to pass arguments by value and also |
||
28 | because Fortran compilers typically expect identifiers to be mangled |
||
29 | <A NAME="IDX291"></A> |
||
30 | somehow for linking. However, if C functions are written in a special |
||
31 | way, they <EM>are</EM> callable from Fortran, and we have employed this |
||
32 | technique to create Fortran-callable "wrapper" functions around the |
||
33 | main FFTW routines. These wrapper functions are included in the FFTW |
||
34 | libraries by default, unless a Fortran compiler isn't found on your |
||
35 | system or <CODE>--disable-fortran</CODE> is included in the <CODE>configure</CODE> |
||
36 | flags. |
||
37 | |||
38 | |||
39 | <P> |
||
40 | As a result, calling FFTW from Fortran requires little more than |
||
41 | appending <SAMP>`<CODE>_f77</CODE>'</SAMP> to the function names and then linking |
||
42 | normally with the FFTW libraries. There are a few wrinkles, however, as |
||
43 | we shall discuss below. |
||
44 | |||
45 | |||
46 | |||
47 | |||
48 | <H2><A NAME="SEC63">Wrapper Routines</A></H2> |
||
49 | |||
50 | <P> |
||
51 | All of the uniprocessor and multi-threaded transform routines have |
||
52 | Fortran-callable wrappers, except for the wisdom import/export functions |
||
53 | (since it is not possible to exchange string and file arguments portably |
||
54 | with Fortran). The name of the wrapper routine is the same as that of |
||
55 | the corresponding C routine, but with <CODE>fftw/fftwnd/rfftw/rfftwnd</CODE> |
||
56 | replaced by <CODE>fftw_f77/fftwnd_f77/rfftw_f77/rfftwnd_f77</CODE>. For |
||
57 | example, in Fortran, instead of calling <CODE>fftw_one</CODE> you would call |
||
58 | <CODE>fftw_f77_one</CODE>.<A NAME="DOCF8" HREF="fftw_foot.html#FOOT8">(8)</A> |
||
59 | <A NAME="IDX292"></A> |
||
60 | For the most part, all of the arguments to the functions are the same, |
||
61 | with the following exceptions: |
||
62 | |||
63 | |||
64 | |||
65 | <UL> |
||
66 | |||
67 | <LI> |
||
68 | |||
69 | Any function that returns a value (e.g. <CODE>fftw_create_plan</CODE>) is |
||
70 | converted into a subroutine. The return value is converted into an |
||
71 | additional (first) parameter of the wrapper subroutine. (The reason for |
||
72 | this is that some Fortran implementations seem to have trouble with C |
||
73 | function return values.) |
||
74 | |||
75 | <LI> |
||
76 | |||
77 | <A NAME="IDX293"></A> |
||
78 | When performing one-dimensional <CODE>FFTW_IN_PLACE</CODE> transforms, you |
||
79 | don't have the option of passing <CODE>NULL</CODE> for the <CODE>out</CODE> argument |
||
80 | (since there is no way to pass <CODE>NULL</CODE> from Fortran). Therefore, |
||
81 | when performing such transforms, you <EM>must</EM> allocate and pass a |
||
82 | contiguous scratch array of the same size as the transform. Note that |
||
83 | for in-place multi-dimensional (<CODE>(r)fftwnd</CODE>) transforms, the |
||
84 | <CODE>out</CODE> argument is ignored, so you can pass anything for that |
||
85 | parameter. |
||
86 | |||
87 | <LI> |
||
88 | |||
89 | <A NAME="IDX294"></A> |
||
90 | The wrapper routines expect multi-dimensional arrays to be in |
||
91 | column-major order, which is the ordinary format of Fortran arrays. |
||
92 | They do this transparently and costlessly simply by reversing the order |
||
93 | of the dimensions passed to FFTW, but this has one important consequence |
||
94 | for multi-dimensional real-complex transforms, discussed below. |
||
95 | |||
96 | <LI> |
||
97 | |||
98 | <CODE>plan</CODE> variables (what would be of type <CODE>fftw_plan</CODE>, |
||
99 | <CODE>rfftwnd_plan</CODE>, etcetera, in C), must be declared as a type that is |
||
100 | the same size as a pointer (address) on your machine. (Fortran has no |
||
101 | generic pointer type.) The Fortran <CODE>integer</CODE> type is usually the |
||
102 | same size as a pointer, but you need to be wary (especially on 64-bit |
||
103 | machines). (You could also use <CODE>integer*4</CODE> on a 32-bit machine and |
||
104 | <CODE>integer*8</CODE> on a 64-bit machine.) Ugh. (<CODE>g77</CODE> has a special |
||
105 | type, <CODE>integer(kind=7)</CODE>, that is defined to be the same size as a |
||
106 | pointer.) |
||
107 | |||
108 | </UL> |
||
109 | |||
110 | <P> |
||
111 | <A NAME="IDX295"></A> |
||
112 | In general, you should take care to use Fortran data types that |
||
113 | correspond to (i.e. are the same size as) the C types used by FFTW. If |
||
114 | your C and Fortran compilers are made by the same vendor, the |
||
115 | correspondence is usually straightforward (i.e. <CODE>integer</CODE> |
||
116 | corresponds to <CODE>int</CODE>, <CODE>real</CODE> corresponds to <CODE>float</CODE>, |
||
117 | etcetera). Such simple correspondences are assumed in the examples |
||
118 | below. The examples also assume that FFTW was compiled in |
||
119 | double precision (the default). |
||
120 | |||
121 | |||
122 | |||
123 | |||
124 | <H2><A NAME="SEC64">FFTW Constants in Fortran</A></H2> |
||
125 | |||
126 | <P> |
||
127 | When creating plans in FFTW, a number of constants are used to specify |
||
128 | options, such as <CODE>FFTW_FORWARD</CODE> or <CODE>FFTW_USE_WISDOM</CODE>. The |
||
129 | same constants must be used with the wrapper routines, but of course the |
||
130 | C header files where the constants are defined can't be incorporated |
||
131 | directly into Fortran code. |
||
132 | |||
133 | |||
134 | <P> |
||
135 | Instead, we have placed Fortran equivalents of the FFTW constant |
||
136 | definitions in the file <CODE>fortran/fftw_f77.i</CODE> of the FFTW package. |
||
137 | If your Fortran compiler supports a preprocessor, you can use that to |
||
138 | incorporate this file into your code whenever you need to call FFTW. |
||
139 | Otherwise, you will have to paste the constant definitions in directly. |
||
140 | They are: |
||
141 | |||
142 | |||
143 | |||
144 | <PRE> |
||
145 | integer FFTW_FORWARD,FFTW_BACKWARD |
||
146 | parameter (FFTW_FORWARD=-1,FFTW_BACKWARD=1) |
||
147 | |||
148 | integer FFTW_REAL_TO_COMPLEX,FFTW_COMPLEX_TO_REAL |
||
149 | parameter (FFTW_REAL_TO_COMPLEX=-1,FFTW_COMPLEX_TO_REAL=1) |
||
150 | |||
151 | integer FFTW_ESTIMATE,FFTW_MEASURE |
||
152 | parameter (FFTW_ESTIMATE=0,FFTW_MEASURE=1) |
||
153 | |||
154 | integer FFTW_OUT_OF_PLACE,FFTW_IN_PLACE,FFTW_USE_WISDOM |
||
155 | parameter (FFTW_OUT_OF_PLACE=0) |
||
156 | parameter (FFTW_IN_PLACE=8,FFTW_USE_WISDOM=16) |
||
157 | |||
158 | integer FFTW_THREADSAFE |
||
159 | parameter (FFTW_THREADSAFE=128) |
||
160 | </PRE> |
||
161 | |||
162 | <P> |
||
163 | <A NAME="IDX296"></A> |
||
164 | In C, you combine different flags (like <CODE>FFTW_USE_WISDOM</CODE> and |
||
165 | <CODE>FFTW_MEASURE</CODE>) using the <SAMP>`<CODE>|</CODE>'</SAMP> operator; in Fortran you |
||
166 | should just use <SAMP>`<CODE>+</CODE>'</SAMP>. |
||
167 | |||
168 | |||
169 | |||
170 | |||
171 | <H2><A NAME="SEC65">Fortran Examples</A></H2> |
||
172 | |||
173 | <P> |
||
174 | In C you might have something like the following to transform a |
||
175 | one-dimensional complex array: |
||
176 | |||
177 | |||
178 | |||
179 | <PRE> |
||
180 | fftw_complex in[N], *out[N]; |
||
181 | fftw_plan plan; |
||
182 | |||
183 | plan = fftw_create_plan(N,FFTW_FORWARD,FFTW_ESTIMATE); |
||
184 | fftw_one(plan,in,out); |
||
185 | fftw_destroy_plan(plan); |
||
186 | </PRE> |
||
187 | |||
188 | <P> |
||
189 | In Fortran, you use the following to accomplish the same thing: |
||
190 | |||
191 | |||
192 | |||
193 | <PRE> |
||
194 | double complex in, out |
||
195 | dimension in(N), out(N) |
||
196 | integer plan |
||
197 | |||
198 | call fftw_f77_create_plan(plan,N,FFTW_FORWARD,FFTW_ESTIMATE) |
||
199 | call fftw_f77_one(plan,in,out) |
||
200 | call fftw_f77_destroy_plan(plan) |
||
201 | </PRE> |
||
202 | |||
203 | <P> |
||
204 | <A NAME="IDX297"></A> |
||
205 | <A NAME="IDX298"></A> |
||
206 | <A NAME="IDX299"></A> |
||
207 | |||
208 | |||
209 | <P> |
||
210 | Notice how all routines are called as Fortran subroutines, and the plan |
||
211 | is returned via the first argument to <CODE>fftw_f77_create_plan</CODE>. To |
||
212 | do the same thing, but using 8 threads in parallel |
||
213 | (see Section <A HREF="fftw_4.html#SEC48">Multi-threaded FFTW</A>), you would simply replace the call to |
||
214 | <CODE>fftw_f77_one</CODE> with: |
||
215 | |||
216 | |||
217 | |||
218 | <PRE> |
||
219 | call fftw_f77_threads_one(8,plan,in,out) |
||
220 | </PRE> |
||
221 | |||
222 | <P> |
||
223 | <A NAME="IDX300"></A> |
||
224 | |||
225 | |||
226 | <P> |
||
227 | To transform a three-dimensional array in-place with C, you might do: |
||
228 | |||
229 | |||
230 | |||
231 | <PRE> |
||
232 | fftw_complex arr[L][M][N]; |
||
233 | fftwnd_plan plan; |
||
234 | int n[3] = {L,M,N}; |
||
235 | |||
236 | plan = fftwnd_create_plan(3,n,FFTW_FORWARD, |
||
237 | FFTW_ESTIMATE | FFTW_IN_PLACE); |
||
238 | fftwnd_one(plan, arr, 0); |
||
239 | fftwnd_destroy_plan(plan); |
||
240 | </PRE> |
||
241 | |||
242 | <P> |
||
243 | In Fortran, you would use this instead: |
||
244 | |||
245 | |||
246 | |||
247 | <PRE> |
||
248 | double complex arr |
||
249 | dimension arr(L,M,N) |
||
250 | integer n |
||
251 | dimension n(3) |
||
252 | integer plan |
||
253 | |||
254 | n(1) = L |
||
255 | n(2) = M |
||
256 | n(3) = N |
||
257 | call fftwnd_f77_create_plan(plan,3,n,FFTW_FORWARD, |
||
258 | + FFTW_ESTIMATE + FFTW_IN_PLACE) |
||
259 | call fftwnd_f77_one(plan, arr, 0) |
||
260 | call fftwnd_f77_destroy_plan(plan) |
||
261 | </PRE> |
||
262 | |||
263 | <P> |
||
264 | <A NAME="IDX301"></A> |
||
265 | <A NAME="IDX302"></A> |
||
266 | <A NAME="IDX303"></A> |
||
267 | |||
268 | |||
269 | <P> |
||
270 | Instead of calling <CODE>fftwnd_f77_create_plan(plan,3,n,...)</CODE>, we could |
||
271 | also have called <CODE>fftw3d_f77_create_plan(plan,L,M,N,...)</CODE>. |
||
272 | <A NAME="IDX304"></A> |
||
273 | |||
274 | |||
275 | <P> |
||
276 | Note that we pass the array dimensions in the "natural" order; also |
||
277 | note that the last argument to <CODE>fftwnd_f77</CODE> is ignored since the |
||
278 | transform is <CODE>FFTW_IN_PLACE</CODE>. |
||
279 | |||
280 | |||
281 | <P> |
||
282 | To transform a one-dimensional real array in Fortran, you might do: |
||
283 | |||
284 | |||
285 | |||
286 | <PRE> |
||
287 | double precision in, out |
||
288 | dimension in(N), out(N) |
||
289 | integer plan |
||
290 | |||
291 | call rfftw_f77_create_plan(plan,N,FFTW_REAL_TO_COMPLEX, |
||
292 | + FFTW_ESTIMATE) |
||
293 | call rfftw_f77_one(plan,in,out) |
||
294 | call rfftw_f77_destroy_plan(plan) |
||
295 | </PRE> |
||
296 | |||
297 | <P> |
||
298 | <A NAME="IDX305"></A> |
||
299 | <A NAME="IDX306"></A> |
||
300 | <A NAME="IDX307"></A> |
||
301 | |||
302 | |||
303 | <P> |
||
304 | To transform a two-dimensional real array, out of place, you might use |
||
305 | the following: |
||
306 | |||
307 | |||
308 | |||
309 | <PRE> |
||
310 | double precision in |
||
311 | double complex out |
||
312 | dimension in(M,N), out(M/2 + 1, N) |
||
313 | integer plan |
||
314 | |||
315 | call rfftw2d_f77_create_plan(plan,M,N,FFTW_REAL_TO_COMPLEX, |
||
316 | + FFTW_ESTIMATE) |
||
317 | call rfftwnd_f77_one_real_to_complex(plan, in, out) |
||
318 | call rfftwnd_f77_destroy_plan(plan) |
||
319 | </PRE> |
||
320 | |||
321 | <P> |
||
322 | <A NAME="IDX308"></A> |
||
323 | <A NAME="IDX309"></A> |
||
324 | <A NAME="IDX310"></A> |
||
325 | |||
326 | |||
327 | <P> |
||
328 | <B>Important:</B> Notice that it is the <EM>first</EM> dimension of the |
||
329 | complex output array that is cut in half in Fortran, rather than the |
||
330 | last dimension as in C. This is a consequence of the wrapper routines |
||
331 | reversing the order of the array dimensions passed to FFTW so that the |
||
332 | Fortran program can use its ordinary column-major order. |
||
333 | <A NAME="IDX311"></A> |
||
334 | <A NAME="IDX312"></A> |
||
335 | |||
336 | |||
337 | <P><HR><P> |
||
338 | Go to the <A HREF="fftw_1.html">first</A>, <A HREF="fftw_4.html">previous</A>, <A HREF="fftw_6.html">next</A>, <A HREF="fftw_10.html">last</A> section, <A HREF="fftw_toc.html">table of contents</A>. |
||
339 | </BODY> |
||
340 | </HTML> |