xref: /petsc/src/mat/tests/ex59.c (revision 98d129c30f3ee9fdddc40fdbc5a989b7be64f888)
1 static char help[] = "Tests MatCreateSubmatrix() in parallel.";
2 
3 #include <petscmat.h>
4 
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;
13 
14   PetscFunctionBeginUser;
15   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
16   PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_matmatmult", &test_matmatmult, NULL));
17 
18   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
19   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
20   n = 2 * size;
21 
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));
26 
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   }
38 
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));
43 
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));
53 
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));
57 
58   if (test_matmatmult) {
59     PetscCall(MatDestroy(&C));
60     PetscCall(MatMatMult(A, A, MAT_INITIAL_MATRIX, PETSC_DEFAULT, &C));
61     PetscCall(MatView(C, PETSC_VIEWER_STDOUT_WORLD));
62   }
63 
64   PetscCall(ISDestroy(&isrow));
65   PetscCall(ISDestroy(&iscol));
66   PetscCall(MatDestroy(&A));
67   PetscCall(MatDestroy(&C));
68   PetscCall(PetscFinalize());
69   return 0;
70 }
71 
72 /*TEST
73 
74    test:
75 
76    test:
77       suffix: 2
78       nsize: 3
79 
80    test:
81       suffix: 2_baij
82       nsize: 3
83       args: -mat_type baij
84 
85    test:
86       suffix: 2_sbaij
87       nsize: 3
88       args: -mat_type sbaij
89 
90    test:
91       suffix: baij
92       args: -mat_type baij
93       output_file: output/ex59_1_baij.out
94 
95    test:
96       suffix: sbaij
97       args: -mat_type sbaij
98       output_file: output/ex59_1_sbaij.out
99 
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
107 
108 TEST*/
109