Actual source code: ex2.c


  2: static char help[] = "Tests shared memory subcommunicators\n\n";
  3: #include <petscsys.h>
  4: #include <petscvec.h>

  6: /*
  7:    One can use petscmpiexec -n 3 -hosts localhost,Barrys-MacBook-Pro.local ./ex2 -info to mimic
  8:   having two nodes that do not share common memory
  9: */

 11: int main(int argc, char **args)
 12: {
 13:   PetscCommShared scomm;
 14:   MPI_Comm        comm;
 15:   PetscMPIInt     lrank, rank, size, i;
 16:   Vec             x, y;
 17:   VecScatter      vscat;
 18:   IS              isstride, isblock;
 19:   PetscViewer     singleton;
 20:   PetscInt        indices[] = {0, 1, 2};

 22:   PetscInitialize(&argc, &args, (char *)0, help);
 23:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
 24:   MPI_Comm_size(PETSC_COMM_WORLD, &size);

 27:   PetscCommDuplicate(PETSC_COMM_WORLD, &comm, NULL);
 28:   PetscCommSharedGet(comm, &scomm);

 30:   for (i = 0; i < size; i++) {
 31:     PetscCommSharedGlobalToLocal(scomm, i, &lrank);
 32:     PetscSynchronizedPrintf(PETSC_COMM_WORLD, "[%d] Global rank %d shared memory comm rank %d\n", rank, i, lrank);
 33:   }
 34:   PetscSynchronizedFlush(PETSC_COMM_WORLD, stdout);
 35:   PetscCommDestroy(&comm);

 37:   VecCreateMPI(PETSC_COMM_WORLD, 2, PETSC_DETERMINE, &x);
 38:   VecSetBlockSize(x, 2);
 39:   VecSetValue(x, 2 * rank, (PetscScalar)(2 * rank + 10), INSERT_VALUES);
 40:   VecSetValue(x, 2 * rank + 1, (PetscScalar)(2 * rank + 1 + 10), INSERT_VALUES);
 41:   VecAssemblyBegin(x);
 42:   VecAssemblyEnd(x);
 43:   VecView(x, PETSC_VIEWER_STDOUT_WORLD);

 45:   VecCreateSeq(PETSC_COMM_SELF, 6, &y);
 46:   VecSetBlockSize(y, 2);
 47:   ISCreateStride(PETSC_COMM_SELF, 6, 0, 1, &isstride);
 48:   ISCreateBlock(PETSC_COMM_SELF, 2, 3, indices, PETSC_COPY_VALUES, &isblock);
 49:   VecScatterCreate(x, isblock, y, isstride, &vscat);
 50:   VecScatterBegin(vscat, x, y, INSERT_VALUES, SCATTER_FORWARD);
 51:   VecScatterEnd(vscat, x, y, INSERT_VALUES, SCATTER_FORWARD);
 52:   VecScatterDestroy(&vscat);
 53:   PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &singleton);
 54:   VecView(y, singleton);
 55:   PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &singleton);

 57:   ISDestroy(&isstride);
 58:   ISDestroy(&isblock);
 59:   VecDestroy(&x);
 60:   VecDestroy(&y);
 61:   PetscFinalize();
 62:   return 0;
 63: }