C++和MPI如何将部分代码编写为并行？

Question

我一直在使用PETSc库编写一些代码，现在我要改变它的一部分以并行方式运行。大多数我想要并行化的东西是矩阵初始化和我生成和计算大量值的部分。无论如何，我的问题在于，如果我运行的代码超过1核心，出于某种原因，代码的所有部分的运行次数与我使用的核心数量一样多。

这是我测试的PETSc和MPI

int main(int argc, char** argv) 
{ 
    time_t rawtime; 
    time (&rawtime); 
    string sta = ctime (&rawtime); 
    cout << "Solving began..." << endl; 

PetscInitialize(&argc, &argv, 0, 0); 

    Mat   A;  /* linear system matrix */ 
    PetscInt  i,j,Ii,J,Istart,Iend,m = 120000,n = 3,its; 
    PetscErrorCode ierr; 
    PetscBool  flg = PETSC_FALSE; 
    PetscScalar v; 
#if defined(PETSC_USE_LOG) 
    PetscLogStage stage; 
#endif 

    /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 
     Compute the matrix and right-hand-side vector that define 
     the linear system, Ax = b. 
    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ 
    /* 
    Create parallel matrix, specifying only its global dimensions. 
    When using MatCreate(), the matrix format can be specified at 
    runtime. Also, the parallel partitioning of the matrix is 
    determined by PETSc at runtime. 

    Performance tuning note: For problems of substantial size, 
    preallocation of matrix memory is crucial for attaining good 
    performance. See the matrix chapter of the users manual for details. 
    */ 
    ierr = MatCreate(PETSC_COMM_WORLD,&A);CHKERRQ(ierr); 
    ierr = MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,m,n);CHKERRQ(ierr); 
    ierr = MatSetFromOptions(A);CHKERRQ(ierr); 
    ierr = MatMPIAIJSetPreallocation(A,5,PETSC_NULL,5,PETSC_NULL);CHKERRQ(ierr); 
    ierr = MatSeqAIJSetPreallocation(A,5,PETSC_NULL);CHKERRQ(ierr); 
    ierr = MatSetUp(A);CHKERRQ(ierr); 

    /* 
    Currently, all PETSc parallel matrix formats are partitioned by 
    contiguous chunks of rows across the processors. Determine which 
    rows of the matrix are locally owned. 
    */ 
    ierr = MatGetOwnershipRange(A,&Istart,&Iend);CHKERRQ(ierr); 

    /* 
    Set matrix elements for the 2-D, five-point stencil in parallel. 
     - Each processor needs to insert only elements that it owns 
     locally (but any non-local elements will be sent to the 
     appropriate processor during matrix assembly). 
     - Always specify global rows and columns of matrix entries. 

    Note: this uses the less common natural ordering that orders first 
    all the unknowns for x = h then for x = 2h etc; Hence you see J = Ii +- n 
    instead of J = I +- m as you might expect. The more standard ordering 
    would first do all variables for y = h, then y = 2h etc. 

    */ 
PetscMPIInt rank;  // processor rank 
PetscMPIInt size;  // size of communicator 
MPI_Comm_rank(PETSC_COMM_WORLD,&rank); 
MPI_Comm_size(PETSC_COMM_WORLD,&size); 

cout << "Rank = " << rank << endl; 
cout << "Size = " << size << endl; 

cout << "Generating 2D-Array" << endl; 

double temp2D[120000][3]; 
for (Ii=Istart; Ii<Iend; Ii++) { 
    for(J=0; J<n;J++){ 
     temp2D[Ii][J] = 1; 
    } 
    } 
    cout << "Processor " << rank << " set values : " << Istart << " - " << Iend << " into 2D-Array" << endl; 

    v = -1.0; 
    for (Ii=Istart; Ii<Iend; Ii++) { 
    for(J=0; J<n;J++){ 
     MatSetValues(A,1,&Ii,1,&J,&v,INSERT_VALUES);CHKERRQ(ierr); 
    } 
    } 
    cout << "Ii = " << Ii << " processor " << rank << " and it owns: " << Istart << " - " << Iend << endl; 

    /* 
    Assemble matrix, using the 2-step process: 
     MatAssemblyBegin(), MatAssemblyEnd() 
    Computations can be done while messages are in transition 
    by placing code between these two statements. 
    */ 
    ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 
    ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 

    MPI_Finalize(); 
