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