1c4762a1bSJed Brown 2c4762a1bSJed Brown static char help[] = "Tests MatIncreaseOverlap() - the parallel case. This example\n\ 3c4762a1bSJed Brown is similar to ex40.c; here the index sets used are random. Input arguments are:\n\ 4c4762a1bSJed Brown -f <input_file> : file to load. For example see $PETSC_DIR/share/petsc/datafiles/matrices\n\ 5c4762a1bSJed Brown -nd <size> : > 0 no of domains per processor \n\ 6c4762a1bSJed Brown -ov <overlap> : >=0 amount of overlap between domains\n\n"; 7c4762a1bSJed Brown 8c4762a1bSJed Brown #include <petscmat.h> 9c4762a1bSJed Brown 10c4762a1bSJed Brown int main(int argc,char **args) 11c4762a1bSJed Brown { 12c4762a1bSJed Brown PetscInt nd = 2,ov=1,i,j,m,n,*idx,lsize; 13c4762a1bSJed Brown PetscMPIInt rank; 14c4762a1bSJed Brown PetscBool flg; 15c4762a1bSJed Brown Mat A,B; 16c4762a1bSJed Brown char file[PETSC_MAX_PATH_LEN]; 17c4762a1bSJed Brown PetscViewer fd; 18c4762a1bSJed Brown IS *is1,*is2; 19c4762a1bSJed Brown PetscRandom r; 20c4762a1bSJed Brown PetscScalar rand; 21c4762a1bSJed Brown 22*b122ec5aSJacob Faibussowitsch CHKERRQ(PetscInitialize(&argc,&args,(char*)0,help)); 235f80ce2aSJacob Faibussowitsch CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 245f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),NULL)); 255f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-nd",&nd,NULL)); 265f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-ov",&ov,NULL)); 27c4762a1bSJed Brown 28c4762a1bSJed Brown /* Read matrix and RHS */ 295f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&fd)); 305f80ce2aSJacob Faibussowitsch CHKERRQ(MatCreate(PETSC_COMM_WORLD,&A)); 315f80ce2aSJacob Faibussowitsch CHKERRQ(MatSetType(A,MATMPIAIJ)); 325f80ce2aSJacob Faibussowitsch CHKERRQ(MatLoad(A,fd)); 335f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerDestroy(&fd)); 34c4762a1bSJed Brown 35c4762a1bSJed Brown /* Read the matrix again as a seq matrix */ 365f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerBinaryOpen(PETSC_COMM_SELF,file,FILE_MODE_READ,&fd)); 375f80ce2aSJacob Faibussowitsch CHKERRQ(MatCreate(PETSC_COMM_SELF,&B)); 385f80ce2aSJacob Faibussowitsch CHKERRQ(MatSetType(B,MATSEQAIJ)); 395f80ce2aSJacob Faibussowitsch CHKERRQ(MatLoad(B,fd)); 405f80ce2aSJacob Faibussowitsch CHKERRQ(PetscViewerDestroy(&fd)); 41c4762a1bSJed Brown 42c4762a1bSJed Brown /* Create the Random no generator */ 435f80ce2aSJacob Faibussowitsch CHKERRQ(MatGetSize(A,&m,&n)); 445f80ce2aSJacob Faibussowitsch CHKERRQ(PetscRandomCreate(PETSC_COMM_SELF,&r)); 455f80ce2aSJacob Faibussowitsch CHKERRQ(PetscRandomSetFromOptions(r)); 46c4762a1bSJed Brown 47c4762a1bSJed Brown /* Create the IS corresponding to subdomains */ 485f80ce2aSJacob Faibussowitsch CHKERRQ(PetscMalloc1(nd,&is1)); 495f80ce2aSJacob Faibussowitsch CHKERRQ(PetscMalloc1(nd,&is2)); 505f80ce2aSJacob Faibussowitsch CHKERRQ(PetscMalloc1(m ,&idx)); 51c4762a1bSJed Brown 52c4762a1bSJed Brown /* Create the random Index Sets */ 53c4762a1bSJed Brown for (i=0; i<nd; i++) { 54c4762a1bSJed Brown for (j=0; j<rank; j++) { 555f80ce2aSJacob Faibussowitsch CHKERRQ(PetscRandomGetValue(r,&rand)); 56c4762a1bSJed Brown } 575f80ce2aSJacob Faibussowitsch CHKERRQ(PetscRandomGetValue(r,&rand)); 58c4762a1bSJed Brown lsize = (PetscInt)(rand*m); 59c4762a1bSJed Brown for (j=0; j<lsize; j++) { 605f80ce2aSJacob Faibussowitsch CHKERRQ(PetscRandomGetValue(r,&rand)); 61c4762a1bSJed Brown idx[j] = (PetscInt)(rand*m); 62c4762a1bSJed Brown } 635f80ce2aSJacob Faibussowitsch CHKERRQ(ISCreateGeneral(PETSC_COMM_SELF,lsize,idx,PETSC_COPY_VALUES,is1+i)); 645f80ce2aSJacob Faibussowitsch CHKERRQ(ISCreateGeneral(PETSC_COMM_SELF,lsize,idx,PETSC_COPY_VALUES,is2+i)); 65c4762a1bSJed Brown } 66c4762a1bSJed Brown 675f80ce2aSJacob Faibussowitsch CHKERRQ(MatIncreaseOverlap(A,nd,is1,ov)); 685f80ce2aSJacob Faibussowitsch CHKERRQ(MatIncreaseOverlap(B,nd,is2,ov)); 69c4762a1bSJed Brown 70c4762a1bSJed Brown /* Now see if the serial and parallel case have the same answers */ 71c4762a1bSJed Brown for (i=0; i<nd; ++i) { 72c4762a1bSJed Brown PetscInt sz1,sz2; 735f80ce2aSJacob Faibussowitsch CHKERRQ(ISEqual(is1[i],is2[i],&flg)); 745f80ce2aSJacob Faibussowitsch CHKERRQ(ISGetSize(is1[i],&sz1)); 755f80ce2aSJacob Faibussowitsch CHKERRQ(ISGetSize(is2[i],&sz2)); 7628b400f6SJacob Faibussowitsch PetscCheck(flg,PETSC_COMM_SELF,PETSC_ERR_PLIB,"proc:[%d], i=%" PetscInt_FMT ", flg =%d sz1 = %" PetscInt_FMT " sz2 = %" PetscInt_FMT,rank,i,(int)flg,sz1,sz2); 77c4762a1bSJed Brown } 78c4762a1bSJed Brown 79c4762a1bSJed Brown /* Free Allocated Memory */ 80c4762a1bSJed Brown for (i=0; i<nd; ++i) { 815f80ce2aSJacob Faibussowitsch CHKERRQ(ISDestroy(&is1[i])); 825f80ce2aSJacob Faibussowitsch CHKERRQ(ISDestroy(&is2[i])); 83c4762a1bSJed Brown } 845f80ce2aSJacob Faibussowitsch CHKERRQ(PetscRandomDestroy(&r)); 855f80ce2aSJacob Faibussowitsch CHKERRQ(PetscFree(is1)); 865f80ce2aSJacob Faibussowitsch CHKERRQ(PetscFree(is2)); 875f80ce2aSJacob Faibussowitsch CHKERRQ(MatDestroy(&A)); 885f80ce2aSJacob Faibussowitsch CHKERRQ(MatDestroy(&B)); 895f80ce2aSJacob Faibussowitsch CHKERRQ(PetscFree(idx)); 90*b122ec5aSJacob Faibussowitsch CHKERRQ(PetscFinalize()); 91*b122ec5aSJacob Faibussowitsch return 0; 92c4762a1bSJed Brown } 93c4762a1bSJed Brown 94c4762a1bSJed Brown /*TEST 95c4762a1bSJed Brown 96c4762a1bSJed Brown build: 97c4762a1bSJed Brown requires: !complex 98c4762a1bSJed Brown 99c4762a1bSJed Brown test: 100c4762a1bSJed Brown nsize: 3 101dfd57a17SPierre Jolivet requires: datafilespath double !defined(PETSC_USE_64BIT_INDICES) !complex 102c4762a1bSJed Brown args: -f ${DATAFILESPATH}/matrices/arco1 -nd 3 -ov 1 103c4762a1bSJed Brown 104c4762a1bSJed Brown TEST*/ 105