- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

I would like to perform inplace real-to-complex and inverse transforms! The function has a signature

fftw_plan fftw_mpi_plan_dft_r2c_2d(ptrdiff_t n0, ptrdiff_t n1, double *in, fftw_complex *out, MPI_Comm comm, unsigned flags);

For an inplace complex-to-complex transform *in,*out are pointers of same type. However, here the pointer type differs. In C you can possibly typecast

a double pointer to a fftw_complex. I am confused how to do it in fortran? Any suggestions?

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

Navdeep R. wrote:

..

program fftw3_demo .. !! compiles and runs correctly till here. pfor = fftw_mpi_plan_dft_r2c_3d(N(3),N(2),N(1),u_r,u_c,MPI_COMM_WORLD,FFTW_FORWARD,FFTW_ESTIMATE) ..

According to this **site**, the interface for fftw_mpi_plan_dft_r2c_3d is

type(C_PTR) function fftw_mpi_plan_dft_r2c_3d(n0,n1,n2,in,out,comm,flags) bind(C, name='fftw_mpi_plan_dft_r2c_3d_f03') import integer(C_INTPTR_T), value :: n0 integer(C_INTPTR_T), value :: n1 integer(C_INTPTR_T), value :: n2 real(C_DOUBLE), dimension(*), intent(out) :: in complex(C_DOUBLE_COMPLEX), dimension(*), intent(out) :: out integer(C_MPI_FINT), value :: comm integer(C_INT), value :: flags end function fftw_mpi_plan_dft_r2c_3d

which has 7 dummy arguments, however your code has 8 which would explain, "`Error: More actual than formal arguments `

`in`

`procedure call at (1)`

". from gfortran.

Is one of FFTW_FORWARD or FFTW_ESTIMATE superfluous?

Link Copied

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

From what I read, fftw_complex is the same as COMPLEX(C_DOUBLE) (or COMPLEX(8)) in Intel Fortran. It is not a double pointer. http://www.fftw.org/doc/Complex-numbers.html

Note that COMPLEX(8) is not the same as COMPLEX*8!

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

Steve Lionel (Ret.) wrote:

From what I read, fftw_complex is the same as COMPLEX(C_DOUBLE) (or COMPLEX(8)) in Intel Fortran. It is not a double pointer. http://www.fftw.org/doc/Complex-numbers.html

Note that COMPLEX(8) is not the same as COMPLEX*8!

This I understand but an inplace transform means the input (double *in) and output (fftw_complex *out) are the same. In C you can do (fftw_complex *) double_pointer to achieve it. How do I do it in fortran?

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

Why not just pass the same variable as both arguments? Unless my understanding of C is much worse than I think it is (always a possibility), you're not being asked to pass a pointer to a pointer, unless I am seriously misreading http://www.fftw.org/doc/Real_002ddata-DFTs.html

The problem you'll run into in Fortran is type checking if you have declared an explicit interface for the FFTW routine that has the arguments of different types. The easy solution is to give both arguments the same type in your interface. A more complicated solution would be to declare a pointer to an array of double and use C_LOC and C_F_POINTER to "cast" the complex to double (or vice-versa.) The easiest solution is to not use an explicit interface at all, though I find that distasteful.

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

fftw_complex is defined as

typedef double fftw_complex[2]

so It is not a pointer to a pointer as you mentioned.

Yes I have an explicit interface and I plan to keep it. As I understand it, I have to define two pointers

real(C_DOUBLE), pointer, dimension(:,:,:) :: u_r complex(C_DOUBLE_COMPLEX), pointer, dimension(:,:,:) :: u_c

one for real and other for complex data. I allocate my data using fftw_alloc_real i.e.

type(C_PTR) :: ptr_u ptr_u = fftw_alloc_real(2*alloc_size)

and then typecast as

call c_f_pointer(ptr_u,u_r,[2*(N(1)/2+1),N(2),local_n3]) call c_f_pointer(ptr_u,u_c,[(N(1)/2+1),N(2),local_n3])

the shapes come out correctly, which I checked but when I go on creating a plan i.e

pfor = fftw_mpi_plan_dft_r2c_3d(N(3),N(2),N(1),u_r,u_c,MPI_COMM_WORLD,FFTW_FORWARD,FFTW_ESTIMATE)

It spits out an error i.e.

Error: More actual than formal arguments in procedure call at (1)

Any ideas on this?

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

That's not an ifort error (gfortran?), and I don't know what your interface looks like. It does seem a simple error to resolve if you compare the call with the interface.

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

Navdeep R. wrote:

.. Any ideas on this?

Assuming the function prototype is

fftw_plan fftw_mpi_plan_dft_r2c_3d(ptrdiff_t n0, ptrdiff_t n1, ptrdiff_t n2, double *in, fftw_complex *out, MPI_Comm comm, unsigned flags);

and with typedef fftw_complex as you indicate, your code is inconsistent when it comes to 5th parameter in this function i.e., with respect to fftw_complex *out:

`complex`

`(C_DOUBLE_COMPLEX), pointer is`

**not**the same as fftw_complex * because- fftw_complex * essentially reduces to double * given the above typedef

You can just have both your u_r and u_c objects in Fortran declared as `real`

`(C_DOUBLE), `

`pointer`

`, `

`dimension`

`(:,:,:) and i think your code will work ok.`

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

FortranFan wrote:

Quote:You can just have both your u_r and u_c objects in Fortran declared as real(C_DOUBLE), pointer, dimension(:,:,:) and i think your code will work ok.

