Actual source code: ex2.c

petsc-3.8.4 2018-03-24
Report Typos and Errors
  1: static const char help[] = "Test overlapped communication on a single star forest (PetscSF)\n\n";

  3:  #include <petscvec.h>
  4:  #include <petscsf.h>
  5:  #include <petscviewer.h>

  7: int main(int argc, char **argv)
  8: {
  9:   PetscInt    ierr;
 10:   PetscSF     sf;
 11:   Vec         A,Aout;
 12:   Vec         B,Bout;
 13:   PetscScalar *bufA;
 14:   PetscScalar *bufAout;
 15:   PetscScalar *bufB;
 16:   PetscScalar *bufBout;
 17:   PetscMPIInt rank, size;
 18:   PetscInt    nroots, nleaves;
 19:   PetscInt    i;
 20:   PetscInt    *ilocal;
 21:   PetscSFNode *iremote;

 23:   PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr;
 24:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 25:   MPI_Comm_size(PETSC_COMM_WORLD,&size);

 27:   if (size != 2) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_USER, "Only coded for two MPI processes\n");

 29:   PetscSFCreate(PETSC_COMM_WORLD,&sf);
 30:   PetscSFSetFromOptions(sf);

 32:   nleaves = 2;
 33:   nroots = 1;
 34:   PetscMalloc1(nleaves,&ilocal);

 36:   for (i = 0; i<nleaves; i++) {
 37:     ilocal[i] = i;
 38:   }

 40:   PetscMalloc1(nleaves,&iremote);
 41:   if (rank == 0) {
 42:     iremote[0].rank = 0;
 43:     iremote[0].index = 0;
 44:     iremote[1].rank = 1;
 45:     iremote[1].index = 0;
 46:   } else {
 47:     iremote[0].rank = 1;
 48:     iremote[0].index = 0;
 49:     iremote[1].rank = 0;
 50:     iremote[1].index = 0;
 51:   }
 52:   PetscSFSetGraph(sf,nroots,nleaves,ilocal,PETSC_OWN_POINTER,iremote,PETSC_OWN_POINTER);
 53:   PetscSFSetUp(sf);
 54:   PetscSFView(sf,PETSC_VIEWER_STDOUT_WORLD);
 55:   VecCreate(PETSC_COMM_WORLD,&A);
 56:   VecSetSizes(A,2,PETSC_DETERMINE);
 57:   VecSetFromOptions(A);
 58:   VecSetUp(A);

 60:   VecDuplicate(A,&B);
 61:   VecDuplicate(A,&Aout);
 62:   VecDuplicate(A,&Bout);
 63:   VecGetArray(A,&bufA);
 64:   VecGetArray(B,&bufB);
 65:   for (i=0; i<2; i++) {
 66:     bufA[i] = (PetscScalar)rank;
 67:     bufB[i] = (PetscScalar)(rank) + 10.0;
 68:   }
 69:   VecRestoreArray(A,&bufA);
 70:   VecRestoreArray(B,&bufB);

 72:   VecGetArrayRead(A,(const PetscScalar**)&bufA);
 73:   VecGetArrayRead(B,(const PetscScalar**)&bufB);
 74:   VecGetArray(Aout,&bufAout);
 75:   VecGetArray(Bout,&bufBout);
 76:   PetscSFBcastBegin(sf,MPIU_SCALAR,(const void*)bufA,(void *)bufAout);
 77:   PetscSFBcastBegin(sf,MPIU_SCALAR,(const void*)bufB,(void *)bufBout);
 78:   PetscSFBcastEnd(sf,MPIU_SCALAR,(const void*)bufA,(void *)bufAout);
 79:   PetscSFBcastEnd(sf,MPIU_SCALAR,(const void*)bufB,(void *)bufBout);
 80:   VecRestoreArrayRead(A,(const PetscScalar**)&bufA);
 81:   VecRestoreArrayRead(B,(const PetscScalar**)&bufB);
 82:   VecRestoreArray(Aout,&bufAout);
 83:   VecRestoreArray(Bout,&bufBout);

 85:   VecView(Aout,PETSC_VIEWER_STDOUT_WORLD);
 86:   VecView(Bout,PETSC_VIEWER_STDOUT_WORLD);
 87:   VecDestroy(&A);
 88:   VecDestroy(&B);
 89:   VecDestroy(&Aout);
 90:   VecDestroy(&Bout);
 91:   PetscSFDestroy(&sf);

 93:   PetscFinalize();
 94:   return ierr;
 95: }