1 2 static char help[] = "Creates a matrix, inserts some values, and tests MatCreateSubMatrices() and MatZeroEntries().\n\n"; 3 4 #include <petscmat.h> 5 6 int main(int argc,char **argv) 7 { 8 Mat mat,submat,submat1,*submatrices; 9 PetscInt m = 10,n = 10,i = 4,tmp,rstart,rend; 10 IS irow,icol; 11 PetscScalar value = 1.0; 12 PetscViewer sviewer; 13 PetscBool allA = PETSC_FALSE; 14 15 CHKERRQ(PetscInitialize(&argc,&argv,(char*)0,help)); 16 CHKERRQ(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_COMMON)); 17 CHKERRQ(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_SELF,PETSC_VIEWER_ASCII_COMMON)); 18 19 CHKERRQ(MatCreate(PETSC_COMM_WORLD,&mat)); 20 CHKERRQ(MatSetSizes(mat,PETSC_DECIDE,PETSC_DECIDE,m,n)); 21 CHKERRQ(MatSetFromOptions(mat)); 22 CHKERRQ(MatSetUp(mat)); 23 CHKERRQ(MatGetOwnershipRange(mat,&rstart,&rend)); 24 for (i=rstart; i<rend; i++) { 25 value = (PetscReal)i+1; tmp = i % 5; 26 CHKERRQ(MatSetValues(mat,1,&tmp,1,&i,&value,INSERT_VALUES)); 27 } 28 CHKERRQ(MatAssemblyBegin(mat,MAT_FINAL_ASSEMBLY)); 29 CHKERRQ(MatAssemblyEnd(mat,MAT_FINAL_ASSEMBLY)); 30 CHKERRQ(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Original matrix\n")); 31 CHKERRQ(MatView(mat,PETSC_VIEWER_STDOUT_WORLD)); 32 33 /* Test MatCreateSubMatrix_XXX_All(), i.e., submatrix = A */ 34 CHKERRQ(PetscOptionsGetBool(NULL,NULL,"-test_all",&allA,NULL)); 35 if (allA) { 36 CHKERRQ(ISCreateStride(PETSC_COMM_SELF,m,0,1,&irow)); 37 CHKERRQ(ISCreateStride(PETSC_COMM_SELF,n,0,1,&icol)); 38 CHKERRQ(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_INITIAL_MATRIX,&submatrices)); 39 CHKERRQ(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_REUSE_MATRIX,&submatrices)); 40 submat = *submatrices; 41 42 /* sviewer will cause the submatrices (one per processor) to be printed in the correct order */ 43 CHKERRQ(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"\nSubmatrices with all\n")); 44 CHKERRQ(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"--------------------\n")); 45 CHKERRQ(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 46 CHKERRQ(MatView(submat,sviewer)); 47 CHKERRQ(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 48 CHKERRQ(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 49 50 CHKERRQ(ISDestroy(&irow)); 51 CHKERRQ(ISDestroy(&icol)); 52 53 /* test getting a reference on a submat */ 54 CHKERRQ(PetscObjectReference((PetscObject)submat)); 55 CHKERRQ(MatDestroySubMatrices(1,&submatrices)); 56 CHKERRQ(MatDestroy(&submat)); 57 } 58 59 /* Form submatrix with rows 2-4 and columns 4-8 */ 60 CHKERRQ(ISCreateStride(PETSC_COMM_SELF,3,2,1,&irow)); 61 CHKERRQ(ISCreateStride(PETSC_COMM_SELF,5,4,1,&icol)); 62 CHKERRQ(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_INITIAL_MATRIX,&submatrices)); 63 submat = *submatrices; 64 65 /* Test reuse submatrices */ 66 CHKERRQ(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_REUSE_MATRIX,&submatrices)); 67 68 /* sviewer will cause the submatrices (one per processor) to be printed in the correct order */ 69 CHKERRQ(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"\nSubmatrices\n")); 70 CHKERRQ(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 71 CHKERRQ(MatView(submat,sviewer)); 72 CHKERRQ(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 73 CHKERRQ(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 74 CHKERRQ(PetscObjectReference((PetscObject)submat)); 75 CHKERRQ(MatDestroySubMatrices(1,&submatrices)); 76 CHKERRQ(MatDestroy(&submat)); 77 78 /* Form submatrix with rows 2-4 and all columns */ 79 CHKERRQ(ISDestroy(&icol)); 80 CHKERRQ(ISCreateStride(PETSC_COMM_SELF,10,0,1,&icol)); 81 CHKERRQ(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_INITIAL_MATRIX,&submatrices)); 82 CHKERRQ(MatCreateSubMatrices(mat,1,&irow,&icol,MAT_REUSE_MATRIX,&submatrices)); 83 submat = *submatrices; 84 85 CHKERRQ(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"\nSubmatrices with allcolumns\n")); 86 CHKERRQ(PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 87 CHKERRQ(MatView(submat,sviewer)); 88 CHKERRQ(PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD,PETSC_COMM_SELF,&sviewer)); 89 CHKERRQ(PetscViewerFlush(PETSC_VIEWER_STDOUT_WORLD)); 90 91 /* Test MatDuplicate */ 92 CHKERRQ(MatDuplicate(submat,MAT_COPY_VALUES,&submat1)); 93 CHKERRQ(MatDestroy(&submat1)); 94 95 /* Zero the original matrix */ 96 CHKERRQ(PetscViewerASCIIPrintf(PETSC_VIEWER_STDOUT_WORLD,"Original zeroed matrix\n")); 97 CHKERRQ(MatZeroEntries(mat)); 98 CHKERRQ(MatView(mat,PETSC_VIEWER_STDOUT_WORLD)); 99 100 CHKERRQ(ISDestroy(&irow)); 101 CHKERRQ(ISDestroy(&icol)); 102 CHKERRQ(PetscObjectReference((PetscObject)submat)); 103 CHKERRQ(MatDestroySubMatrices(1,&submatrices)); 104 CHKERRQ(MatDestroy(&submat)); 105 CHKERRQ(MatDestroy(&mat)); 106 CHKERRQ(PetscFinalize()); 107 return 0; 108 } 109 110 /*TEST 111 112 test: 113 args: -mat_type aij 114 115 test: 116 suffix: 2 117 args: -mat_type dense 118 119 test: 120 suffix: 3 121 nsize: 3 122 args: -mat_type aij 123 124 test: 125 suffix: 4 126 nsize: 3 127 args: -mat_type dense 128 129 test: 130 suffix: 5 131 nsize: 3 132 args: -mat_type aij -test_all 133 134 TEST*/ 135