Actual source code: ex35.c
2: static char help[] = "Tests MatCreateSubMatrices().\n\n";
4: #include <petscmat.h>
6: int main(int argc, char **args)
7: {
8: Mat A, B, *Bsub;
9: PetscInt i, j, m = 6, n = 6, N = 36, Ii, J;
10: PetscScalar v;
11: IS isrow;
13: PetscFunctionBeginUser;
14: PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
15: PetscCall(MatCreateSeqAIJ(PETSC_COMM_WORLD, N, N, 5, NULL, &A));
16: for (i = 0; i < m; i++) {
17: for (j = 0; j < n; j++) {
18: v = -1.0;
19: Ii = j + n * i;
20: if (i > 0) {
21: J = Ii - n;
22: PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES));
23: }
24: if (i < m - 1) {
25: J = Ii + n;
26: PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES));
27: }
28: if (j > 0) {
29: J = Ii - 1;
30: PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES));
31: }
32: if (j < n - 1) {
33: J = Ii + 1;
34: PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES));
35: }
36: v = 4.0;
37: PetscCall(MatSetValues(A, 1, &Ii, 1, &Ii, &v, INSERT_VALUES));
38: }
39: }
40: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
41: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
42: PetscCall(MatView(A, PETSC_VIEWER_STDOUT_SELF));
44: /* take the first diagonal block */
45: PetscCall(ISCreateStride(PETSC_COMM_WORLD, m, 0, 1, &isrow));
46: PetscCall(MatCreateSubMatrices(A, 1, &isrow, &isrow, MAT_INITIAL_MATRIX, &Bsub));
47: B = *Bsub;
48: PetscCall(ISDestroy(&isrow));
49: PetscCall(MatView(B, PETSC_VIEWER_STDOUT_SELF));
50: PetscCall(MatDestroySubMatrices(1, &Bsub));
52: /* take a strided block */
53: PetscCall(ISCreateStride(PETSC_COMM_WORLD, m, 0, 2, &isrow));
54: PetscCall(MatCreateSubMatrices(A, 1, &isrow, &isrow, MAT_INITIAL_MATRIX, &Bsub));
55: B = *Bsub;
56: PetscCall(ISDestroy(&isrow));
57: PetscCall(MatView(B, PETSC_VIEWER_STDOUT_SELF));
58: PetscCall(MatDestroySubMatrices(1, &Bsub));
60: /* take the last block */
61: PetscCall(ISCreateStride(PETSC_COMM_WORLD, m, N - m - 1, 1, &isrow));
62: PetscCall(MatCreateSubMatrices(A, 1, &isrow, &isrow, MAT_INITIAL_MATRIX, &Bsub));
63: B = *Bsub;
64: PetscCall(ISDestroy(&isrow));
65: PetscCall(MatView(B, PETSC_VIEWER_STDOUT_SELF));
67: PetscCall(MatDestroySubMatrices(1, &Bsub));
68: PetscCall(MatDestroy(&A));
70: PetscCall(PetscFinalize());
71: return 0;
72: }
74: /*TEST
76: test:
78: TEST*/