2: static char help[] = "Creates a matrix, inserts some values, and tests MatCreateSubMatrices() and MatZeroEntries().\n\n";

  4:  #include <petscmat.h>

  6: int main(int argc,char **argv)
  7: {
  8:   Mat            mat,submat,*submatrices;
  9:   PetscInt       m = 10,n = 10,i = 4,tmp,rstart,rend;
 11:   IS             irow,icol;
 12:   PetscScalar    value = 1.0;
 13:   PetscViewer    sviewer;
 14:   PetscBool      allA = PETSC_FALSE;

 16:   PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
 17:   PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_COMMON);
 18:   PetscViewerPushFormat(PETSC_VIEWER_STDOUT_SELF,PETSC_VIEWER_ASCII_COMMON);

 20:   MatCreate(PETSC_COMM_WORLD,&mat);
 21:   MatSetSizes(mat,PETSC_DECIDE,PETSC_DECIDE,m,n);
 22:   MatSetFromOptions(mat);
 23:   MatSetUp(mat);
 24:   MatGetOwnershipRange(mat,&rstart,&rend);
 25:   for (i=rstart; i<rend; i++) {
 26:     value = (PetscReal)i+1; tmp = i % 5;
 27:     MatSetValues(mat,1,&tmp,1,&i,&value,INSERT_VALUES);
 28:   }
 29:   MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY);
 30:   MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY);
 31:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Original matrix\n");
 32:   MatView(mat,PETSC_VIEWER_STDOUT_WORLD);

 34:   /* Test MatCreateSubMatrix_XXX_All(), i.e., submatrix = A */
 35:   PetscOptionsGetBool(NULL,NULL,"-test_all",&allA,NULL);
 36:   if (allA) {
 37:     ISCreateStride(PETSC_COMM_SELF,m,0,1,&irow);
 38:     ISCreateStride(PETSC_COMM_SELF,n,0,1,&icol);
 39:     MatCreateSubMatrices(mat,1,&irow,&icol,MAT_INITIAL_MATRIX,&submatrices);
 40:     MatCreateSubMatrices(mat,1,&irow,&icol,MAT_REUSE_MATRIX,&submatrices);
 41:     submat = *submatrices;

 43:     /* sviewer will cause the submatrices (one per processor) to be printed in the correct order */
 44:     PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"\nSubmatrices with all\n");
 45:     PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"--------------------\n");
 46:     PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer);
 47:     MatView(submat,sviewer);
 48:     PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer);
 49:     PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);

 51:     ISDestroy(&irow);
 52:     ISDestroy(&icol);

 54:     /* test getting a reference on a submat */
 55:     PetscObjectReference((PetscObject)submat);
 56:     MatDestroySubMatrices(1,&submatrices);
 57:     MatDestroy(&submat);
 58:   }

 60:   /* Form submatrix with rows 2-4 and columns 4-8 */
 61:   ISCreateStride(PETSC_COMM_SELF,3,2,1,&irow);
 62:   ISCreateStride(PETSC_COMM_SELF,5,4,1,&icol);
 63:   MatCreateSubMatrices(mat,1,&irow,&icol,MAT_INITIAL_MATRIX,&submatrices);
 64:   submat = *submatrices;

 66:   /* Test reuse submatrices */
 67:   MatCreateSubMatrices(mat,1,&irow,&icol,MAT_REUSE_MATRIX,&submatrices);

 69:   /* sviewer will cause the submatrices (one per processor) to be printed in the correct order */
 70:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"\nSubmatrices\n");
 71:   PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer);
 72:   MatView(submat,sviewer);
 73:   PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer);
 74:   PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);
 75:   PetscObjectReference((PetscObject)submat);
 76:   MatDestroySubMatrices(1,&submatrices);
 77:   MatDestroy(&submat);

 79:   /* Form submatrix with rows 2-4 and all columns */
 80:   ISDestroy(&icol);
 81:   ISCreateStride(PETSC_COMM_SELF,10,0,1,&icol);
 82:   MatCreateSubMatrices(mat,1,&irow,&icol,MAT_INITIAL_MATRIX,&submatrices);
 83:   MatCreateSubMatrices(mat,1,&irow,&icol,MAT_REUSE_MATRIX,&submatrices);
 84:   submat = *submatrices;

 86:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"\nSubmatrices with allcolumns\n");
 87:   PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer);
 88:   MatView(submat,sviewer);
 89:   PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer);
 90:   PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD);

 92:   /* Zero the original matrix */
 93:   PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Original zeroed matrix\n");
 94:   MatZeroEntries(mat);
 95:   MatView(mat,PETSC_VIEWER_STDOUT_WORLD);

 97:   ISDestroy(&irow);
 98:   ISDestroy(&icol);
 99:   PetscObjectReference((PetscObject)submat);
100:   MatDestroySubMatrices(1,&submatrices);
101:   MatDestroy(&submat);
102:   MatDestroy(&mat);
103:   PetscFinalize();
104:   return ierr;
105: }