Multithreaded FFTW can be implemented as in this from FFTW homepage. Instead, we want to call the serial FFTW routines within OpenMP parallel environment using multiple threads. We want variable and fourier_variable to be thread-safe. This could be done by using PURE subroutines and declaring variable and fourier_variable inside it. The question here is related to calling FFTW routines like fftw_execute_dft_r2c from within a PURE subroutine.
A stripped-down version of the code is presented here just for reference (the full code is an optimisation solver involving many FFTW calls).
PROGRAM main
USE fft_call
REAL(8), DIMENSION(1:N, 1:N) :: variable
COMPLEX(C_DOUBLE_COMPLEX), DIMENSION(1:N/2+1, 1:N) :: fourier_variable
INTEGER :: JJ
!$OMP PARALLEL
!$OMP DO
DO JJ = 1, 5
call fourier_to_physical(fourier_variable, variable)
END DO
!$OMP END DO
!$OMP END PARALLEL
END PROGRAM main
MODULE fft_call
contains
PURE SUBROUTINE fourier_to_physical( fourier_variable, variable)
IMPLICIT NONE
REAL(8), DIMENSION(1:N, 1:N) :: variable
COMPLEX(C_DOUBLE_COMPLEX), DIMENSION(1:N/2+1, 1:N), INTENT(OUT) :: fourier_variable
CALL fftw_execute_dft_r2c (plan_fwd, variable, fourier_variable)
END SUBROUTINE fourier_to_physical
END MODULE fft_call
The error while calling fftw_plan_dft_r2c_2d from the PURE subroutine fourier_to_physical:
Error: Function reference to 'fftw_plan_dft_r2c' at (1) is to a non-PURE procedure within a PURE procedure
The question: is there a way to call FFTW routines like fftw_execute_dft_r2c from within a PURE subroutine in Fortran90?
Or, in other words, are their PURE versions of fftw_execute_dft_r2c such that we can call them from PURE procedures? We are beginners to OpenMP.