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

FFTW MPI高级接口1D实复变换段错误问题求助

avatar
云原生阿岚
2026-6-14
FFTW MPI实-复/复-实批量变换段错误问题排查

问题描述

使用FFTW高级接口结合MPI批量执行1D傅里叶变换，复-复变换运行正常，但切换为实-复（r2c）和复-实（c2r）变换后出现段错误。已调整变换计划和执行函数，但问题依旧。怀疑复数数组元素分配数量有误（尝试过N/2和N/2+1两种方式），同时注意到错误信息中出现fftw_execute_r2r，但实际调用的是fftw_execute_dft_r2c，对此存在疑问。

代码与错误信息

代码 

#include <stdio.h>
#include <stdlib.h>
#include <mpi.h>
#include <unistd.h>
#include <complex.h>
#include <fftw3-mpi.h>

int main(int argc, char *argv[]) {

/* --- MPI Init --------------------------------------------------- */
    int rank, size;
    int Nproc_x;
    int Nx, Ny;
    int nx = 3;

    MPI_Init(&argc, &argv);
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);
    MPI_Comm_size(MPI_COMM_WORLD, &Nproc_x);

    Nx = Nproc_x*nx;
    Ny = 40;


/* --- FFTW INITIALISATION ----------------------------------------------- */

    int i, N;
    int fftw_rank;
    ptrdiff_t fftw_dims[1];
    ptrdiff_t howmany, block0;
    /* local data size */
    ptrdiff_t fftw_alloc_local;
    ptrdiff_t fftw_local_nx;
    ptrdiff_t fftw_local_x_start;

    /* Complex arrays */
    double * rho_real;
    fftw_complex * rho_cmplx;

    /* Plans for FFTW transforms */
    static fftw_plan r2c_plan;
    static fftw_plan c2r_plan;

    /* Other variables */
    static double FFT_norm;

    fftw_mpi_init();

    fftw_rank    = 1;
    fftw_dims[0] = Nx;
    howmany      = Ny;
    block0 =  nx; /* Array size in radial dir on current processor */
    FFT_norm = (double) 1/Nx;

    /* get local data size and allocate data array */
    fftw_alloc_local = fftw_mpi_local_size_many(fftw_rank,
                                                fftw_dims,
                                                howmany,
                                                block0,
                                                MPI_COMM_WORLD,  /* Transform in azimuth */
                                                &fftw_local_nx,
                                                &fftw_local_x_start);

    N = fftw_alloc_local;

    rho_real      = fftw_alloc_real(N);
    rho_cmplx     = fftw_alloc_complex(N/2+1);


    /* plan for forward transformation: real density to complex density */
    r2c_plan = fftw_mpi_plan_many_dft_r2c(fftw_rank,
                                          fftw_dims,
                                          howmany,
                                          block0,
                                          FFTW_MPI_DEFAULT_BLOCK,
                                          rho_real,
                                          rho_cmplx,
                                          MPI_COMM_WORLD,
                                          FFTW_ESTIMATE);

    c2r_plan = fftw_mpi_plan_many_dft_c2r(fftw_rank,
                                          fftw_dims,
                                          howmany,
                                          FFTW_MPI_DEFAULT_BLOCK,
                                          block0,
                                          rho_cmplx,
                                          rho_real,
                                          MPI_COMM_WORLD,
                                          FFTW_ESTIMATE);

    /* Populate real array */
    for (int i = 0; i < fftw_local_nx; ++i) {
        for (int j = 0; j < Ny; ++j) {
            rho_real[i * Ny + j] = i * Ny + j ;
        }
    }


    fftw_mpi_execute_dft_r2c(r2c_plan, rho_real, rho_cmplx);
    fftw_mpi_execute_dft_c2r(c2r_plan, rho_cmplx, rho_real);

    printf("INFO(gravBessel): fftw execute success ! \n");



    MPI_Finalize();

    return 0;
}

错误信息 

