1 static const char help[] = "Test PetscSF with MPI large count (more than 2 billion elements in messages)\n\n"; 2 3 #include <petscsys.h> 4 #include <petscsf.h> 5 6 int main(int argc,char **argv) 7 { 8 PetscErrorCode ierr; 9 PetscSF sf; 10 PetscInt i,nroots,nleaves; 11 PetscInt n = (1ULL<<31) + 1024; /* a little over 2G elements */ 12 PetscSFNode *iremote = NULL; 13 PetscMPIInt rank,size; 14 char *rootdata=NULL,*leafdata=NULL; 15 Vec x,y; 16 VecScatter vscat; 17 PetscInt rstart,rend; 18 IS ix; 19 const PetscScalar *xv; 20 21 ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 22 ierr = MPI_Comm_rank(PETSC_COMM_WORLD,&rank);CHKERRMPI(ierr); 23 ierr = MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRMPI(ierr); 24 PetscCheckFalse(size != 2,PETSC_COMM_WORLD,PETSC_ERR_WRONG_MPI_SIZE,"The test can only run with two MPI ranks"); 25 26 /* Test PetscSF */ 27 ierr = PetscSFCreate(PETSC_COMM_WORLD,&sf);CHKERRQ(ierr); 28 ierr = PetscSFSetFromOptions(sf);CHKERRQ(ierr); 29 30 if (!rank) { 31 nroots = n; 32 nleaves = 0; 33 } else { 34 nroots = 0; 35 nleaves = n; 36 ierr = PetscMalloc1(nleaves,&iremote);CHKERRQ(ierr); 37 for (i=0; i<nleaves; i++) { 38 iremote[i].rank = 0; 39 iremote[i].index = i; 40 } 41 } 42 ierr = PetscSFSetGraph(sf,nroots,nleaves,NULL,PETSC_COPY_VALUES,iremote,PETSC_COPY_VALUES);CHKERRQ(ierr); 43 ierr = PetscMalloc2(nroots,&rootdata,nleaves,&leafdata);CHKERRQ(ierr); 44 if (!rank) { 45 memset(rootdata,11,nroots); 46 rootdata[nroots-1] = 12; /* Use a different value at the end */ 47 } 48 49 ierr = PetscSFBcastBegin(sf,MPI_SIGNED_CHAR,rootdata,leafdata,MPI_REPLACE);CHKERRQ(ierr); /* rank 0->1, bcast rootdata to leafdata */ 50 ierr = PetscSFBcastEnd(sf,MPI_SIGNED_CHAR,rootdata,leafdata,MPI_REPLACE);CHKERRQ(ierr); 51 ierr = PetscSFReduceBegin(sf,MPI_SIGNED_CHAR,leafdata,rootdata,MPI_SUM);CHKERRQ(ierr); /* rank 1->0, add leafdata to rootdata */ 52 ierr = PetscSFReduceEnd(sf,MPI_SIGNED_CHAR,leafdata,rootdata,MPI_SUM);CHKERRQ(ierr); 53 if (!rank) { 54 PetscCheckFalse(rootdata[0] != 22 || rootdata[nroots-1] != 24,PETSC_COMM_SELF,PETSC_ERR_PLIB,"SF: wrong results"); 55 } 56 57 ierr = PetscFree2(rootdata,leafdata);CHKERRQ(ierr); 58 ierr = PetscFree(iremote);CHKERRQ(ierr); 59 ierr = PetscSFDestroy(&sf);CHKERRQ(ierr); 60 61 /* Test VecScatter */ 62 ierr = VecCreate(PETSC_COMM_WORLD,&x);CHKERRQ(ierr); 63 ierr = VecCreate(PETSC_COMM_WORLD,&y);CHKERRQ(ierr); 64 ierr = VecSetSizes(x,rank==0? n : 64,PETSC_DECIDE);CHKERRQ(ierr); 65 ierr = VecSetSizes(y,rank==0? 64 : n,PETSC_DECIDE);CHKERRQ(ierr); 66 ierr = VecSetFromOptions(x);CHKERRQ(ierr); 67 ierr = VecSetFromOptions(y);CHKERRQ(ierr); 68 69 ierr = VecGetOwnershipRange(x,&rstart,&rend);CHKERRQ(ierr); 70 ierr = ISCreateStride(PETSC_COMM_SELF,rend-rstart,rstart,1,&ix);CHKERRQ(ierr); 71 ierr = VecScatterCreate(x,ix,y,ix,&vscat);CHKERRQ(ierr); 72 73 ierr = VecSet(x,3.0);CHKERRQ(ierr); 74 ierr = VecScatterBegin(vscat,x,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 75 ierr = VecScatterEnd(vscat,x,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 76 77 ierr = VecScatterBegin(vscat,y,x,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 78 ierr = VecScatterEnd(vscat,y,x,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 79 80 ierr = VecGetArrayRead(x,&xv);CHKERRQ(ierr); 81 PetscCheckFalse(xv[0] != 6.0,PETSC_COMM_SELF,PETSC_ERR_PLIB,"VecScatter: wrong results"); 82 ierr = VecRestoreArrayRead(x,&xv);CHKERRQ(ierr); 83 ierr = VecDestroy(&x);CHKERRQ(ierr); 84 ierr = VecDestroy(&y);CHKERRQ(ierr); 85 ierr = VecScatterDestroy(&vscat);CHKERRQ(ierr); 86 ierr = ISDestroy(&ix);CHKERRQ(ierr); 87 88 ierr = PetscFinalize(); 89 return ierr; 90 } 91 92 /**TEST 93 test: 94 requires: defined(PETSC_HAVE_MPI_LARGE_COUNT) defined(PETSC_USE_64BIT_INDICES) 95 TODO: need a machine with big memory (~150GB) to run the test 96 nsize: 2 97 args: -sf_type {{basic neighbor}} 98 99 TEST**/ 100 101