通过MPI传递一个Armadillo C++矩阵

Question

我需要通过MPI传递由Armadillo C++ Matrix Library定义的矩阵或复杂矩阵类型。有什么好的方法吗？我想尝试：

1. 将矩阵写入某种类型的数组，然后发送该数组的行/列，并使用在MPI_send/recv两侧重新构造数组的方法。 2. 使用类似MPI_BYTE类型的内容？

谢谢。

更新

因此，我正在尝试在一个节点上为简单示例实现另一种方案，即发送和接收。

translate.cpp

    #include <mpi.h>
    #include <armadillo>
    #include <vector> 
    #include <cstdlib>

    using namespace std; 
    using namespace arma; 
    using std::vector; 

    class ArmadilloMPI
    {
        public:
            ArmadilloMPI(int nRows, int nCols)
            {
                this->nRows = nRows;
                this->nCols = nCols; 
                realArray = (double **)malloc(nCols * nRows * sizeof(double*));
                imArray = (double **)malloc(nCols * nRows * sizeof(double*));
            }

            ~ArmadilloMPI()
            {
                free(realArray[0]);
                free(realArray);
                free(imArray[0]);
                free(imArray);
            }

            double **realArray; 
            double **imArray; 
            int nCols; 
            int nRows; 

            cx_mat matConstructRecv(int src, int tag)
            {
                cx_mat A(nRows, nCols); 
                MPI_Recv(&(imArray[0][0]),  nRows * nCols, MPI_DOUBLE, src, tag, MPI_COMM_WORLD,0);
                MPI_Recv(&(realArray[0][0]),nRows * nCols, MPI_DOUBLE, src, tag, MPI_COMM_WORLD,0);

                for(int  i = 0; i < nRows; ++i )
                {
                    for(int j = 0; i < nCols; ++j)
                    {
                        real(A(i,j)) = *realArray[i * nRows + j]; 
                        imag(A(i,j)) = *imArray[i * nRows + j];
                    }
                }
                return A; 
            }

            void matDestroySend(cx_mat &A, int dest, int tag)
            {
                for(int  i = 0; i < nRows; ++i )
                {
                    for(int j = 0; i < nCols; ++j)
                    {
                        realArray[i * nRows + j]  = &real(A(i,j)); 
                        imArray[i * nRows + j] = &imag(A(i,j)); 
                    }
                }
                MPI_Send(&(realArray[0][0]), nRows * nCols, MPI_DOUBLE, dest, tag, MPI_COMM_WORLD);
                MPI_Send(&(imArray[0][0]), nRows * nCols, MPI_DOUBLE, dest, tag, MPI_COMM_WORLD);
            }
    };

    int main(int argc, char** argv)
    {
        MPI::Init(argc, argv);

        int size = MPI::COMM_WORLD.Get_size();
        int rank = MPI::COMM_WORLD.Get_rank();

        cout << "test"<<endl; 
        vector<cx_mat> world; 
        for(int i = 0; i < size; ++i )
        {
            world.push_back(randu<cx_mat>(4,4));
        }
        cx_mat A;
        A = randu<cx_mat>(4,4);

        ArmadilloMPI* armaMPI = new ArmadilloMPI(4,4); 

        if(rank==0)
        {

            for(int i = 1; i < size; i++)
            {   
                cout << "A is now " << A << endl; 
                A += armaMPI->matConstructRecv(i, 0);
            }
        }
        else
        {
            armaMPI->matDestroySend(world[rank], 1, 0);
        }

        cout << A << endl; 
        delete armaMPI;
        MPI::Finalize();
    }

但是我们遇到了段错误。

*** Process received signal *** 
Signal: Segmentation fault: 11 (11) 
Signal     code:     (0) 
Failing at address: 0x0 translate(1032,0x7fff747ad310) malloc: ***   error for object     0x41434d5f49504d4f: pointer being freed was not allocated

你的想法？

- user3819317

---

我尝试了 MPI_Bcast( &AAA(0,0), sizex*sizey*64, MPI_BYTE, 0, MPI_COMM_WORLD);，但是一旦 AAA 大于16x16左右，就会触发分段错误。 我会选择第一个选项...您可能会对PETSc库的MATMPIDENSE矩阵类型感兴趣。这里有一个例子。 - francis

Answer 1

这里有几个问题:

