You need to enable JavaScript to run this app.
优惠活动
大模型
产品
解决方案
定价
更多

FFTW3高级与基础接口2D C2R变换结果不一致问题排查

FFTW3 Fortran接口多平面C2R变换异常问题

我正在维护一套采用伪谱求解器的遗留CFD代码,计划通过替换更高效的FFT算法实现升级,目前正在学习FFTW3的Fortran接口。我有DNS代码输出的多种格式3D速度数组,用于验证FFT变换的正确性。

当前需求是对一个第一维度为物理空间、后两维度为谱空间的3D数组执行2D逆傅里叶变换(C2R),得到全物理空间的3D数组,但代码运行出现问题:

  • 最初用fftw_plan_dft_c2r成功完成单平面2D变换,但遍历所有平面时结果异常;尝试保留输入数据时触发段错误
  • 转而使用高级接口fftw_plan_many_dft_c2r批量处理多平面DFT,先测试单平面场景:预期输出为mz×mx的约150常量数组,基础接口fftw_plan_dft_c2r符合预期,但高级接口输出异常数值
  • 尝试调整istride/ostride、idist/odist、inembed/onembed等参数后问题仍未解决;更新代码明确网格尺寸和谱数组定义后,原本正常的2D变换也失效了

附测试代码:

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

内容的提问来源于stack exchange,提问作者Ryan Kelly

相关产品推荐
方舟 Agent Plan

超全模态模型 × Harness 升级,最新支持 Deepseek-V4.1-Flash、GLM-5.3 系列、Doubao-Seedream-5.0-pro、Kimi-K3 (部分), 限时 9.9 元起

最近更新时间:2026.06.24 14:19:51