Actual source code: ex26.c
1: static char help[] = "Tests MatGetRowIJ for SeqAIJ, SeqBAIJ and SeqSBAIJ\n\n";
3: #include <petscmat.h>
5: PetscErrorCode DumpCSR(Mat A, PetscInt shift, PetscBool symmetric, PetscBool compressed)
6: {
7: MatType type;
8: PetscInt i, j, nr, bs = 1;
9: const PetscInt *ia, *ja;
10: PetscBool done, isseqbaij, isseqsbaij;
12: PetscFunctionBeginUser;
13: PetscCall(PetscObjectBaseTypeCompare((PetscObject)A, MATSEQBAIJ, &isseqbaij));
14: PetscCall(PetscObjectBaseTypeCompare((PetscObject)A, MATSEQSBAIJ, &isseqsbaij));
15: if (isseqbaij || isseqsbaij) PetscCall(MatGetBlockSize(A, &bs));
16: PetscCall(MatGetType(A, &type));
17: PetscCall(MatGetRowIJ(A, shift, symmetric, compressed, &nr, &ia, &ja, &done));
18: PetscCall(PetscPrintf(PETSC_COMM_SELF, "===========================================================\n"));
19: PetscCall(PetscPrintf(PETSC_COMM_SELF, "CSR for %s: shift %" PetscInt_FMT " symmetric %" PetscInt_FMT " compressed %" PetscInt_FMT "\n", type, shift, (PetscInt)symmetric, (PetscInt)compressed));
20: for (i = 0; i < nr; i++) {
21: PetscCall(PetscPrintf(PETSC_COMM_SELF, "%" PetscInt_FMT ":", i + shift));
22: for (j = ia[i]; j < ia[i + 1]; j++) PetscCall(PetscPrintf(PETSC_COMM_SELF, " %" PetscInt_FMT, ja[j - shift]));
23: PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n"));
24: }
25: PetscCall(MatRestoreRowIJ(A, shift, symmetric, compressed, &nr, &ia, &ja, &done));
26: PetscFunctionReturn(PETSC_SUCCESS);
27: }
29: int main(int argc, char **args)
30: {
31: Mat A, B, C;
32: PetscInt i, j, k, m = 3, n = 3, bs = 1;
33: PetscMPIInt size;
35: PetscFunctionBeginUser;
36: PetscCall(PetscInitialize(&argc, &args, NULL, help));
37: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
38: PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!");
39: PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL));
40: PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
41: PetscCall(PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL));
42: /* adjust sizes by block size */
43: if (m % bs) m += bs - m % bs;
44: if (n % bs) n += bs - n % bs;
46: PetscCall(MatCreate(PETSC_COMM_SELF, &A));
47: PetscCall(MatSetSizes(A, m * n, m * n, PETSC_DECIDE, PETSC_DECIDE));
48: PetscCall(MatSetBlockSize(A, bs));
49: PetscCall(MatSetType(A, MATSEQAIJ));
50: PetscCall(MatSetUp(A));
51: PetscCall(MatCreate(PETSC_COMM_SELF, &B));
52: PetscCall(MatSetSizes(B, m * n, m * n, PETSC_DECIDE, PETSC_DECIDE));
53: PetscCall(MatSetBlockSize(B, bs));
54: PetscCall(MatSetType(B, MATSEQBAIJ));
55: PetscCall(MatSetUp(B));
56: PetscCall(MatCreate(PETSC_COMM_SELF, &C));
57: PetscCall(MatSetSizes(C, m * n, m * n, PETSC_DECIDE, PETSC_DECIDE));
58: PetscCall(MatSetBlockSize(C, bs));
59: PetscCall(MatSetType(C, MATSEQSBAIJ));
60: PetscCall(MatSetUp(C));
61: PetscCall(MatSetOption(C, MAT_IGNORE_LOWER_TRIANGULAR, PETSC_TRUE));
63: for (i = 0; i < m; i++) {
64: for (j = 0; j < n; j++) {
65: PetscScalar v = -1.0;
66: PetscInt Ii = j + n * i, J;
67: J = Ii - n;
68: if (J >= 0) {
69: PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES));
70: PetscCall(MatSetValues(B, 1, &Ii, 1, &J, &v, INSERT_VALUES));
71: PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES));
72: }
73: J = Ii + n;
74: if (J < m * n) {
75: PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES));
76: PetscCall(MatSetValues(B, 1, &Ii, 1, &J, &v, INSERT_VALUES));
77: PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES));
78: }
79: J = Ii - 1;
80: if (J >= 0) {
81: PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES));
82: PetscCall(MatSetValues(B, 1, &Ii, 1, &J, &v, INSERT_VALUES));
83: PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES));
84: }
85: J = Ii + 1;
86: if (J < m * n) {
87: PetscCall(MatSetValues(A, 1, &Ii, 1, &J, &v, INSERT_VALUES));
88: PetscCall(MatSetValues(B, 1, &Ii, 1, &J, &v, INSERT_VALUES));
89: PetscCall(MatSetValues(C, 1, &Ii, 1, &J, &v, INSERT_VALUES));
90: }
91: v = 4.0;
92: PetscCall(MatSetValues(A, 1, &Ii, 1, &Ii, &v, INSERT_VALUES));
93: PetscCall(MatSetValues(B, 1, &Ii, 1, &Ii, &v, INSERT_VALUES));
94: PetscCall(MatSetValues(C, 1, &Ii, 1, &Ii, &v, INSERT_VALUES));
95: }
96: }
97: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
98: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
99: PetscCall(MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY));
100: PetscCall(MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY));
101: PetscCall(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY));
102: PetscCall(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY));
104: /* test MatGetRowIJ for the three Mat types */
105: PetscCall(MatView(A, NULL));
106: PetscCall(MatView(B, NULL));
107: PetscCall(MatView(C, NULL));
108: for (i = 0; i < 2; i++) {
109: PetscInt shift = i;
110: for (j = 0; j < 2; j++) {
111: PetscBool symmetric = ((j > 0) ? PETSC_FALSE : PETSC_TRUE);
112: for (k = 0; k < 2; k++) {
113: PetscBool compressed = ((k > 0) ? PETSC_FALSE : PETSC_TRUE);
114: PetscCall(DumpCSR(A, shift, symmetric, compressed));
115: PetscCall(DumpCSR(B, shift, symmetric, compressed));
116: PetscCall(DumpCSR(C, shift, symmetric, compressed));
117: }
118: }
119: }
120: PetscCall(MatDestroy(&A));
121: PetscCall(MatDestroy(&B));
122: PetscCall(MatDestroy(&C));
123: PetscCall(PetscFinalize());
124: return 0;
125: }
127: /*TEST
129: test:
131: test:
132: suffix: 2
133: args: -bs 2
135: TEST*/