cout << "No more MPI" << endl; 
return 0; 

} 

而我真正的程序有几个不同的.cpp文件只是简单的示例代码。我在主程序中初始化MPI，在另一个.cpp文件中调用了一个函数，我在那里实现了相同类型的矩阵填充，但是在填充矩阵之前程序所做的所有cout的打印次数将与我的核心数一样多。

我可以运行我的测试程序的程序mpiexec -n 4测试，它成功运行，但由于某种原因，我要运行我的真正的程序为我的测试程序的程序mpiexec -n 4 ./myprog

输出是继

Solving began... 
Solving began... 
Solving began... 
Solving began... 
Rank = 0 
Size = 4 
Generating 2D-Array 
Processor 0 set values : 0 - 30000 into 2D-Array 
Rank = 2 
Size = 4 
Generating 2D-Array 
Processor 2 set values : 60000 - 90000 into 2D-Array 
Rank = 3 
Size = 4 
Generating 2D-Array 
Processor 3 set values : 90000 - 120000 into 2D-Array 
Rank = 1 
Size = 4 
Generating 2D-Array 
Processor 1 set values : 30000 - 60000 into 2D-Array 
Ii = 30000 processor 0 and it owns: 0 - 30000 
Ii = 90000 processor 2 and it owns: 60000 - 90000 
Ii = 120000 processor 3 and it owns: 90000 - 120000 
Ii = 60000 processor 1 and it owns: 30000 - 60000 
no more MPI 
no more MPI 
no more MPI 
no more MPI 

编辑后的两点意见： 所以我的目标是拥有20个节点的小集群上运行这一点，每个节点有2个内核。后来这应该在超级计算机上运行，​​所以mpi绝对是我需要去的方式。我目前正在两台不同的机器上测试它，其中一台具有1个处理器/ 4个内核，另一个具有4个处理器/ 16个内核。

Answer 1

MPI是SPMD/MPMD模型（单程序多数据/多程序多数据）的实现。 MPI作业由同时运行的进程组成，它们互相交换消息以便合作解决问题。您不能只并行运行部分代码。你只能有部分代码不相互通信，但仍然可以同时执行。而你应该使用mpirun或mpiexec以并行模式启动你的应用程序。

如果您只想将代码的一部分并行化，并且可以忍受只能在单台机器上运行代码的限制，那么您需要的是OpenMP而不是MPI。或者，您也可以根据PETSc网站使用底层的POSIX线程编程，它支持pthreads。而且OpenMP建立在pthreads之上，所以使用PETSc和OpenMP可能是可能的。

Answer 2

为了增加Hristo的答案，MPI被构建为以分布式方式运行，即完全独立的进程。它们必须是分开的，因为它们应该在不同的物理机器上。您可以在一台机器上运行多个MPI进程，例如每个核心一个进程。这很好，但MPI没有任何工具来利用共享内存上下文。换句话说，您不能让某些MPI等级（进程）对由另一个MPI进程拥有的矩阵进行工作，因为您无法共享该矩阵。

当您启动x MPI流程时，您将获得x个运行相同的确切程序的副本。您需要的代码类似于

if (rank == 0) 
    do something 
else 
    do something else 

使不同的过程做不同的事情。进程可以通过发送消息相互通信，但它们都运行相同的确切二进制文件。 如果你没有代码分歧，那么你会得到同一个程序的x个副本，并给出相同的结果x次。