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