1 static char help[] = "Test MatCreateSubMatrices\n\n"; 2 3 #include <petscis.h> 4 #include <petscmat.h> 5 6 int main(int argc,char **args) 7 { 8 Mat A,*submats,*submats2; 9 IS *irow,*icol; 10 PetscInt i,n; 11 PetscMPIInt rank; 12 PetscViewer matfd,rowfd,colfd; 13 PetscBool same; 14 char matfile[PETSC_MAX_PATH_LEN],rowfile[PETSC_MAX_PATH_LEN],colfile[PETSC_MAX_PATH_LEN]; 15 char rankstr[16]={0}; 16 17 PetscCall(PetscInitialize(&argc,&args,(char*)0,help)); 18 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 19 20 PetscCall(PetscOptionsGetString(NULL,NULL,"-A",matfile,sizeof(matfile),NULL)); 21 PetscCall(PetscOptionsGetString(NULL,NULL,"-row",rowfile,sizeof(rowfile),NULL)); 22 PetscCall(PetscOptionsGetString(NULL,NULL,"-col",colfile,sizeof(colfile),NULL)); 23 24 /* Each rank has its own files for row/col ISes */ 25 PetscCall(PetscSNPrintf(rankstr,16,"-%d",rank)); 26 PetscCall(PetscStrlcat(rowfile,rankstr,PETSC_MAX_PATH_LEN)); 27 PetscCall(PetscStrlcat(colfile,rankstr,PETSC_MAX_PATH_LEN)); 28 29 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,matfile,FILE_MODE_READ,&matfd)); 30 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF,rowfile,FILE_MODE_READ,&rowfd)); 31 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF,colfile,FILE_MODE_READ,&colfd)); 32 33 PetscCall(MatCreate(PETSC_COMM_WORLD,&A)); 34 PetscCall(MatSetFromOptions(A)); 35 PetscCall(MatLoad(A,matfd)); 36 37 /* We stored the number of ISes at the beginning of rowfd */ 38 PetscCall(PetscViewerBinaryRead(rowfd,&n,1,NULL,PETSC_INT)); 39 PetscCall(PetscMalloc2(n,&irow,n,&icol)); 40 for (i=0; i<n; i++) { 41 PetscCall(ISCreate(PETSC_COMM_SELF,&irow[i])); 42 PetscCall(ISCreate(PETSC_COMM_SELF,&icol[i])); 43 PetscCall(ISLoad(irow[i],rowfd)); 44 PetscCall(ISLoad(icol[i],colfd)); 45 } 46 47 PetscCall(PetscViewerDestroy(&matfd)); 48 PetscCall(PetscViewerDestroy(&rowfd)); 49 PetscCall(PetscViewerDestroy(&colfd)); 50 51 /* Create submats for the first time */ 52 PetscCall(MatCreateSubMatrices(A,n,irow,icol,MAT_INITIAL_MATRIX,&submats)); 53 54 /* Dup submats to submats2 for later comparison */ 55 PetscCall(PetscMalloc1(n,&submats2)); 56 for (i=0; i<n; i++) { 57 PetscCall(MatDuplicate(submats[i],MAT_COPY_VALUES,&submats2[i])); 58 } 59 60 /* Create submats again */ 61 PetscCall(MatCreateSubMatrices(A,n,irow,icol,MAT_REUSE_MATRIX,&submats)); 62 63 /* Compare submats and submats2 */ 64 for (i=0; i<n; i++) { 65 PetscCall(MatEqual(submats[i],submats2[i],&same)); 66 PetscCheck(same,PETSC_COMM_SELF,PETSC_ERR_PLIB,"submatrix %" PetscInt_FMT " is not same",i); 67 } 68 69 PetscCall(MatDestroy(&A)); 70 for (i=0; i<n; i++) { 71 PetscCall(ISDestroy(&irow[i])); 72 PetscCall(ISDestroy(&icol[i])); 73 } 74 PetscCall(MatDestroySubMatrices(n,&submats)); 75 PetscCall(MatDestroyMatrices(n,&submats2)); 76 PetscCall(PetscFree2(irow,icol)); 77 PetscCall(PetscFinalize()); 78 return 0; 79 } 80 81 /*TEST 82 83 test: 84 suffix: 1 85 nsize: 2 86 requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES) 87 args: -mat_type {{aij baij}} -A ${DATAFILESPATH}/matrices/CreateSubMatrices/A -row ${DATAFILESPATH}/matrices/CreateSubMatrices/row -col ${DATAFILESPATH}/matrices/CreateSubMatrices/col 88 89 TEST*/ 90