Ошибка вычисления матрицы появляется, когда размеры становятся большими

Я запускаю код, в котором я просто создаю 2 матрицы: одна матрица имеет размеры arows x nsame, а другая имеет размеры nsame x bcols. В результате получается массив размеров arows x bcols. Это довольно просто реализовать с помощью BLAS, и следующий код, похоже, работает так, как задумано, при использовании приведенной ниже модели ведущий-ведомый с OpenMPI:

#include <iostream>
#include <stdio.h>
#include <iostream>
#include <cmath>
#include <mpi.h>
#include <gsl/gsl_blas.h>
using namespace std;`

int main(int argc, char** argv){
    int noprocs, nid;
    MPI_Status status;
    MPI_Init(&argc, &argv);
    MPI_Comm_rank(MPI_COMM_WORLD, &nid);
    MPI_Comm_size(MPI_COMM_WORLD, &noprocs);
    int master = 0;

    const int nsame = 500; //must be same if matrices multiplied together = acols = brows
    const int arows = 500;
    const int bcols = 527; //works for 500 x 500 x 527 and 6000 x 100 x 36
    int rowsent;
    double buff[nsame];
    double b[nsame*bcols];
    double c[arows][bcols];
    double CC[1*bcols]; //here ncols corresponds to numbers of rows for matrix b
    for (int i = 0; i < bcols; i++){
                CC[i] = 0.;
    }; 
    // Master part
    if (nid == master ) { 

        double a [arows][nsame]; //creating identity matrix of dimensions arows x nsame (it is I if arows = nsame)
        for (int i = 0; i < arows; i++){
            for (int j = 0; j < nsame; j++){
                if (i == j)
                    a[i][j] = 1.;
                else
                    a[i][j] = 0.;
            }
        }
        double b[nsame*bcols];//here ncols corresponds to numbers of rows for matrix b
            for (int i = 0; i < (nsame*bcols); i++){
                b[i] = (10.*i + 3.)/(3.*i - 2.) ;
            }; 
        MPI_Bcast(b,nsame*bcols, MPI_DOUBLE_PRECISION, master, MPI_COMM_WORLD);  
        rowsent=0;
        for (int i=1; i < (noprocs); i++) {  
            // Note A is a 2D array so A[rowsent]=&A[rowsent][0]
            MPI_Send(a[rowsent], nsame, MPI_DOUBLE_PRECISION,i,rowsent+1,MPI_COMM_WORLD);
            rowsent++; 
        }

        for (int i=0; i<arows; i++) { 
            MPI_Recv(CC, bcols, MPI_DOUBLE_PRECISION, MPI_ANY_SOURCE, MPI_ANY_TAG,
                     MPI_COMM_WORLD, &status); 
            int sender = status.MPI_SOURCE;
            int anstype = status.MPI_TAG;            //row number+1
            int IND_I = 0;
            while (IND_I < bcols){
                c[anstype - 1][IND_I] = CC[IND_I]; 
                IND_I++;
            }
            if (rowsent < arows) {
                MPI_Send(a[rowsent], nsame,MPI_DOUBLE_PRECISION,sender,rowsent+1,MPI_COMM_WORLD);
                rowsent++; 
            }
            else {       // tell sender no more work to do via a 0 TAG
                MPI_Send(MPI_BOTTOM,0,MPI_DOUBLE_PRECISION,sender,0,MPI_COMM_WORLD);
            }
        }
    }

    // Slave part
    else { 
        MPI_Bcast(b,nsame*bcols, MPI_DOUBLE_PRECISION, master, MPI_COMM_WORLD); 
        MPI_Recv(buff,nsame,MPI_DOUBLE_PRECISION,master,MPI_ANY_TAG,MPI_COMM_WORLD,&status); 
        while(status.MPI_TAG != 0) {
            int crow = status.MPI_TAG; 
            gsl_matrix_view AAAA = gsl_matrix_view_array(buff, 1, nsame);
            gsl_matrix_view BBBB = gsl_matrix_view_array(b, nsame, bcols);
            gsl_matrix_view CCCC = gsl_matrix_view_array(CC, 1, bcols);

            /* Compute C = A B */
            gsl_blas_dgemm (CblasNoTrans, CblasNoTrans, 1.0, &AAAA.matrix, &BBBB.matrix,
                            0.0, &CCCC.matrix); 

            MPI_Send(CC,bcols,MPI_DOUBLE_PRECISION, master, crow, MPI_COMM_WORLD);
            MPI_Recv(buff,nsame,MPI_DOUBLE_PRECISION,master,MPI_ANY_TAG,MPI_COMM_WORLD,&status); 
        }
    }

    // output c here on master node //uncomment the below lines if I wish to see the output
    //    if (nid == master){
//        if (rowsent == arows){
//            //            cout << rowsent;
//            int IND_F = 0;
//            while (IND_F < arows){
//                int IND_K = 0;
//                while (IND_K < bcols){
//                    cout << "[" << IND_F << "]" << "[" << IND_K << "] = " << c[IND_F][IND_K] << " ";
//                    IND_K++;
//                }
//                cout << "\n";
//                IND_F++;
//            }
//        }
//    }
    MPI_Finalize();
    //free any allocated space here
    return 0;
};

Теперь то, что кажется странным, заключается в том, что когда я увеличиваю размер матриц (например, с nsame = 500 до nsame = 501), код больше не работает. Я получаю следующую ошибку:

mpirun noticed that process rank 0 with PID 0 on node Users-MacBook-Air exited on signal 11 (Segmentation fault: 11).

Я пробовал это с другими комбинациями размеров для матриц, и всегда кажется, что существует верхний предел для размера самих матриц (который, кажется, зависит от того, как я изменяю сами разные размеры). Я также пробовал изменять значения самих матриц, хотя это, похоже, ничего не меняет. Я понимаю, что есть альтернативные способы инициализации матриц в моем примере (например, с использованием вектора), но мне просто интересно, почему моя текущая схема умножения матриц произвольного размера, похоже, работает только до определенной степени.