Actual source code: ex59.c
2: static char help[] = "Tests MatCreateSubmatrix() in parallel.";
4: #include <petscmat.h>
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;
15: PetscFunctionBeginUser;
16: PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
17: PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_matmatmult", &test_matmatmult, NULL));
19: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
20: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
21: n = 2 * size;
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));
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: }
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));
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));
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));
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: }
65: PetscCall(ISDestroy(&isrow));
66: PetscCall(ISDestroy(&iscol));
67: PetscCall(MatDestroy(&A));
68: PetscCall(MatDestroy(&C));
69: PetscCall(PetscFinalize());
70: return 0;
71: }
73: /*TEST
75: test:
77: test:
78: suffix: 2
79: nsize: 3
81: test:
82: suffix: 2_baij
83: nsize: 3
84: args: -mat_type baij
86: test:
87: suffix: 2_sbaij
88: nsize: 3
89: args: -mat_type sbaij
91: test:
92: suffix: baij
93: args: -mat_type baij
94: output_file: output/ex59_1_baij.out
96: test:
97: suffix: sbaij
98: args: -mat_type sbaij
99: output_file: output/ex59_1_sbaij.out
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
109: TEST*/