1 static char help[] = "Tests various routines in MatKAIJ format.\n";
2
3 #include <petscmat.h>
4 #define IMAX 15
5
main(int argc,char ** args)6 int main(int argc, char **args)
7 {
8 Mat A, B, TA;
9 PetscScalar *S, *T;
10 PetscViewer fd;
11 char file[PETSC_MAX_PATH_LEN];
12 PetscInt m, n, M, N, p = 1, q = 1, i, j;
13 PetscMPIInt rank, size;
14 PetscBool flg;
15
16 PetscFunctionBeginUser;
17 PetscCall(PetscInitialize(&argc, &args, NULL, help));
18 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
19 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
20
21 /* Load AIJ matrix A */
22 PetscCall(PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &flg));
23 PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_USER, "Must indicate binary file with the -f option");
24 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &fd));
25 PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
26 PetscCall(MatLoad(A, fd));
27 PetscCall(PetscViewerDestroy(&fd));
28
29 /* Get dof, then create S and T */
30 PetscCall(PetscOptionsGetInt(NULL, NULL, "-p", &p, NULL));
31 PetscCall(PetscOptionsGetInt(NULL, NULL, "-q", &q, NULL));
32 PetscCall(PetscMalloc2(p * q, &S, p * q, &T));
33 for (i = 0; i < p * q; i++) S[i] = 0;
34
35 for (i = 0; i < p; i++) {
36 for (j = 0; j < q; j++) {
37 /* Set some random non-zero values */
38 S[i + p * j] = ((PetscReal)((i + 1) * (j + 1))) / ((PetscReal)(p + q));
39 T[i + p * j] = ((PetscReal)((p - i) + j)) / ((PetscReal)(p * q));
40 }
41 }
42
43 /* Test KAIJ when both S & T are not NULL */
44
45 /* Create KAIJ matrix TA */
46 PetscCall(MatCreateKAIJ(A, p, q, S, T, &TA));
47 PetscCall(MatGetLocalSize(A, &m, &n));
48 PetscCall(MatGetSize(A, &M, &N));
49
50 PetscCall(MatConvert(TA, MATAIJ, MAT_INITIAL_MATRIX, &B));
51
52 /* Test MatKAIJGetScaledIdentity() */
53 PetscCall(MatKAIJGetScaledIdentity(TA, &flg));
54 PetscCheck(!flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in Test 1: MatKAIJGetScaledIdentity()");
55 /* Test MatMult() */
56 PetscCall(MatMultEqual(TA, B, 10, &flg));
57 PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_CONV_FAILED, "Error in Test 1: MatMult() for KAIJ matrix");
58 /* Test MatMultAdd() */
59 PetscCall(MatMultAddEqual(TA, B, 10, &flg));
60 PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_CONV_FAILED, "Error in Test 1: MatMultAdd() for KAIJ matrix");
61
62 PetscCall(MatDestroy(&TA));
63 PetscCall(MatDestroy(&B));
64
65 /* Test KAIJ when S is NULL */
66
67 /* Create KAIJ matrix TA */
68 PetscCall(MatCreateKAIJ(A, p, q, NULL, T, &TA));
69 PetscCall(MatGetLocalSize(A, &m, &n));
70 PetscCall(MatGetSize(A, &M, &N));
71
72 PetscCall(MatConvert(TA, MATAIJ, MAT_INITIAL_MATRIX, &B));
73
74 /* Test MatKAIJGetScaledIdentity() */
75 PetscCall(MatKAIJGetScaledIdentity(TA, &flg));
76 PetscCheck(!flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in Test 2: MatKAIJGetScaledIdentity()");
77 /* Test MatMult() */
78 PetscCall(MatMultEqual(TA, B, 10, &flg));
79 PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_CONV_FAILED, "Error in Test 2: MatMult() for KAIJ matrix");
80 /* Test MatMultAdd() */
81 PetscCall(MatMultAddEqual(TA, B, 10, &flg));
82 PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_CONV_FAILED, "Error in Test 2: MatMultAdd() for KAIJ matrix");
83
84 PetscCall(MatDestroy(&TA));
85 PetscCall(MatDestroy(&B));
86
87 /* Test KAIJ when T is NULL */
88
89 /* Create KAIJ matrix TA */
90 PetscCall(MatCreateKAIJ(A, p, q, S, NULL, &TA));
91 PetscCall(MatGetLocalSize(A, &m, &n));
92 PetscCall(MatGetSize(A, &M, &N));
93
94 PetscCall(MatConvert(TA, MATAIJ, MAT_INITIAL_MATRIX, &B));
95
96 /* Test MatKAIJGetScaledIdentity() */
97 PetscCall(MatKAIJGetScaledIdentity(TA, &flg));
98 PetscCheck(!flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in Test 3: MatKAIJGetScaledIdentity()");
99 /* Test MatMult() */
100 PetscCall(MatMultEqual(TA, B, 10, &flg));
101 PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_CONV_FAILED, "Error in Test 3: MatMult() for KAIJ matrix");
102 /* Test MatMultAdd() */
103 PetscCall(MatMultAddEqual(TA, B, 10, &flg));
104 PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_CONV_FAILED, "Error in Test 3: MatMultAdd() for KAIJ matrix");
105
106 PetscCall(MatDestroy(&TA));
107 PetscCall(MatDestroy(&B));
108
109 /* Test KAIJ when T is an identity matrix */
110
111 if (p == q) {
112 for (i = 0; i < p; i++) {
113 for (j = 0; j < q; j++) {
114 if (i == j) T[i + j * p] = 1.0;
115 else T[i + j * p] = 0.0;
116 }
117 }
118
119 /* Create KAIJ matrix TA */
120 PetscCall(MatCreateKAIJ(A, p, q, S, T, &TA));
121 PetscCall(MatGetLocalSize(A, &m, &n));
122 PetscCall(MatGetSize(A, &M, &N));
123
124 PetscCall(MatConvert(TA, MATAIJ, MAT_INITIAL_MATRIX, &B));
125
126 /* Test MatKAIJGetScaledIdentity() */
127 PetscCall(MatKAIJGetScaledIdentity(TA, &flg));
128 PetscCheck(!flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in Test 4: MatKAIJGetScaledIdentity()");
129 /* Test MatMult() */
130 PetscCall(MatMultEqual(TA, B, 10, &flg));
131 PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_CONV_FAILED, "Error in Test 4: MatMult() for KAIJ matrix");
132 /* Test MatMultAdd() */
133 PetscCall(MatMultAddEqual(TA, B, 10, &flg));
134 PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_CONV_FAILED, "Error in Test 4: MatMultAdd() for KAIJ matrix");
135
136 PetscCall(MatDestroy(&TA));
137 PetscCall(MatDestroy(&B));
138
139 PetscCall(MatCreateKAIJ(A, p, q, NULL, T, &TA));
140 /* Test MatKAIJGetScaledIdentity() */
141 PetscCall(MatKAIJGetScaledIdentity(TA, &flg));
142 PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in Test 5: MatKAIJGetScaledIdentity()");
143 PetscCall(MatDestroy(&TA));
144
145 for (i = 0; i < p; i++) {
146 for (j = 0; j < q; j++) {
147 if (i == j) S[i + j * p] = T[i + j * p] = 2.0;
148 else S[i + j * p] = T[i + j * p] = 0.0;
149 }
150 }
151 PetscCall(MatCreateKAIJ(A, p, q, S, T, &TA));
152 /* Test MatKAIJGetScaledIdentity() */
153 PetscCall(MatKAIJGetScaledIdentity(TA, &flg));
154 PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Error in Test 6: MatKAIJGetScaledIdentity()");
155 PetscCall(MatDestroy(&TA));
156 }
157
158 /* Done with all tests */
159
160 PetscCall(PetscFree2(S, T));
161 PetscCall(MatDestroy(&A));
162 PetscCall(PetscFinalize());
163 return 0;
164 }
165
166 /*TEST
167
168 build:
169 requires: !complex
170
171 test:
172 requires: datafilespath !complex double !defined(PETSC_USE_64BIT_INDICES)
173 output_file: output/empty.out
174 nsize: {{1 2 3 4}}
175 args: -f ${DATAFILESPATH}/matrices/small -p {{2 3 7}} -q {{3 7}} -viewer_binary_skip_info
176
177 TEST*/
178