[euler:12536] *** Process received signal ***
[euler:12536] Signal: Segmentation fault (11)
[euler:12536] Signal code: Address not mapped (1)
[euler:12536] Failing at address: (nil)
[euler:12536] [ 0] /lib/x86_64-linux-gnu/libc.so.6(+0x43090)[0x7f2f7ff04090]
[euler:12536] [ 1] /lib/x86_64-linux-gnu/libfftw3.so.3(fftw_execute_r2r+0x4)[0x7f2f802e59d4]
[euler:12536] [ 2] my_program(+0x1509)[0x55da6283a509]
[euler:12536] [ 3] /lib/x86_64-linux-gnu/libc.so.6(__libc_start_main+0xf3)[0x7f2f7fee5083]
[euler:12536] [ 4] my_program(+0x122e)[0x55da6283a22e]
[euler:12536] *** End of error message ***

问题根源与修复方案

1. 复数数组分配错误

fftw_mpi_local_size_many仅适用于复-复变换的本地大小计算，实-复/复-实变换的复数数组需要使用专用函数fftw_mpi_local_size_many_dft_r2c计算正确的本地元素数量——N/2+1是全局复数数组的大小逻辑，不适用于MPI分布式场景，直接推导会导致数组空间不足或冗余，触发越界访问。

修复代码：

// 替换原复数数组分配部分
ptrdiff_t fftw_local_nx_cmplx, fftw_local_x_start_cmplx;
ptrdiff_t fftw_alloc_local_cmplx = fftw_mpi_local_size_many_dft_r2c(
    fftw_rank, fftw_dims, howmany,
    block0, FFTW_MPI_DEFAULT_BLOCK,
    MPI_COMM_WORLD,
    &fftw_local_nx_cmplx, &fftw_local_x_start_cmplx
);
rho_cmplx = fftw_alloc_complex(fftw_alloc_local_cmplx);

2. 数组初始化索引错误

原代码中数组初始化的索引顺序rho_real[i * Ny + j]颠倒了变换维度与批量维度的存储顺序。FFTW的many类变换采用批量优先的存储布局：每个批量（howmany）的完整变换数据连续存储，即第j个批量的第i个本地元素应为rho_real[j * fftw_local_nx + i]。原索引会导致数组越界写入，触发段错误。

修复代码：

/* Populate real array (修正索引顺序) */
for (int j = 0; j < Ny; ++j) {
    for (int i = 0; i < fftw_local_nx; ++i) {
        rho_real[j * fftw_local_nx + i] = j * fftw_local_nx + i;
    }
}

3. 关于fftw_execute_r2r的疑问

错误信息中出现fftw_execute_r2r属于正常现象：FFTW内部会将实-复/复-实变换转换为实-实变换的特殊类型实现，无需额外处理，此并非问题根源。

4. 额外建议：添加计划创建校验

计划创建失败会返回NULL，直接调用执行函数会导致崩溃，建议添加校验逻辑：

if (!r2c_plan || !c2r_plan) {
    fprintf(stderr, "Rank %d: Failed to create FFTW plans\n", rank);
    MPI_Abort(MPI_COMM_WORLD, 1);
}

修复后完整代码示例 

#include <stdio.h>
#include <stdlib.h>
#include <mpi.h>
#include <unistd.h>
#include <complex.h>
#include <fftw3-mpi.h>

int main(int argc, char *argv[]) {

/* --- MPI Init --------------------------------------------------- */
    int rank, size;
    int Nproc_x;
    int Nx, Ny;
    int nx = 3;

    MPI_Init(&argc, &argv);
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);
    MPI_Comm_size(MPI_COMM_WORLD, &Nproc_x);

    Nx = Nproc_x*nx;
    Ny = 40;


