1*97929ea7SJunchao Zhang static char help[]= " Test VecScatter with x, y on different communicators\n\n"; 2*97929ea7SJunchao Zhang 3*97929ea7SJunchao Zhang #include <petscvec.h> 4*97929ea7SJunchao Zhang 5*97929ea7SJunchao Zhang int main(int argc,char **argv) 6*97929ea7SJunchao Zhang { 7*97929ea7SJunchao Zhang PetscErrorCode ierr; 8*97929ea7SJunchao Zhang PetscInt i,n=5,rstart; 9*97929ea7SJunchao Zhang PetscScalar *val; 10*97929ea7SJunchao Zhang const PetscScalar *dat; 11*97929ea7SJunchao Zhang Vec x,y1,y2; 12*97929ea7SJunchao Zhang MPI_Comm newcomm; 13*97929ea7SJunchao Zhang PetscMPIInt nproc,rank; 14*97929ea7SJunchao Zhang IS ix; 15*97929ea7SJunchao Zhang VecScatter vscat1,vscat2; 16*97929ea7SJunchao Zhang 17*97929ea7SJunchao Zhang PetscFunctionBegin; 18*97929ea7SJunchao Zhang ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr; 19*97929ea7SJunchao Zhang ierr = MPI_Comm_size(PETSC_COMM_WORLD,&nproc);CHKERRQ(ierr); 20*97929ea7SJunchao Zhang ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr); 21*97929ea7SJunchao Zhang 22*97929ea7SJunchao Zhang if (nproc != 2) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"This test must run with exactly two MPI ranks\n"); 23*97929ea7SJunchao Zhang 24*97929ea7SJunchao Zhang /* Create MPI vectors x and y, which are on the same comm (i.e., MPI_IDENT) */ 25*97929ea7SJunchao Zhang ierr = VecCreateMPI(PETSC_COMM_WORLD,n,PETSC_DECIDE,&x);CHKERRQ(ierr); 26*97929ea7SJunchao Zhang ierr = VecDuplicate(x,&y1);CHKERRQ(ierr); 27*97929ea7SJunchao Zhang ierr = VecGetOwnershipRange(x,&rstart,NULL);CHKERRQ(ierr); 28*97929ea7SJunchao Zhang 29*97929ea7SJunchao Zhang /* Set x's value locally. x would be {0., 1., 2., ..., 9.} */ 30*97929ea7SJunchao Zhang ierr = VecGetArray(x,&val);CHKERRQ(ierr); 31*97929ea7SJunchao Zhang for (i=0; i<n; i++) val[i] = rstart + i; 32*97929ea7SJunchao Zhang ierr = VecRestoreArray(x,&val);CHKERRQ(ierr); 33*97929ea7SJunchao Zhang 34*97929ea7SJunchao Zhang /* Create index set ix = {0, 1, 2, ..., 9} */ 35*97929ea7SJunchao Zhang ierr = ISCreateStride(PETSC_COMM_WORLD,n,rstart,1,&ix);CHKERRQ(ierr); 36*97929ea7SJunchao Zhang 37*97929ea7SJunchao Zhang /* Create newcomm that reverses processes in x's comm, and then create y2 on it*/ 38*97929ea7SJunchao Zhang ierr = MPI_Comm_split(PETSC_COMM_WORLD,0/*color*/,-rank/*key*/,&newcomm);CHKERRQ(ierr); /* use -rank as key to reverse processes in newcomm */ 39*97929ea7SJunchao Zhang ierr = VecCreateMPI(newcomm,n,PETSC_DECIDE,&y2);CHKERRQ(ierr); 40*97929ea7SJunchao Zhang 41*97929ea7SJunchao Zhang /* It looks vscat1/2 are the same, but actually not. y2 is on a different communicator than x */ 42*97929ea7SJunchao Zhang ierr = VecScatterCreate(x,ix,y1,ix,&vscat1);CHKERRQ(ierr); 43*97929ea7SJunchao Zhang ierr = VecScatterCreate(x,ix,y2,ix,&vscat2);CHKERRQ(ierr); 44*97929ea7SJunchao Zhang 45*97929ea7SJunchao Zhang ierr = VecScatterBegin(vscat1,x,y1,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 46*97929ea7SJunchao Zhang ierr = VecScatterBegin(vscat2,x,y2,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 47*97929ea7SJunchao Zhang ierr = VecScatterEnd (vscat1,x,y1,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 48*97929ea7SJunchao Zhang ierr = VecScatterEnd (vscat2,x,y2,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 49*97929ea7SJunchao Zhang 50*97929ea7SJunchao Zhang /* View on rank 0 of x's comm, which is PETSC_COMM_WORLD */ 51*97929ea7SJunchao Zhang if (!rank) { 52*97929ea7SJunchao Zhang /* Print the part of x on rank 0, which is 0 1 2 3 4 */ 53*97929ea7SJunchao Zhang ierr = PetscPrintf(PETSC_COMM_SELF,"\nOn rank 0 of PETSC_COMM_WORLD, x = ");CHKERRQ(ierr); 54*97929ea7SJunchao Zhang ierr = VecGetArrayRead(x,&dat);CHKERRQ(ierr); 55*97929ea7SJunchao Zhang for (i=0; i<n; i++) {ierr = PetscPrintf(PETSC_COMM_SELF," %.0f",PetscRealPart(dat[i]));CHKERRQ(ierr);} 56*97929ea7SJunchao Zhang ierr = VecRestoreArrayRead(x,&dat);CHKERRQ(ierr); 57*97929ea7SJunchao Zhang 58*97929ea7SJunchao Zhang /* Print the part of y1 on rank 0, which is 0 1 2 3 4 */ 59*97929ea7SJunchao Zhang ierr = PetscPrintf(PETSC_COMM_SELF,"\nOn rank 0 of PETSC_COMM_WORLD, y1 = ");CHKERRQ(ierr); 60*97929ea7SJunchao Zhang ierr = VecGetArrayRead(y1,&dat);CHKERRQ(ierr); 61*97929ea7SJunchao Zhang for (i=0; i<n; i++) {ierr = PetscPrintf(PETSC_COMM_SELF," %.0f",PetscRealPart(dat[i]));CHKERRQ(ierr);} 62*97929ea7SJunchao Zhang ierr = VecRestoreArrayRead(y1,&dat);CHKERRQ(ierr); 63*97929ea7SJunchao Zhang 64*97929ea7SJunchao Zhang /* Print the part of y2 on rank 0, which is 5 6 7 8 9 since y2 swapped the processes of PETSC_COMM_WORLD */ 65*97929ea7SJunchao Zhang ierr = PetscPrintf(PETSC_COMM_SELF,"\nOn rank 0 of PETSC_COMM_WORLD, y2 = ");CHKERRQ(ierr); 66*97929ea7SJunchao Zhang ierr = VecGetArrayRead(y2,&dat);CHKERRQ(ierr); 67*97929ea7SJunchao Zhang for (i=0; i<n; i++) {ierr = PetscPrintf(PETSC_COMM_SELF," %.0f",PetscRealPart(dat[i]));CHKERRQ(ierr);} 68*97929ea7SJunchao Zhang ierr = VecRestoreArrayRead(y2,&dat);CHKERRQ(ierr); 69*97929ea7SJunchao Zhang ierr = PetscPrintf(PETSC_COMM_SELF,"\n");CHKERRQ(ierr); 70*97929ea7SJunchao Zhang } 71*97929ea7SJunchao Zhang 72*97929ea7SJunchao Zhang ierr = ISDestroy(&ix);CHKERRQ(ierr); 73*97929ea7SJunchao Zhang ierr = VecDestroy(&x);CHKERRQ(ierr); 74*97929ea7SJunchao Zhang ierr = VecDestroy(&y1);CHKERRQ(ierr); 75*97929ea7SJunchao Zhang ierr = VecDestroy(&y2);CHKERRQ(ierr); 76*97929ea7SJunchao Zhang ierr = VecScatterDestroy(&vscat1);CHKERRQ(ierr); 77*97929ea7SJunchao Zhang ierr = VecScatterDestroy(&vscat2);CHKERRQ(ierr); 78*97929ea7SJunchao Zhang ierr = MPI_Comm_free(&newcomm);CHKERRQ(ierr); 79*97929ea7SJunchao Zhang ierr = PetscFinalize(); 80*97929ea7SJunchao Zhang return ierr; 81*97929ea7SJunchao Zhang } 82*97929ea7SJunchao Zhang 83*97929ea7SJunchao Zhang /*TEST 84*97929ea7SJunchao Zhang 85*97929ea7SJunchao Zhang test: 86*97929ea7SJunchao Zhang nsize: 2 87*97929ea7SJunchao Zhang requires: double 88*97929ea7SJunchao Zhang TEST*/ 89*97929ea7SJunchao Zhang 90