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