我正在使用使用伪谱求解器的遗留CFD代码，最终我想用更好的FFT算法更新代码，所以我目前正在学习Fortran中的FFTW3接口。我有一个来自 DNS 代码的各种形式的 3D 速度数组，用于测试 FFT 并确保我正确地转换内容。
目前，我有一个 3D 数组，它在第一维是物理的，但在其他两个维度是光谱的。我的目标是对此数组进行 2D 转换以获得 3D 物理值，但我一直在试图找出当前代码的问题所在。
我已经使用
fftw_plan_dft_c2r()我想我可能最好使用高级接口在飞机上进行多个 DFT，使用
fftw_plan_many_dft_c2r()mzmx 9576957.89326171 -6100449.91245205 4941.42749753309
-2029633.19975117 378.904251844637 -1219266.13730242
1476.66846120109 -868929.573548530 -20.9026499204410
-676593.055185564 450.197926313034 -552121.186219578
-425.926031798779 -467648.354118945 -220.256552425897
...
我真的不确定我还能改变什么。出于好奇，我尝试了 {i,o}stride、{i,o}dist 和 {i,o}nembed，但似乎没有什么能让我更接近正确的答案。
代码如下：
program fftw_test
use, intrinsic :: iso_c_binding
implicit none
include 'fftw3.f03'
! Grid sizes
integer, parameter :: nxh = 128, nz = 128
integer, parameter :: mx = 384, mz = 192
! FFTW Variables
type(C_PTR) :: plan2,plan3
complex(C_DOUBLE_COMPLEX), dimension(nz,nxh) :: u_phys
complex(C_DOUBLE_COMPLEX), dimension(nz,nxh) :: u1_2d,u3_2d
real(C_DOUBLE), dimension(mz,mx) :: u2_2d,u4_2d
integer(C_INT), dimension(2) :: msize_2d = (mx,mz)
integer(C_INT) :: rank, howmany, idist, odist, istride, ostride
! ------------------------------------------------------------------------------------ !
! ------------------------------------------------------------------------------------ !
! Create FFT plan (before initializing variable)
print *,'Create FFT Plan(s)'
rank = 2; howmany = 1; ! 1 2D c2r transform on an (mz) x (mx) array
! (msize_2d is flipped because of column-major arrays)
idist = 0; odist = 0; ! Unused since howmany = 1
istride = 1; ostride = 1; ! Array is contiguous in memory (I think?)
plan2 = fftw_plan_dft_c2r (rank, msize_2d,u1_2d,u2_2d,FFTW_ESTIMATE)
plan3 = fftw_plan_many_dft_c2r(rank, msize_2d,howmany, &
u3_2d,msize_2d,istride,idist, &
u4_2d,msize_2d,ostride,odist, FFTW_ESTIMATE)
u_phys = (0.0,0.0)
u_phys(1,1) = (150.569261744544,1.136868377216160E-013)
u1_2d = u_phys ! temp storage
u3_2d = u_phys ! temp storage
! Execute FFTs
print *,'Execute FFTs'
call fftw_execute_dft_c2r(plan2,u1_2d,u2_2d)
write(200,*) u2_2d
call fftw_execute_dft_c2r(plan3,u3_2d,u4_2d)
write(201,*) u4_2d
! Destroy FFT plans
call fftw_destroy_plan(plan2)
call fftw_destroy_plan(plan3)
end program fftw_test
编辑：我更新了代码以显式定义网格大小变量和要转换的样本光谱数组。然而，由于某种原因，尽管我之前的验证没有更改代码中的任何内容，但 2D 变换不再对我起作用。
为了完整起见，我将分享我对此问题的解决方案作为答案。以下是稍作修改的代码，确实给出了预期的输出：
program fftw_test
use,intrinsic :: iso_c_binding
implicit none
include 'fftw3.f03'
! Grid sizes
integer,parameter :: nz = 128, nxh = 128
integer,parameter :: mz = 192, mx = 384
! FFTW Plans
type(C_PTR) :: plan1,plan2
! Complex arrays
complex(C_DOUBLE_COMPLEX), dimension(nz,nxh) :: u_phys ! input array
complex(C_DOUBLE_COMPLEX), dimension(mz,mx) :: u1_spec_2D
complex(C_DOUBLE_COMPLEX), dimension(mz,mx) :: u2_spec_2D
! Real arrays
real(C_DOUBLE), dimension(mz,mx) :: u1_real_2D ! 2D c2r transform target array
real(C_DOUBLE), dimension(mz,mx) :: u2_real_2D ! 2D c2r transform target array
! fft planl vars
integer(C_INT) :: rank, howmany, idist, odist, istride, ostride
integer(C_INT), dimension(2) :: msz = (/mz,mx/)
! ------------------------------------------------------------------------------------ !
! ------------------------------------------------------------------------------------ !
! Create FFT plan
rank = 2; howmany = 1
idist = 0; odist = 0
istride = 1; ostride = 1
plan1 = fftw_plan_dft_c2r(rank,msz,u1_spec_2D,u1_real_2D,FFTW_ESTIMATE)
plan2 = fftw_plan_many_dft_c2r(rank,msz,howmany, &
u2_spec_2D,msz,istride,idist, &
u2_real_2D,msz,ostride,odist,FFTW_ESTIMATE)
u_phys = (0.0,0.0)
u_phys(1,1) = (150.569261744544,1.136868377216160E-013)
u1_spec_2D(1:nz,1:nxh) = u_phys
u2_spec_2D(1:nz,1:nxh) = u_phys
! Execute FFT
call fftw_execute_dft_c2r(plan1,u1_spec_2D,u1_real_2D)
call fftw_execute_dft_c2r(plan2,u2_spec_2D,u2_real_2D)
write(100,*) u1_real_2D
write(200,*) u2_real_2D
! Destroy plans
call fftw_destroy_plan(plan1)
call fftw_destroy_plan(plan2)
end program fftw_test
除了一些名称更改之外，唯一的区别是我在要转换的复杂数组中添加了额外的空间
nznxhmzmx