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