Actual source code: ex2.c
petsc-3.9.4 2018-09-11
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: */
13: int main(int argc,char **args)
14: {
15: PetscErrorCode ierr;
16: PetscCommShared scomm;
17: MPI_Comm comm;
18: PetscMPIInt lrank,rank,size,i;
19: Vec x,y;
20: VecScatter vscat;
21: IS isstride,isblock;
22: PetscViewer singleton;
23: PetscInt indices[] = {0,1,2};
25: PetscInitialize(&argc,&args,(char*)0,help);
26: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
27: MPI_Comm_size(PETSC_COMM_WORLD,&size);
28: if (size != 3) SETERRQ(PETSC_COMM_WORLD,1,"This example only works for 3 processes");
30: PetscCommDuplicate(PETSC_COMM_WORLD,&comm,NULL);
31: PetscCommSharedGet(comm,&scomm);
33: for (i=0; i<size; i++) {
34: PetscCommSharedGlobalToLocal(scomm,i,&lrank);
35: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d] Global rank %d shared memory comm rank %d\n",rank,i,lrank);
36: }
37: PetscSynchronizedFlush(PETSC_COMM_WORLD,stdout);
38: PetscCommDestroy(&comm);
40: VecCreateMPI(PETSC_COMM_WORLD,2,PETSC_DETERMINE,&x);
41: VecSetBlockSize(x,2);
42: VecSetValue(x,2*rank,(PetscScalar)(2*rank+10),INSERT_VALUES);
43: VecSetValue(x,2*rank+1,(PetscScalar)(2*rank+1+10),INSERT_VALUES);
44: VecAssemblyBegin(x);
45: VecAssemblyEnd(x);
46: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
48: VecCreateSeq(PETSC_COMM_SELF,6,&y);
49: VecSetBlockSize(y,2);
50: ISCreateStride(PETSC_COMM_SELF,6,0,1,&isstride);
51: ISCreateBlock(PETSC_COMM_SELF,2,3,indices,PETSC_COPY_VALUES,&isblock);
52: VecScatterCreate(x,isblock,y,isstride,&vscat);
53: VecScatterBegin(vscat,x,y,INSERT_VALUES,SCATTER_FORWARD);
54: VecScatterEnd(vscat,x,y,INSERT_VALUES,SCATTER_FORWARD);
55: VecScatterDestroy(&vscat);
56: PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&singleton);
57: VecView(y,singleton);
58: PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&singleton);
60: ISDestroy(&isstride);
61: ISDestroy(&isblock);
62: VecDestroy(&x);
63: VecDestroy(&y);
64: PetscFinalize();
65: return 0;
66: }