MPI矩阵向量乘法程序部分维度下出现段错误问题求助

问题根因 #

你触发段错误的核心原因是栈溢出，除此之外代码还存在多处逻辑bug和内存浪费问题。

具体问题点 #

• 栈溢出触发段错误
  matmul函数中定义的float matrix [n][n]是GCC扩展支持的变长数组（VLA），默认分配在程序栈空间。通常操作系统默认栈大小仅为8MB左右，当n=2880时，该数组占用内存为2880*2880*4B = 33.18MB，直接超出栈上限触发段错误；小维度n=720时占用720*720*4B=2.07MB未超出栈上限，所以可以正常运行。
• 内存严重浪费
  按行拆分矩阵的场景下，每个进程只需要存储n/size行的局部矩阵，你现在给每个进程都分配了完整的n*n大小的矩阵，内存占用是实际需要的size倍，大维度下极易触发内存不足。同时main函数中所有进程都分配了完整的全局矩阵A，只有rank0进程需要持有全局矩阵，其他进程完全不需要分配这块内存。
• 内存泄漏
  matmul函数中new申请的localresult和局部矩阵没有对应的delete释放，运行时间长后会持续占用内存。
• 兼容性问题
  变长数组不是C++标准语法，仅GCC等少数编译器支持，代码可移植性极差。
• 边界逻辑bug
  当前MPI_Scatter和MPI_Gather的计数都是硬编码的n/size，如果n不能被进程数整除会出现数据截断、结果错误。

修复方案 #

核心修改点 #

1. 把栈上分配的变长数组改成堆上分配，且只分配进程实际需要的局部矩阵大小
2. 释放所有new申请的堆内存避免泄漏
3. 仅在rank0进程分配和初始化全局矩阵A，其他进程不分配这块内存
4. （可选）改用MPI_Scatterv/MPI_Gatherv兼容n不能被进程数整除的场景

修正后的核心代码示例 #

#include <iostream>
#include <cmath>
#include <cstdlib>
#include <chrono>
#include <mpi.h>
float genA (int row, int col) {
    if (row > col)
        return 1.;
    else
        return 0.;
}

float genx0 (int i) {
    return 1.;
}

void checkx (int iter, long i, float xval) {
    if (iter == 1) {
        float shouldbe = i;
        if (fabs(xval/shouldbe) > 1.01 || fabs(xval/shouldbe) < .99 )
            std::cout<<"incorrect : x["<<i<<"] at iteration "<<iter<<" should be "<<shouldbe<<" not "<<xval<<std::endl;
    }

    if (iter == 2) {
        float shouldbe =(i-1)*i/2;
        if (fabs(xval/shouldbe) > 1.01 || fabs(xval/shouldbe) < .99)
            std::cout<<"incorrect : x["<<i<<"] at iteration "<<iter<<" should be "<<shouldbe<<" not "<<xval<<std::endl;
    }
}

//perform dense y=Ax on an n \times n matrix
void matmul(float*A, float*x, float*y, long n) {
    int rank;
    int size;
    MPI_Comm_rank (MPI_COMM_WORLD, &rank);
    MPI_Comm_size (MPI_COMM_WORLD, &size); 
    long local_rows = n / size;
    float* localresult = new float[local_rows]{};
    // 堆上分配仅需要的局部矩阵大小，不再用栈上VLA
    float* matrix = new float[local_rows * n];
    MPI_Scatter(A, local_rows * n, MPI_FLOAT, matrix, local_rows * n, MPI_FLOAT, 0, MPI_COMM_WORLD);
    MPI_Bcast(x, n, MPI_FLOAT, 0, MPI_COMM_WORLD);
    for (long row = 0; row < local_rows; row++) {
        for (long col = 0; col < n; ++col) {
            localresult[row] += x[col] * matrix[row * n + col];
        }
    }
    MPI_Gather(localresult, local_rows, MPI_FLOAT, y, local_rows, MPI_FLOAT, 0, MPI_COMM_WORLD);
    // 释放堆内存
    delete[] localresult;
    delete[] matrix;
}

int main (int argc, char*argv[]) {
    if (argc < 3) {
        std::cout<<"usage: "<<argv[0]<<" <n> <iteration>"<<std::endl;
        return 1; // 缺少参数直接退出，避免后续逻辑崩溃
    }
    bool check = true;
    long n = atol(argv[1]);
    long iter = atol(argv[2]);
    int my_rank;
    int nprocs;
    MPI_Init(&argc,&argv);
    MPI_Comm_size(MPI_COMM_WORLD,&nprocs);
    MPI_Comm_rank(MPI_COMM_WORLD,&my_rank); 
    
    float* A = nullptr;
    // 仅rank0分配和初始化全局矩阵A
    if (my_rank == 0) {
        A = new float[n*n];
        for (long row = 0; row<n; row++) {
            for (long col=0; col<n; col++) {
                A[row*n+col] = genA(row, col);
            }
        }
    }
    
    float* x = new float[n];
    for (long i=0; i<n; ++i)
        x[i] = genx0(i);
    float* y = new float[n];
    
    std::chrono::time_point<std::chrono::system_clock> start = std::chrono::system_clock::now();

    for (int it = 0; it<iter; ++it) {
        matmul(A, x, y, n);
        {
            float*t = x;
            x=y;
            y=t;
        }
        if (check && my_rank == 0) // 仅rank0做校验，避免重复输出
            for (long i = 0; i<n; ++i)
                checkx (it+1, i, x[i]);
    }
    
    if (my_rank == 0) {
        std::chrono::time_point<std::chrono::system_clock> end = std::chrono::system_clock::now();
        std::chrono::duration<double> elapsed_seconds = end-start;
        std::cerr<<elapsed_seconds.count()<<std::endl;
        delete[] A;
    }
    delete[] x;
    delete[] y;
    MPI_Finalize();
    return 0;
}

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