Как сделать «транспонированное» двумерное преобразование fftw3 MPI, если это вообще возможно?
Рассмотрим двумерное преобразование в форме L x M (настройка основного столбца) из сложного массиваЦСИ в реальный массивTGT, Или, на фортранском,
complex(C_DOUBLE_COMPLEX), pointer :: src(:,:)
real(8), pointer :: tgt(:,:) .
Соответствующие указатели
type(C_PTR) :: csrc,ctgt .
Я бы выделил их следующим образом:
! The complex array first
alloc_local = fftw_mpi_local_size_2d(M,L/2+1,MPI_COMM_WORLD,local_M,local_offset1)
csrc = fftw_alloc_complex(alloc_local)
call c_f_pointer(csrc, src, [L/2,local_M])
! Now the real array
alloc_local = fftw_mpi_local_size_2d(2*(L/2+1),M, &
MPI_COMM_WORLD,local_L,local_offset2)
ctgt = fftw_alloc_real(alloc_local)
call c_f_pointer(ctgt, tgt, [M,local_L])
Теперь план будет создан как:
! Create c-->r transform with one transposition left out
plan = fftw_mpi_plan_dft_c2r_2d(M,L,src,tgt, MPI_COMM_WORLD, &
ior(FFTW_MEASURE,FFTW_MPI_TRANSPOSED_OUT))
Наконец, преобразование будет выполнено как:
call fftw_mpi_execute_dft_c2r(plan, src, tgt)
Однако этот рецепт не работает. Последний вызов вызывает ошибку сегментации. Сначала я подумал, что это может быть связано с тем, как яЦСИ а такжеTGT массивы, но игра с другим объемом памяти, выделенной дляTGT не дал никакого результата. Итак, я либо делаю что-то действительно глупое, либо это невозможно вообще.
РЕДАКТИРОВАТЬ: МИНИМАЛИСТИЧЕСКИЙ СОВМЕСТИМЫЙ ПРИМЕР
program trashingfftw
use, intrinsic :: iso_c_binding
use MPI
implicit none
include 'fftw3-mpi.f03'
integer(C_INTPTR_T), parameter :: L = 256
integer(C_INTPTR_T), parameter :: M = 256
type(C_PTR) :: plan, ctgt, csrc
complex(C_DOUBLE_COMPLEX), pointer :: src(:,:)
real(8), pointer :: tgt(:,:)
integer(C_INTPTR_T) :: alloc_local, local_M, &
& local_L,local_offset1,local_offset2
integer :: ierr,id
call mpi_init(ierr)
call mpi_comm_rank(MPI_COMM_WORLD,id,ierr)
call fftw_mpi_init()
alloc_local = fftw_mpi_local_size_2d(M,L/2+1, MPI_COMM_WORLD, &
local_M, local_offset1)
csrc = fftw_alloc_complex(alloc_local)
call c_f_pointer(csrc, src, [L/2,local_M])
alloc_local = fftw_mpi_local_size_2d(2*(L/2+1),M, MPI_COMM_WORLD, &
& local_L, local_offset2)
ctgt = fftw_alloc_real(alloc_local)
call c_f_pointer(ctgt, tgt, [M,local_L])
plan = fftw_mpi_plan_dft_c2r_2d(M,L,src,tgt, MPI_COMM_WORLD, &
ior(FFTW_MEASURE, FFTW_MPI_TRANSPOSED_OUT))
call fftw_mpi_execute_dft_c2r(plan, src, tgt)
call mpi_finalize(ierr)
end program trashingfftw