• In c and c++, array and vector start at 0, not 1. So the following code will fail :

   vector<cx_mat> world; 
   world.resize(1);
   world[1] = randu<cx_mat>(4,4); //problem to come !
  

  You may change for :

  vector<cx_mat> world;
  world.push_back(randu<cx_mat>(4,4));
  

• Dynamic allocation of 2D array with contiguous memory. You need one new for an array of double, and another new for array of pointers to double. Then set each pointer to point to the first item of the row.

  double *data=new double[nCols * nRows ];
  realArray = new double*[( nRows )];
  for(int i=0;i<nRows;i++){
       realArray[i]=&data[i*nCols];
  }
  

• You could guess this one...Why don't compilers warn about this kind of stuff ? Because it could make sense, but not here.

  for(int j = 0; i < nCols; ++j)
  

• You may add a different tag to each message to avoid switching the real part and the imaginary part

  MPI_Send(&(realArray[0][0]), nRows * nCols, MPI_DOUBLE, dest, tag, MPI_COMM_WORLD);
  MPI_Send(&(imArray[0][0]), nRows * nCols, MPI_DOUBLE, dest, tag+1, MPI_COMM_WORLD);
  

代码变为：

#include <mpi.h>
#include <armadillo>
#include <vector>
#include <iostream>
#include <cstdlib>

using namespace std;
using namespace arma;
using std::vector;

class ArmadilloMPI
{
public:
    ArmadilloMPI(int nRows, int nCols)
    {
        this->nRows = nRows;
        this->nCols = nCols;
        double *data=new double[nCols * nRows ];
        realArray = new double*[( nRows )];
        for(int i=0;i<nRows;i++){
            realArray[i]=&data[i*nCols];
        }
        double *datai=new double[(nCols * nRows )];
        imArray =new double*[( nRows )];
        for(int i=0;i<nRows;i++){
            imArray[i]=&datai[i*nCols];
        }

    }

    ~ArmadilloMPI()
    {
        delete[] realArray[0];
        delete[] realArray;
        delete[] imArray[0];
        delete[] imArray;
    }

    double **realArray;
    double **imArray;
    int nCols;
    int nRows;

    cx_mat matConstructRecv(int tag, int src)
    {
        cx_mat A(nRows, nCols);
        MPI_Recv(&(imArray[0][0]),  nRows * nCols, MPI_DOUBLE, src, tag+1, MPI_COMM_WORLD,0);
        MPI_Recv(&(realArray[0][0]),nRows * nCols, MPI_DOUBLE, src, tag, MPI_COMM_WORLD,0);

        for(int  i = 0; i < nRows; ++i )
        {
            for(int j = 0; j < nCols; ++j)
            {
                real(A(i,j)) = realArray[i][j];
                imag(A(i,j)) = imArray[i][j];
            }
        }
        return A;
    }

    void matDestroySend(cx_mat &A, int dest, int tag)
    {
        for(int  i = 0; i < nRows; ++i )
        {
            for(int j = 0; j < nCols; ++j)
            {
                realArray[i][j] = real((A(i,j)));
                imArray[i][j] = imag((A(i,j)));
            }
        }
        MPI_Send(&(realArray[0][0]), nRows * nCols, MPI_DOUBLE, dest, tag, MPI_COMM_WORLD);
        MPI_Send(&(imArray[0][0]), nRows * nCols, MPI_DOUBLE, dest, tag+1, MPI_COMM_WORLD);
    }
};

int main(int argc, char **argv)
{
    int rank;
    int size;
    MPI_Init(&argc, &argv);
    MPI_Comm_size(MPI_COMM_WORLD, &size);
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);

    srand (time(NULL)+rank);

    vector<cx_mat> world;
    world.push_back(randu<cx_mat>(4,4));

    cx_mat A;
    ArmadilloMPI* armaMPI = new ArmadilloMPI(4,4);
    if(rank==0)
    {
        world[0].print("world[0] on 0:");

        armaMPI->matDestroySend(world[0], 1, 0);
    }
    if(rank==1){
        A = armaMPI->matConstructRecv(0, 0);
        A.print("A on 1:");
    }

    delete armaMPI;
    MPI_Finalize();
}

编译代码的步骤如下：

 mpiCC -O2 -o main main.cpp -larmadillo -llapack -lblas -Wall

运行：

mpiexec -np 2 main

- francis

---

非常感谢。嗯，有趣的是，直到发布后你才会注意到所有你没有注意到的事情。 - user3819317