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