1*97929ea7SJunchao Zhang static char help[]= "Test vecscatter of different block sizes across processes\n\n"; 2*97929ea7SJunchao Zhang 3*97929ea7SJunchao Zhang #include <petscvec.h> 4*97929ea7SJunchao Zhang int main(int argc,char **argv) 5*97929ea7SJunchao Zhang { 6*97929ea7SJunchao Zhang PetscErrorCode ierr; 7*97929ea7SJunchao Zhang PetscInt i,bs,n,low,high; 8*97929ea7SJunchao Zhang PetscMPIInt nproc,rank; 9*97929ea7SJunchao Zhang Vec x,y,z; 10*97929ea7SJunchao Zhang IS ix,iy; 11*97929ea7SJunchao Zhang VecScatter vscat; 12*97929ea7SJunchao Zhang const PetscScalar *yv; 13*97929ea7SJunchao Zhang 14*97929ea7SJunchao Zhang ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr; 15*97929ea7SJunchao Zhang ierr = MPI_Comm_size(PETSC_COMM_WORLD,&nproc);CHKERRQ(ierr); 16*97929ea7SJunchao Zhang ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRQ(ierr); 17*97929ea7SJunchao Zhang if (nproc != 2) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"This test can only run on two MPI ranks"); 18*97929ea7SJunchao Zhang 19*97929ea7SJunchao Zhang /* Create an MPI vector x of size 12 on two processes, and set x = {0, 1, 2, .., 11} */ 20*97929ea7SJunchao Zhang ierr = VecCreateMPI(PETSC_COMM_WORLD,6,PETSC_DECIDE,&x);CHKERRQ(ierr); 21*97929ea7SJunchao Zhang ierr = VecGetOwnershipRange(x,&low,&high);CHKERRQ(ierr); 22*97929ea7SJunchao Zhang for (i=low; i<high; i++) {ierr = VecSetValue(x,i,(PetscScalar)i,INSERT_VALUES);CHKERRQ(ierr);} 23*97929ea7SJunchao Zhang ierr = VecAssemblyBegin(x);CHKERRQ(ierr); 24*97929ea7SJunchao Zhang ierr = VecAssemblyEnd(x);CHKERRQ(ierr); 25*97929ea7SJunchao Zhang 26*97929ea7SJunchao Zhang /* Create a seq vector y, and a parallel to sequential (PtoS) vecscatter to scatter x to y */ 27*97929ea7SJunchao Zhang if (!rank) { 28*97929ea7SJunchao Zhang /* On rank 0, seq y is of size 6. We will scatter x[0,1,2,6,7,8] to y[0,1,2,3,4,5] using IS with bs=3 */ 29*97929ea7SJunchao Zhang PetscInt idx[2]={0,2}; 30*97929ea7SJunchao Zhang PetscInt idy[2]={0,1}; 31*97929ea7SJunchao Zhang n = 6; 32*97929ea7SJunchao Zhang bs = 3; 33*97929ea7SJunchao Zhang ierr = VecCreateSeq(PETSC_COMM_SELF,n,&y);CHKERRQ(ierr); 34*97929ea7SJunchao Zhang ierr = ISCreateBlock(PETSC_COMM_SELF,bs,2,idx,PETSC_COPY_VALUES,&ix);CHKERRQ(ierr); 35*97929ea7SJunchao Zhang ierr = ISCreateBlock(PETSC_COMM_SELF,bs,2,idy,PETSC_COPY_VALUES,&iy);CHKERRQ(ierr); 36*97929ea7SJunchao Zhang } else { 37*97929ea7SJunchao Zhang /* On rank 1, seq y is of size 4. We will scatter x[4,5,10,11] to y[0,1,2,3] using IS with bs=2 */ 38*97929ea7SJunchao Zhang PetscInt idx[2]= {2,5}; 39*97929ea7SJunchao Zhang PetscInt idy[2]= {0,1}; 40*97929ea7SJunchao Zhang n = 4; 41*97929ea7SJunchao Zhang bs = 2; 42*97929ea7SJunchao Zhang ierr = VecCreateSeq(PETSC_COMM_SELF,n,&y);CHKERRQ(ierr); 43*97929ea7SJunchao Zhang ierr = ISCreateBlock(PETSC_COMM_SELF,bs,2,idx,PETSC_COPY_VALUES,&ix);CHKERRQ(ierr); 44*97929ea7SJunchao Zhang ierr = ISCreateBlock(PETSC_COMM_SELF,bs,2,idy,PETSC_COPY_VALUES,&iy);CHKERRQ(ierr); 45*97929ea7SJunchao Zhang } 46*97929ea7SJunchao Zhang ierr = VecScatterCreate(x,ix,y,iy,&vscat);CHKERRQ(ierr); 47*97929ea7SJunchao Zhang 48*97929ea7SJunchao Zhang /* Do the vecscatter */ 49*97929ea7SJunchao Zhang ierr = VecScatterBegin(vscat,x,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 50*97929ea7SJunchao Zhang ierr = VecScatterEnd(vscat,x,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 51*97929ea7SJunchao Zhang 52*97929ea7SJunchao Zhang /* Print y. Since y is sequential, we put y in a parallel z to print its value on both ranks */ 53*97929ea7SJunchao Zhang ierr = VecGetArrayRead(y,&yv);CHKERRQ(ierr); 54*97929ea7SJunchao Zhang ierr = VecCreateMPIWithArray(PETSC_COMM_WORLD,1,n,PETSC_DECIDE,yv,&z);CHKERRQ(ierr); 55*97929ea7SJunchao Zhang ierr = VecView(z,PETSC_VIEWER_STDOUT_WORLD);CHKERRQ(ierr); 56*97929ea7SJunchao Zhang ierr = VecRestoreArrayRead(y,&yv);CHKERRQ(ierr); 57*97929ea7SJunchao Zhang 58*97929ea7SJunchao Zhang ierr = ISDestroy(&ix);CHKERRQ(ierr); 59*97929ea7SJunchao Zhang ierr = ISDestroy(&iy);CHKERRQ(ierr); 60*97929ea7SJunchao Zhang ierr = VecDestroy(&x);CHKERRQ(ierr); 61*97929ea7SJunchao Zhang ierr = VecDestroy(&y);CHKERRQ(ierr); 62*97929ea7SJunchao Zhang ierr = VecDestroy(&z);CHKERRQ(ierr); 63*97929ea7SJunchao Zhang ierr = VecScatterDestroy(&vscat);CHKERRQ(ierr); 64*97929ea7SJunchao Zhang 65*97929ea7SJunchao Zhang ierr = PetscFinalize(); 66*97929ea7SJunchao Zhang return ierr; 67*97929ea7SJunchao Zhang } 68*97929ea7SJunchao Zhang 69*97929ea7SJunchao Zhang /*TEST 70*97929ea7SJunchao Zhang 71*97929ea7SJunchao Zhang test: 72*97929ea7SJunchao Zhang nsize: 2 73*97929ea7SJunchao Zhang args: 74*97929ea7SJunchao Zhang requires: 75*97929ea7SJunchao Zhang TEST*/ 76*97929ea7SJunchao Zhang 77