/* --- FFTW INITIALISATION ----------------------------------------------- */

    int i;
    int fftw_rank;
    ptrdiff_t fftw_dims[1];
    ptrdiff_t howmany, block0;
    /* local data size */
    ptrdiff_t fftw_alloc_local;
    ptrdiff_t fftw_local_nx;
    ptrdiff_t fftw_local_x_start;

    /* Complex arrays */
    double * rho_real;
    fftw_complex * rho_cmplx;

    /* Plans for FFTW transforms */
    fftw_plan r2c_plan;
    fftw_plan c2r_plan;

    /* Other variables */
    double FFT_norm;

    fftw_mpi_init();

    fftw_rank    = 1;
    fftw_dims[0] = Nx;
    howmany      = Ny;
    block0 =  nx; /* Array size in radial dir on current processor */
    FFT_norm = (double) 1/Nx;

    /* get local data size for real array */
    fftw_alloc_local = fftw_mpi_local_size_many(fftw_rank,
                                                fftw_dims,
                                                howmany,
                                                block0,
                                                MPI_COMM_WORLD,
                                                &fftw_local_nx,
                                                &fftw_local_x_start);

    rho_real = fftw_alloc_real(fftw_alloc_local);

    /* get local data size for complex array (r2c专用) */
    ptrdiff_t fftw_alloc_local_cmplx;
    ptrdiff_t fftw_local_nx_cmplx;
    ptrdiff_t fftw_local_x_start_cmplx;
    fftw_alloc_local_cmplx = fftw_mpi_local_size_many_dft_r2c(
        fftw_rank, fftw_dims, howmany,
        block0, FFTW_MPI_DEFAULT_BLOCK,
        MPI_COMM_WORLD,
        &fftw_local_nx_cmplx, &fftw_local_x_start_cmplx
    );
    rho_cmplx = fftw_alloc_complex(fftw_alloc_local_cmplx);


    /* plan for forward transformation: real density to complex density */
    r2c_plan = fftw_mpi_plan_many_dft_r2c(fftw_rank,
                                          fftw_dims,
                                          howmany,
                                          block0,
                                          FFTW_MPI_DEFAULT_BLOCK,
                                          rho_real,
                                          rho_cmplx,
                                          MPI_COMM_WORLD,
                                          FFTW_ESTIMATE);

    c2r_plan = fftw_mpi_plan_many_dft_c2r(fftw_rank,
                                          fftw_dims,
                                          howmany,
                                          FFTW_MPI_DEFAULT_BLOCK,
                                          block0,
                                          rho_cmplx,
                                          rho_real,
                                          MPI_COMM_WORLD,
                                          FFTW_ESTIMATE);

    /* 校验计划是否创建成功 */
    if (!r2c_plan || !c2r_plan) {
        fprintf(stderr, "Rank %d: Failed to create FFTW plans\n", rank);
        MPI_Abort(MPI_COMM_WORLD, 1);
    }

    /* Populate real array (修正索引顺序) */
    for (int j = 0; j < Ny; ++j) {
        for (int i = 0; i < fftw_local_nx; ++i) {
            rho_real[j * fftw_local_nx + i] = j * fftw_local_nx + i;
        }
    }


    fftw_mpi_execute_dft_r2c(r2c_plan, rho_real, rho_cmplx);
    fftw_mpi_execute_dft_c2r(c2r_plan, rho_cmplx, rho_real);

    printf("INFO(gravBessel): fftw execute success ! Rank %d\n", rank);

    /* 释放资源 */
    fftw_destroy_plan(r2c_plan);
    fftw_destroy_plan(c2r_plan);
    fftw_free(rho_real);
    fftw_free(rho_cmplx);

    MPI_Finalize();

    return 0;
}

内容的提问来源于stack exchange，提问作者steven

相关产品推荐
方舟 Agent Plan

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

产品介绍相关文档
相关产品推荐
方舟 Agent Plan

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

产品介绍相关文档
火山引擎 最新活动
全面开放
Doubao-Seedance-2.5全新模型升级

模型再升级,30秒超长叙事, 模态参考扩容

HOT
方舟Codingplan

模型自由，工具不限，最新支持 Deepseek-V4 系列、GLM-5.3 系列

HOT
方舟Agent plan

超全模态模型 × Harness 升级，最新支持 Deepseek-V4.1-Flash、GLM-5.3 系列

HOT
新用户特惠专场

大模型19元起，Al应用9.9元畅享，新人首购爆款尽享优惠

最近更新时间：2026.06.14 06:05:53