1 2 static char help[] = "Tests MatCreateSubmatrix() in parallel."; 3 4 #include <petscmat.h> 5 6 int main(int argc, char **args) 7 { 8 Mat C, A; 9 PetscInt i, j, m = 3, n = 2, rstart, rend; 10 PetscMPIInt size, rank; 11 PetscScalar v; 12 IS isrow, iscol; 13 PetscBool test_matmatmult = PETSC_FALSE; 14 15 PetscFunctionBeginUser; 16 PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 17 PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_matmatmult", &test_matmatmult, NULL)); 18 19 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 20 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 21 n = 2 * size; 22 23 PetscCall(MatCreate(PETSC_COMM_WORLD, &C)); 24 PetscCall(MatSetSizes(C, PETSC_DECIDE, PETSC_DECIDE, m * n, m * n)); 25 PetscCall(MatSetFromOptions(C)); 26 PetscCall(MatSetUp(C)); 27 28 /* 29 This is JUST to generate a nice test matrix, all processors fill up 30 the entire matrix. This is not something one would ever do in practice. 31 */ 32 PetscCall(MatGetOwnershipRange(C, &rstart, &rend)); 33 for (i = rstart; i < rend; i++) { 34 for (j = 0; j < m * n; j++) { 35 v = i + j + 1; 36 PetscCall(MatSetValues(C, 1, &i, 1, &j, &v, INSERT_VALUES)); 37 } 38 } 39 40 PetscCall(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY)); 41 PetscCall(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY)); 42 PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_ASCII_COMMON)); 43 PetscCall(MatView(C, PETSC_VIEWER_STDOUT_WORLD)); 44 45 /* 46 Generate a new matrix consisting of every second row and column of 47 the original matrix 48 */ 49 PetscCall(MatGetOwnershipRange(C, &rstart, &rend)); 50 /* Create parallel IS with the rows we want on THIS processor */ 51 PetscCall(ISCreateStride(PETSC_COMM_WORLD, (rend - rstart) / 2, rstart, 2, &isrow)); 52 /* Create parallel IS with the rows we want on THIS processor (same as rows for now) */ 53 PetscCall(ISCreateStride(PETSC_COMM_WORLD, (rend - rstart) / 2, rstart, 2, &iscol)); 54 55 PetscCall(MatCreateSubMatrix(C, isrow, iscol, MAT_INITIAL_MATRIX, &A)); 56 PetscCall(MatCreateSubMatrix(C, isrow, iscol, MAT_REUSE_MATRIX, &A)); 57 PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 58 59 if (test_matmatmult) { 60 PetscCall(MatDestroy(&C)); 61 PetscCall(MatMatMult(A, A, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C)); 62 PetscCall(MatView(C, PETSC_VIEWER_STDOUT_WORLD)); 63 } 64 65 PetscCall(ISDestroy(&isrow)); 66 PetscCall(ISDestroy(&iscol)); 67 PetscCall(MatDestroy(&A)); 68 PetscCall(MatDestroy(&C)); 69 PetscCall(PetscFinalize()); 70 return 0; 71 } 72 73 /*TEST 74 75 test: 76 77 test: 78 suffix: 2 79 nsize: 3 80 81 test: 82 suffix: 2_baij 83 nsize: 3 84 args: -mat_type baij 85 86 test: 87 suffix: 2_sbaij 88 nsize: 3 89 args: -mat_type sbaij 90 91 test: 92 suffix: baij 93 args: -mat_type baij 94 output_file: output/ex59_1_baij.out 95 96 test: 97 suffix: sbaij 98 args: -mat_type sbaij 99 output_file: output/ex59_1_sbaij.out 100 101 test: 102 suffix: kok 103 nsize: 3 104 requires: kokkos_kernels 105 args: -mat_type aijkokkos -test_matmatmult 106 filter: grep -v -i type 107 output_file: output/ex59_kok.out 108 109 TEST*/ 110