1 2 static char help[] = "Tests MatCreateSubMatrices() for SBAIJ matrices\n\n"; 3 4 #include <petscmat.h> 5 6 int main(int argc,char **args) 7 { 8 Mat BAIJ,SBAIJ,*subBAIJ,*subSBAIJ; 9 PetscViewer viewer; 10 char file[PETSC_MAX_PATH_LEN]; 11 PetscBool flg; 12 PetscErrorCode ierr; 13 PetscInt n = 2,issize,M,N; 14 PetscMPIInt rank; 15 IS isrow,iscol,irow[n],icol[n]; 16 17 ierr = PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr; 18 CHKERRQ(PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&flg)); 19 CHKERRQ(PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&viewer)); 20 21 CHKERRQ(MatCreate(PETSC_COMM_WORLD,&BAIJ)); 22 CHKERRQ(MatSetType(BAIJ,MATMPIBAIJ)); 23 CHKERRQ(MatLoad(BAIJ,viewer)); 24 CHKERRQ(PetscViewerDestroy(&viewer)); 25 26 CHKERRQ(PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&viewer)); 27 CHKERRQ(MatCreate(PETSC_COMM_WORLD,&SBAIJ)); 28 CHKERRQ(MatSetType(SBAIJ,MATMPISBAIJ)); 29 CHKERRQ(MatLoad(SBAIJ,viewer)); 30 CHKERRQ(PetscViewerDestroy(&viewer)); 31 32 CHKERRQ(MatGetSize(BAIJ,&M,&N)); 33 issize = M/4; 34 CHKERRQ(ISCreateStride(PETSC_COMM_SELF,issize,0,1,&isrow)); 35 irow[0] = irow[1] = isrow; 36 issize = N/2; 37 CHKERRQ(ISCreateStride(PETSC_COMM_SELF,issize,0,1,&iscol)); 38 icol[0] = icol[1] = iscol; 39 CHKERRQ(MatCreateSubMatrices(BAIJ,n,irow,icol,MAT_INITIAL_MATRIX,&subBAIJ)); 40 CHKERRQ(MatCreateSubMatrices(BAIJ,n,irow,icol,MAT_REUSE_MATRIX,&subBAIJ)); 41 42 /* irow and icol must be same for SBAIJ matrices! */ 43 icol[0] = icol[1] = isrow; 44 CHKERRQ(MatCreateSubMatrices(SBAIJ,n,irow,icol,MAT_INITIAL_MATRIX,&subSBAIJ)); 45 CHKERRQ(MatCreateSubMatrices(SBAIJ,n,irow,icol,MAT_REUSE_MATRIX,&subSBAIJ)); 46 47 CHKERRMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 48 if (rank == 0) { 49 CHKERRQ(MatView(subBAIJ[0],PETSC_VIEWER_STDOUT_SELF)); 50 CHKERRQ(MatView(subSBAIJ[0],PETSC_VIEWER_STDOUT_SELF)); 51 } 52 53 /* Free data structures */ 54 CHKERRQ(ISDestroy(&isrow)); 55 CHKERRQ(ISDestroy(&iscol)); 56 CHKERRQ(MatDestroySubMatrices(n,&subBAIJ)); 57 CHKERRQ(MatDestroySubMatrices(n,&subSBAIJ)); 58 CHKERRQ(MatDestroy(&BAIJ)); 59 CHKERRQ(MatDestroy(&SBAIJ)); 60 61 ierr = PetscFinalize(); 62 return ierr; 63 } 64