Actual source code: ex87.c


  2: static char help[] = "Tests MatCreateSubMatrices() for SBAIJ matrices\n\n";

  4: #include <petscmat.h>

  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:   PetscInt       n = 2,issize,M,N;
 13:   PetscMPIInt    rank;
 14:   IS             isrow,iscol,irow[n],icol[n];

 16:   PetscInitialize(&argc,&args,(char*)0,help);
 17:   PetscOptionsGetString(NULL,NULL,"-f",file,sizeof(file),&flg);
 18:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&viewer);

 20:   MatCreate(PETSC_COMM_WORLD,&BAIJ);
 21:   MatSetType(BAIJ,MATMPIBAIJ);
 22:   MatLoad(BAIJ,viewer);
 23:   PetscViewerDestroy(&viewer);

 25:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,file,FILE_MODE_READ,&viewer);
 26:   MatCreate(PETSC_COMM_WORLD,&SBAIJ);
 27:   MatSetType(SBAIJ,MATMPISBAIJ);
 28:   MatLoad(SBAIJ,viewer);
 29:   PetscViewerDestroy(&viewer);

 31:   MatGetSize(BAIJ,&M,&N);
 32:   issize = M/4;
 33:   ISCreateStride(PETSC_COMM_SELF,issize,0,1,&isrow);
 34:   irow[0] = irow[1] = isrow;
 35:   issize = N/2;
 36:   ISCreateStride(PETSC_COMM_SELF,issize,0,1,&iscol);
 37:   icol[0] = icol[1] = iscol;
 38:   MatCreateSubMatrices(BAIJ,n,irow,icol,MAT_INITIAL_MATRIX,&subBAIJ);
 39:   MatCreateSubMatrices(BAIJ,n,irow,icol,MAT_REUSE_MATRIX,&subBAIJ);

 41:   /* irow and icol must be same for SBAIJ matrices! */
 42:   icol[0] = icol[1] = isrow;
 43:   MatCreateSubMatrices(SBAIJ,n,irow,icol,MAT_INITIAL_MATRIX,&subSBAIJ);
 44:   MatCreateSubMatrices(SBAIJ,n,irow,icol,MAT_REUSE_MATRIX,&subSBAIJ);

 46:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 47:   if (rank == 0) {
 48:     MatView(subBAIJ[0],PETSC_VIEWER_STDOUT_SELF);
 49:     MatView(subSBAIJ[0],PETSC_VIEWER_STDOUT_SELF);
 50:   }

 52:   /* Free data structures */
 53:   ISDestroy(&isrow);
 54:   ISDestroy(&iscol);
 55:   MatDestroySubMatrices(n,&subBAIJ);
 56:   MatDestroySubMatrices(n,&subSBAIJ);
 57:   MatDestroy(&BAIJ);
 58:   MatDestroy(&SBAIJ);

 60:   PetscFinalize();
 61:   return 0;
 62: }