That's the first thing I tried, but it gives a type mismatch error.

Error: Type mismatch in argument 'out' at (1); passed REAL(8) to COMPLEX(8)

Hence the question!

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

I'm just curious why you are asking in this forum when you're not using ifort, not that we're unwilling to help.

Please show your explicit interface.

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

I have a student license for ifort and related tools, which is installed on my main machine.

For the time being I am on my laptop which only has gfortran.

I inherited a pseudo-spectral CFD code which uses FFTW2 interface, and I would like to switch to FFTW3

program fftw3_demo use mpi use, intrinsic :: iso_c_binding implicit none include "fftw3-mpi.f03" integer(C_INTPTR_T), dimension(3) :: N = [32,32,32] real(C_DOUBLE), pointer, dimension(:,:,:) :: u_r complex(C_DOUBLE_COMPLEX), pointer, dimension(:,:,:) :: u_c integer(C_INTPTR_T) :: alloc_size,local_n3,local_n3_start type(C_PTR) :: ptr_u,pfor,pinv integer :: mpi_error_flag,mpi_rank,mpi_procs real(kind=8) :: u_sum_old call mpi_init(mpi_error_flag) call mpi_comm_rank(MPI_COMM_WORLD, mpi_rank,mpi_error_flag) call mpi_comm_size(MPI_COMM_WORLD, mpi_procs, mpi_error_flag) call fftw_mpi_init() alloc_size = fftw_mpi_local_size_3d(N(3),N(2),N(1)/2+1,MPI_COMM_WORLD,local_n3,local_n3_start) ptr_u = fftw_alloc_real(2*alloc_size) call c_f_pointer(ptr_u,u_r,[2*(N(1)/2+1),N(2),local_n3]) call c_f_pointer(ptr_u,u_c,[(N(1)/2+1),N(2),local_n3]) write(*,*) shape(u_r), shape(u_c) !! compiles and runs correctly till here. pfor = fftw_mpi_plan_dft_r2c_3d(N(3),N(2),N(1),u_r,u_c,MPI_COMM_WORLD,FFTW_FORWARD,FFTW_ESTIMATE) pinv = fftw_mpi_plan_dft_c2r_3d(N(3),N(2),N(1),u_c,u_r,MPI_COMM_WORLD,FFTW_BACKWARD,FFTW_ESTIMATE) !!u_r(:,:,:) = 0.d0 !!u_r(:N(1),:,:) = 1.d0 !!u_sum_old = sum(u_r(:N(1),:,:)) !!call fftw_execute(pfor) !!call fftw_execute(pinv) !!write(*,*) u_sum_old,sum(u_r(:N(1),:,:))/(N(1)*N(2)*N(3)) !!call fftw_free(ptr_u) call fftw_mpi_cleanup() call mpi_finalize(mpi_error_flag) end program fftw3_demo

.

About the interface, right now it's just single program.

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

Navdeep R. wrote:

..

program fftw3_demo .. !! compiles and runs correctly till here. pfor = fftw_mpi_plan_dft_r2c_3d(N(3),N(2),N(1),u_r,u_c,MPI_COMM_WORLD,FFTW_FORWARD,FFTW_ESTIMATE) ..

According to this **site**, the interface for fftw_mpi_plan_dft_r2c_3d is

type(C_PTR) function fftw_mpi_plan_dft_r2c_3d(n0,n1,n2,in,out,comm,flags) bind(C, name='fftw_mpi_plan_dft_r2c_3d_f03') import integer(C_INTPTR_T), value :: n0 integer(C_INTPTR_T), value :: n1 integer(C_INTPTR_T), value :: n2 real(C_DOUBLE), dimension(*), intent(out) :: in complex(C_DOUBLE_COMPLEX), dimension(*), intent(out) :: out integer(C_MPI_FINT), value :: comm integer(C_INT), value :: flags end function fftw_mpi_plan_dft_r2c_3d

which has 7 dummy arguments, however your code has 8 which would explain, "`Error: More actual than formal arguments `

`in`

`procedure call at (1)`

". from gfortran.

Is one of FFTW_FORWARD or FFTW_ESTIMATE superfluous?

- Mark as New
- Bookmark
- Subscribe
- Mute
- Subscribe to RSS Feed
- Permalink
- Report Inappropriate Content

FortranFan wrote:

According to this

site, the interface for fftw_mpi_plan_dft_r2c_3d istype(C_PTR) function fftw_mpi_plan_dft_r2c_3d(n0,n1,n2,in,out,comm,flags) bind(C, name='fftw_mpi_plan_dft_r2c_3d_f03') import integer(C_INTPTR_T), value :: n0 integer(C_INTPTR_T), value :: n1 integer(C_INTPTR_T), value :: n2 real(C_DOUBLE), dimension(*), intent(out) :: in complex(C_DOUBLE_COMPLEX), dimension(*), intent(out) :: out integer(C_MPI_FINT), value :: comm integer(C_INT), value :: flags end function fftw_mpi_plan_dft_r2c_3dwhich has 7 dummy arguments, however your code has 8 which would explain, "Error: More actual than formal arguments in procedure call at (1)". from gfortran.

Is one of FFTW_FORWARD or FFTW_ESTIMATE superfluous?

Yes indeed. When I wrote the line, I pulled in FFTW_FORWARD/ BACKWARD from memory, which is not required here.

I will tweak the code and see if it works. Thanks a ton.

**Update : It works. So the solution is to use two different pointers, one real and other complex.**

- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Printer Friendly Page