1 static char help[] = "3D, tri-linear quadrilateral (Q1), displacement finite element formulation\n\
2 of linear elasticity. E=1.0, nu=0.25.\n\
3 Unit square domain with Dirichelet boundary condition on the y=0 side only.\n\
4 Load of 1.0 in x + 2y direction on all nodes (not a true uniform load).\n\
5 -ne <size> : number of (square) quadrilateral elements in each dimension\n\
6 -alpha <v> : scaling of material coefficient in embedded circle\n\n";
7
8 #include <petscksp.h>
9
10 static PetscBool log_stages = PETSC_TRUE;
11
MaybeLogStagePush(PetscLogStage stage)12 static PetscErrorCode MaybeLogStagePush(PetscLogStage stage)
13 {
14 return log_stages ? PetscLogStagePush(stage) : PETSC_SUCCESS;
15 }
16
MaybeLogStagePop(void)17 static PetscErrorCode MaybeLogStagePop(void)
18 {
19 return log_stages ? PetscLogStagePop() : PETSC_SUCCESS;
20 }
21
22 PetscErrorCode elem_3d_elast_v_25(PetscScalar *);
23
main(int argc,char ** args)24 int main(int argc, char **args)
25 {
26 Mat Amat;
27 PetscInt m, nn, M, Istart, Iend, i, j, k, ii, jj, kk, ic, ne = 4, id;
28 PetscReal x, y, z, h, *coords, soft_alpha = 1.e-3;
29 PetscBool two_solves = PETSC_FALSE, test_nonzero_cols = PETSC_FALSE, use_nearnullspace = PETSC_FALSE, test_late_bs = PETSC_FALSE, test_rap_bs = PETSC_FALSE;
30 Vec xx, bb;
31 KSP ksp;
32 MPI_Comm comm;
33 PetscMPIInt npe, mype;
34 PetscScalar DD[24][24], DD2[24][24];
35 PetscLogStage stage[6];
36 PetscScalar DD1[24][24];
37
38 PetscFunctionBeginUser;
39 PetscCall(PetscInitialize(&argc, &args, NULL, help));
40 comm = PETSC_COMM_WORLD;
41 PetscCallMPI(MPI_Comm_rank(comm, &mype));
42 PetscCallMPI(MPI_Comm_size(comm, &npe));
43
44 PetscOptionsBegin(comm, NULL, "3D bilinear Q1 elasticity options", "");
45 {
46 char nestring[256];
47 PetscCall(PetscSNPrintf(nestring, sizeof nestring, "number of elements in each direction, ne+1 must be a multiple of %" PetscInt_FMT " (sizes^{1/3})", (PetscInt)(PetscPowReal((PetscReal)npe, 1. / 3.) + .5)));
48 PetscCall(PetscOptionsInt("-ne", nestring, "", ne, &ne, NULL));
49 PetscCall(PetscOptionsBool("-log_stages", "Log stages of solve separately", "", log_stages, &log_stages, NULL));
50 PetscCall(PetscOptionsReal("-alpha", "material coefficient inside circle", "", soft_alpha, &soft_alpha, NULL));
51 PetscCall(PetscOptionsBool("-two_solves", "solve additional variant of the problem", "", two_solves, &two_solves, NULL));
52 PetscCall(PetscOptionsBool("-test_nonzero_cols", "nonzero test", "", test_nonzero_cols, &test_nonzero_cols, NULL));
53 PetscCall(PetscOptionsBool("-use_mat_nearnullspace", "MatNearNullSpace API test", "", use_nearnullspace, &use_nearnullspace, NULL));
54 PetscCall(PetscOptionsBool("-test_late_bs", "", "", test_late_bs, &test_late_bs, NULL));
55 PetscCall(PetscOptionsBool("-test_rap_bs", "", "", test_rap_bs, &test_rap_bs, NULL));
56 }
57 PetscOptionsEnd();
58
59 if (log_stages) {
60 PetscCall(PetscLogStageRegister("Setup", &stage[0]));
61 PetscCall(PetscLogStageRegister("Solve", &stage[1]));
62 PetscCall(PetscLogStageRegister("2nd Setup", &stage[2]));
63 PetscCall(PetscLogStageRegister("2nd Solve", &stage[3]));
64 PetscCall(PetscLogStageRegister("3rd Setup", &stage[4]));
65 PetscCall(PetscLogStageRegister("3rd Solve", &stage[5]));
66 } else {
67 for (i = 0; i < (PetscInt)PETSC_STATIC_ARRAY_LENGTH(stage); i++) stage[i] = -1;
68 }
69
70 h = 1. / ne;
71 nn = ne + 1;
72 /* ne*ne; number of global elements */
73 M = 3 * nn * nn * nn; /* global number of equations */
74 if (npe == 2) {
75 if (mype == 1) m = 0;
76 else m = nn * nn * nn;
77 npe = 1;
78 } else {
79 m = nn * nn * nn / npe;
80 if (mype == npe - 1) m = nn * nn * nn - (npe - 1) * m;
81 }
82 m *= 3; /* number of equations local*/
83 /* Setup solver */
84 PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
85 PetscCall(KSPSetComputeSingularValues(ksp, PETSC_TRUE));
86 PetscCall(KSPSetFromOptions(ksp));
87 {
88 /* configuration */
89 const PetscInt NP = (PetscInt)(PetscPowReal((PetscReal)npe, 1. / 3.) + .5);
90 const PetscInt ipx = mype % NP, ipy = (mype % (NP * NP)) / NP, ipz = mype / (NP * NP);
91 const PetscInt Ni0 = ipx * (nn / NP), Nj0 = ipy * (nn / NP), Nk0 = ipz * (nn / NP);
92 const PetscInt Ni1 = Ni0 + (m > 0 ? (nn / NP) : 0), Nj1 = Nj0 + (nn / NP), Nk1 = Nk0 + (nn / NP);
93 const PetscInt NN = nn / NP, id0 = ipz * nn * nn * NN + ipy * nn * NN * NN + ipx * NN * NN * NN;
94 PetscInt *d_nnz, *o_nnz, osz[4] = {0, 9, 15, 19}, nbc;
95 PetscScalar vv[24], v2[24];
96 PetscCheck(npe == NP * NP * NP, comm, PETSC_ERR_ARG_WRONG, "npe=%d: npe^{1/3} must be integer", npe);
97 PetscCheck(nn == NP * (nn / NP), comm, PETSC_ERR_ARG_WRONG, "-ne %" PetscInt_FMT ": (ne+1)%%(npe^{1/3}) must equal zero", ne);
98
99 /* count nnz */
100 PetscCall(PetscMalloc1(m + 1, &d_nnz));
101 PetscCall(PetscMalloc1(m + 1, &o_nnz));
102 for (i = Ni0, ic = 0; i < Ni1; i++) {
103 for (j = Nj0; j < Nj1; j++) {
104 for (k = Nk0; k < Nk1; k++) {
105 nbc = 0;
106 if (i == Ni0 || i == Ni1 - 1) nbc++;
107 if (j == Nj0 || j == Nj1 - 1) nbc++;
108 if (k == Nk0 || k == Nk1 - 1) nbc++;
109 for (jj = 0; jj < 3; jj++, ic++) {
110 d_nnz[ic] = 3 * (27 - osz[nbc]);
111 o_nnz[ic] = 3 * osz[nbc];
112 }
113 }
114 }
115 }
116 PetscCheck(ic == m, PETSC_COMM_SELF, PETSC_ERR_PLIB, "ic %" PetscInt_FMT " does not equal m %" PetscInt_FMT, ic, m);
117
118 /* create stiffness matrix */
119 PetscCall(MatCreate(comm, &Amat));
120 PetscCall(MatSetSizes(Amat, m, m, M, M));
121 if (!test_late_bs) PetscCall(MatSetBlockSize(Amat, 3));
122 PetscCall(MatSetType(Amat, MATAIJ));
123 PetscCall(MatSetOption(Amat, MAT_SPD, PETSC_TRUE));
124 PetscCall(MatSetOption(Amat, MAT_SPD_ETERNAL, PETSC_TRUE)); // this keeps CG after switch to negative
125 PetscCall(MatSetFromOptions(Amat));
126 PetscCall(MatSeqAIJSetPreallocation(Amat, 0, d_nnz));
127 PetscCall(MatMPIAIJSetPreallocation(Amat, 0, d_nnz, 0, o_nnz));
128
129 PetscCall(PetscFree(d_nnz));
130 PetscCall(PetscFree(o_nnz));
131 PetscCall(MatCreateVecs(Amat, &bb, &xx));
132
133 PetscCall(MatGetOwnershipRange(Amat, &Istart, &Iend));
134
135 PetscCheck(m == Iend - Istart, PETSC_COMM_SELF, PETSC_ERR_PLIB, "m %" PetscInt_FMT " does not equal Iend %" PetscInt_FMT " - Istart %" PetscInt_FMT, m, Iend, Istart);
136 /* generate element matrices */
137 {
138 PetscBool hasData = PETSC_TRUE;
139 if (!hasData) {
140 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\t No data is provided\n"));
141 for (i = 0; i < 24; i++) {
142 for (j = 0; j < 24; j++) {
143 if (i == j) DD1[i][j] = 1.0;
144 else DD1[i][j] = -.25;
145 }
146 }
147 } else {
148 /* Get array data */
149 PetscCall(elem_3d_elast_v_25((PetscScalar *)DD1));
150 }
151
152 /* BC version of element */
153 for (i = 0; i < 24; i++) {
154 for (j = 0; j < 24; j++) {
155 if (i < 12 || (j < 12 && !test_nonzero_cols)) {
156 if (i == j) DD2[i][j] = 0.1 * DD1[i][j];
157 else DD2[i][j] = 0.0;
158 } else DD2[i][j] = DD1[i][j];
159 }
160 }
161 /* element residual/load vector */
162 for (i = 0; i < 24; i++) {
163 if (i % 3 == 0) vv[i] = h * h;
164 else if (i % 3 == 1) vv[i] = 2.0 * h * h;
165 else vv[i] = .0;
166 }
167 for (i = 0; i < 24; i++) {
168 if (i % 3 == 0 && i >= 12) v2[i] = h * h;
169 else if (i % 3 == 1 && i >= 12) v2[i] = 2.0 * h * h;
170 else v2[i] = .0;
171 }
172 }
173
174 PetscCall(PetscMalloc1(m + 1, &coords));
175 coords[m] = -99.0;
176
177 /* forms the element stiffness and coordinates */
178 for (i = Ni0, ic = 0, ii = 0; i < Ni1; i++, ii++) {
179 for (j = Nj0, jj = 0; j < Nj1; j++, jj++) {
180 for (k = Nk0, kk = 0; k < Nk1; k++, kk++, ic++) {
181 /* coords */
182 x = coords[3 * ic] = h * (PetscReal)i;
183 y = coords[3 * ic + 1] = h * (PetscReal)j;
184 z = coords[3 * ic + 2] = h * (PetscReal)k;
185 /* matrix */
186 id = id0 + ii + NN * jj + NN * NN * kk;
187 if (i < ne && j < ne && k < ne) {
188 /* radius */
189 PetscReal radius = PetscSqrtReal((x - .5 + h / 2) * (x - .5 + h / 2) + (y - .5 + h / 2) * (y - .5 + h / 2) + (z - .5 + h / 2) * (z - .5 + h / 2));
190 PetscReal alpha = 1.0;
191 PetscInt jx, ix, idx[8], idx3[24];
192 idx[0] = id;
193 idx[1] = id + 1;
194 idx[2] = id + NN + 1;
195 idx[3] = id + NN;
196 idx[4] = id + NN * NN;
197 idx[5] = id + 1 + NN * NN;
198 idx[6] = id + NN + 1 + NN * NN;
199 idx[7] = id + NN + NN * NN;
200
201 /* correct indices */
202 if (i == Ni1 - 1 && Ni1 != nn) {
203 idx[1] += NN * (NN * NN - 1);
204 idx[2] += NN * (NN * NN - 1);
205 idx[5] += NN * (NN * NN - 1);
206 idx[6] += NN * (NN * NN - 1);
207 }
208 if (j == Nj1 - 1 && Nj1 != nn) {
209 idx[2] += NN * NN * (nn - 1);
210 idx[3] += NN * NN * (nn - 1);
211 idx[6] += NN * NN * (nn - 1);
212 idx[7] += NN * NN * (nn - 1);
213 }
214 if (k == Nk1 - 1 && Nk1 != nn) {
215 idx[4] += NN * (nn * nn - NN * NN);
216 idx[5] += NN * (nn * nn - NN * NN);
217 idx[6] += NN * (nn * nn - NN * NN);
218 idx[7] += NN * (nn * nn - NN * NN);
219 }
220
221 if (radius < 0.25) alpha = soft_alpha;
222
223 for (ix = 0; ix < 24; ix++) {
224 for (jx = 0; jx < 24; jx++) DD[ix][jx] = alpha * DD1[ix][jx];
225 }
226 if (k > 0) {
227 if (!test_late_bs) {
228 PetscCall(MatSetValuesBlocked(Amat, 8, idx, 8, idx, (const PetscScalar *)DD, ADD_VALUES));
229 PetscCall(VecSetValuesBlocked(bb, 8, idx, (const PetscScalar *)vv, ADD_VALUES));
230 } else {
231 for (ix = 0; ix < 8; ix++) {
232 idx3[3 * ix] = 3 * idx[ix];
233 idx3[3 * ix + 1] = 3 * idx[ix] + 1;
234 idx3[3 * ix + 2] = 3 * idx[ix] + 2;
235 }
236 PetscCall(MatSetValues(Amat, 24, idx3, 24, idx3, (const PetscScalar *)DD, ADD_VALUES));
237 PetscCall(VecSetValues(bb, 24, idx3, (const PetscScalar *)vv, ADD_VALUES));
238 }
239 } else {
240 /* a BC */
241 for (ix = 0; ix < 24; ix++) {
242 for (jx = 0; jx < 24; jx++) DD[ix][jx] = alpha * DD2[ix][jx];
243 }
244 if (!test_late_bs) {
245 PetscCall(MatSetValuesBlocked(Amat, 8, idx, 8, idx, (const PetscScalar *)DD, ADD_VALUES));
246 PetscCall(VecSetValuesBlocked(bb, 8, idx, (const PetscScalar *)v2, ADD_VALUES));
247 } else {
248 for (ix = 0; ix < 8; ix++) {
249 idx3[3 * ix] = 3 * idx[ix];
250 idx3[3 * ix + 1] = 3 * idx[ix] + 1;
251 idx3[3 * ix + 2] = 3 * idx[ix] + 2;
252 }
253 PetscCall(MatSetValues(Amat, 24, idx3, 24, idx3, (const PetscScalar *)DD, ADD_VALUES));
254 PetscCall(VecSetValues(bb, 24, idx3, (const PetscScalar *)v2, ADD_VALUES));
255 }
256 }
257 }
258 }
259 }
260 }
261 PetscCall(MatAssemblyBegin(Amat, MAT_FINAL_ASSEMBLY));
262 PetscCall(MatAssemblyEnd(Amat, MAT_FINAL_ASSEMBLY));
263 PetscCall(VecAssemblyBegin(bb));
264 PetscCall(VecAssemblyEnd(bb));
265 }
266 PetscCall(MatAssemblyBegin(Amat, MAT_FINAL_ASSEMBLY));
267 PetscCall(MatAssemblyEnd(Amat, MAT_FINAL_ASSEMBLY));
268 PetscCall(VecAssemblyBegin(bb));
269 PetscCall(VecAssemblyEnd(bb));
270 if (test_late_bs) {
271 PetscCall(VecSetBlockSize(xx, 3));
272 PetscCall(VecSetBlockSize(bb, 3));
273 PetscCall(MatSetBlockSize(Amat, 3));
274 }
275
276 if (!PETSC_TRUE) {
277 PetscViewer viewer;
278 PetscCall(PetscViewerASCIIOpen(comm, "Amat.m", &viewer));
279 PetscCall(PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_MATLAB));
280 PetscCall(MatView(Amat, viewer));
281 PetscCall(PetscViewerPopFormat(viewer));
282 PetscCall(PetscViewerDestroy(&viewer));
283 }
284
285 /* finish KSP/PC setup */
286 PetscCall(KSPSetOperators(ksp, Amat, Amat));
287 if (use_nearnullspace) {
288 MatNullSpace matnull;
289 Vec vec_coords;
290 PetscScalar *c;
291 PC pc;
292 PetscCall(VecCreate(MPI_COMM_WORLD, &vec_coords));
293 PetscCall(VecSetBlockSize(vec_coords, 3));
294 PetscCall(VecSetSizes(vec_coords, m, PETSC_DECIDE));
295 PetscCall(VecSetUp(vec_coords));
296 PetscCall(VecGetArray(vec_coords, &c));
297 for (i = 0; i < m; i++) c[i] = coords[i]; /* Copy since Scalar type might be Complex */
298 PetscCall(VecRestoreArray(vec_coords, &c));
299 PetscCall(MatNullSpaceCreateRigidBody(vec_coords, &matnull));
300 PetscCall(MatSetNearNullSpace(Amat, matnull));
301 PetscCall(MatNullSpaceDestroy(&matnull));
302 PetscCall(VecDestroy(&vec_coords));
303 PetscCall(KSPGetPC(ksp, &pc));
304 PetscCall(PCJacobiSetRowl1Scale(pc, 0.5));
305 } else {
306 PC pc;
307 PetscInt idx[] = {1, 2};
308 PetscCall(KSPGetPC(ksp, &pc));
309 PetscCall(PCSetCoordinates(pc, 3, m / 3, coords));
310 PetscCall(PCGAMGSetUseSAEstEig(pc, PETSC_FALSE));
311 PetscCall(PCGAMGSetLowMemoryFilter(pc, PETSC_TRUE));
312 PetscCall(PCGAMGMISkSetMinDegreeOrdering(pc, PETSC_TRUE));
313 PetscCall(PCGAMGSetAggressiveSquareGraph(pc, PETSC_FALSE));
314 PetscCall(PCGAMGSetInjectionIndex(pc, 2, idx)); // code coverage, same as command line
315 }
316
317 PetscCall(MaybeLogStagePush(stage[0]));
318
319 /* PC setup basically */
320 PetscCall(KSPSetUp(ksp));
321
322 PetscCall(MaybeLogStagePop());
323
324 if (test_rap_bs) {
325 PC pc;
326 Mat P, cmat;
327 KSP ksp2, cksp;
328 PetscCall(KSPGetPC(ksp, &pc));
329 PetscCall(PCMGGetLevels(pc, &i));
330 PetscCall(PCMGGetInterpolation(pc, i - 1, &P));
331 PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp2));
332 PetscCall(KSPSetOptionsPrefix(ksp2, "rap_"));
333 PetscCall(KSPSetFromOptions(ksp2));
334 PetscCall(KSPGetPC(ksp2, &pc));
335 PetscCall(PCSetType(pc, PCMG));
336 PetscCall(KSPSetOperators(ksp2, Amat, Amat));
337 PetscCall(PCMGSetLevels(pc, 2, NULL));
338 PetscCall(PCMGSetGalerkin(pc, PC_MG_GALERKIN_PMAT));
339 PetscCall(PCMGSetInterpolation(pc, 1, P));
340 PetscCall(VecSet(bb, 1.0));
341 PetscCall(KSPSolve(ksp2, bb, xx));
342 PetscCall(PCMGGetCoarseSolve(pc, &cksp));
343 PetscCall(KSPGetOperators(cksp, &cmat, &cmat));
344 PetscCall(MatViewFromOptions(cmat, NULL, "-rap_mat_view"));
345 /* Free work space and exit */
346 PetscCall(KSPDestroy(&ksp));
347 PetscCall(KSPDestroy(&ksp2));
348 PetscCall(VecDestroy(&xx));
349 PetscCall(VecDestroy(&bb));
350 PetscCall(MatDestroy(&Amat));
351 PetscCall(PetscFree(coords));
352 PetscCall(PetscFinalize());
353 return 0;
354 }
355
356 PetscCall(MaybeLogStagePush(stage[1]));
357
358 /* test BCs */
359 if (test_nonzero_cols) {
360 PetscCall(VecZeroEntries(xx));
361 if (mype == 0) PetscCall(VecSetValue(xx, 0, 1.0, INSERT_VALUES));
362 PetscCall(VecAssemblyBegin(xx));
363 PetscCall(VecAssemblyEnd(xx));
364 PetscCall(KSPSetInitialGuessNonzero(ksp, PETSC_TRUE));
365 }
366
367 /* 1st solve */
368 PetscCall(KSPSolve(ksp, bb, xx));
369
370 PetscCall(MaybeLogStagePop());
371
372 /* 2nd solve */
373 if (two_solves) {
374 PetscReal emax, emin;
375 PetscReal norm, norm2;
376 Vec res;
377
378 PetscCall(MaybeLogStagePush(stage[2]));
379 /* PC setup basically */
380 PetscCall(MatScale(Amat, -100000.0));
381 PetscCall(MatSetOption(Amat, MAT_SPD, PETSC_FALSE));
382 PetscCall(KSPSetOperators(ksp, Amat, Amat));
383 PetscCall(KSPSetUp(ksp));
384
385 PetscCall(MaybeLogStagePop());
386 PetscCall(MaybeLogStagePush(stage[3]));
387 PetscCall(KSPSolve(ksp, bb, xx));
388 PetscCall(KSPComputeExtremeSingularValues(ksp, &emax, &emin));
389
390 PetscCall(MaybeLogStagePop());
391 PetscCall(MaybeLogStagePush(stage[4]));
392
393 PetscCall(MaybeLogStagePop());
394 PetscCall(MaybeLogStagePush(stage[5]));
395
396 /* 3rd solve, test KSP[PC]SetReusePreconditioner */
397 PetscCall(MatScale(Amat, 1. + PETSC_MACHINE_EPSILON));
398 PetscCall(KSPSetReusePreconditioner(ksp, PETSC_TRUE));
399 PetscCall(KSPSolve(ksp, bb, xx));
400
401 PetscCall(MaybeLogStagePop());
402
403 PetscCall(VecNorm(bb, NORM_2, &norm2));
404
405 PetscCall(VecDuplicate(xx, &res));
406 PetscCall(MatMult(Amat, xx, res));
407 PetscCall(VecAXPY(bb, -1.0, res));
408 PetscCall(VecDestroy(&res));
409 PetscCall(VecNorm(bb, NORM_2, &norm));
410 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "[%d]%s |b-Ax|/|b|=%e, |b|=%e, emax=%e\n", 0, PETSC_FUNCTION_NAME, (double)(norm / norm2), (double)norm2, (double)emax));
411 }
412
413 /* Free work space */
414 PetscCall(KSPDestroy(&ksp));
415 PetscCall(VecDestroy(&xx));
416 PetscCall(VecDestroy(&bb));
417 PetscCall(MatDestroy(&Amat));
418 PetscCall(PetscFree(coords));
419
420 PetscCall(PetscFinalize());
421 return 0;
422 }
423
424 /* Data was previously provided in the file data/elem_3d_elast_v_25.tx */
elem_3d_elast_v_25(PetscScalar * dd)425 PetscErrorCode elem_3d_elast_v_25(PetscScalar *dd)
426 {
427 PetscScalar DD[] = {
428 0.18981481481481474, 5.27777777777777568E-002, 5.27777777777777568E-002, -5.64814814814814659E-002, -1.38888888888889072E-002, -1.38888888888889089E-002, -8.24074074074073876E-002, -5.27777777777777429E-002, 1.38888888888888725E-002,
429 4.90740740740740339E-002, 1.38888888888889124E-002, 4.72222222222222071E-002, 4.90740740740740339E-002, 4.72222222222221932E-002, 1.38888888888888968E-002, -8.24074074074073876E-002, 1.38888888888888673E-002, -5.27777777777777429E-002,
430 -7.87037037037036785E-002, -4.72222222222221932E-002, -4.72222222222222071E-002, 1.20370370370370180E-002, -1.38888888888888742E-002, -1.38888888888888829E-002, 5.27777777777777568E-002, 0.18981481481481474, 5.27777777777777568E-002,
431 1.38888888888889124E-002, 4.90740740740740269E-002, 4.72222222222221932E-002, -5.27777777777777637E-002, -8.24074074074073876E-002, 1.38888888888888725E-002, -1.38888888888889037E-002, -5.64814814814814728E-002, -1.38888888888888985E-002,
432 4.72222222222221932E-002, 4.90740740740740478E-002, 1.38888888888888968E-002, -1.38888888888888673E-002, 1.20370370370370058E-002, -1.38888888888888742E-002, -4.72222222222221932E-002, -7.87037037037036785E-002, -4.72222222222222002E-002,
433 1.38888888888888742E-002, -8.24074074074073598E-002, -5.27777777777777568E-002, 5.27777777777777568E-002, 5.27777777777777568E-002, 0.18981481481481474, 1.38888888888889055E-002, 4.72222222222222002E-002, 4.90740740740740269E-002,
434 -1.38888888888888829E-002, -1.38888888888888829E-002, 1.20370370370370180E-002, 4.72222222222222002E-002, 1.38888888888888985E-002, 4.90740740740740339E-002, -1.38888888888888985E-002, -1.38888888888888968E-002, -5.64814814814814520E-002,
435 -5.27777777777777568E-002, 1.38888888888888777E-002, -8.24074074074073876E-002, -4.72222222222222002E-002, -4.72222222222221932E-002, -7.87037037037036646E-002, 1.38888888888888794E-002, -5.27777777777777568E-002, -8.24074074074073598E-002,
436 -5.64814814814814659E-002, 1.38888888888889124E-002, 1.38888888888889055E-002, 0.18981481481481474, -5.27777777777777568E-002, -5.27777777777777499E-002, 4.90740740740740269E-002, -1.38888888888889072E-002, -4.72222222222221932E-002,
437 -8.24074074074073876E-002, 5.27777777777777568E-002, -1.38888888888888812E-002, -8.24074074074073876E-002, -1.38888888888888742E-002, 5.27777777777777499E-002, 4.90740740740740269E-002, -4.72222222222221863E-002, -1.38888888888889089E-002,
438 1.20370370370370162E-002, 1.38888888888888673E-002, 1.38888888888888742E-002, -7.87037037037036785E-002, 4.72222222222222002E-002, 4.72222222222222071E-002, -1.38888888888889072E-002, 4.90740740740740269E-002, 4.72222222222222002E-002,
439 -5.27777777777777568E-002, 0.18981481481481480, 5.27777777777777568E-002, 1.38888888888889020E-002, -5.64814814814814728E-002, -1.38888888888888951E-002, 5.27777777777777637E-002, -8.24074074074073876E-002, 1.38888888888888881E-002,
440 1.38888888888888742E-002, 1.20370370370370232E-002, -1.38888888888888812E-002, -4.72222222222221863E-002, 4.90740740740740339E-002, 1.38888888888888933E-002, -1.38888888888888812E-002, -8.24074074074073876E-002, -5.27777777777777568E-002,
441 4.72222222222222071E-002, -7.87037037037036924E-002, -4.72222222222222140E-002, -1.38888888888889089E-002, 4.72222222222221932E-002, 4.90740740740740269E-002, -5.27777777777777499E-002, 5.27777777777777568E-002, 0.18981481481481477,
442 -4.72222222222222071E-002, 1.38888888888888968E-002, 4.90740740740740131E-002, 1.38888888888888812E-002, -1.38888888888888708E-002, 1.20370370370370267E-002, 5.27777777777777568E-002, 1.38888888888888812E-002, -8.24074074074073876E-002,
443 1.38888888888889124E-002, -1.38888888888889055E-002, -5.64814814814814589E-002, -1.38888888888888812E-002, -5.27777777777777568E-002, -8.24074074074073737E-002, 4.72222222222222002E-002, -4.72222222222222002E-002, -7.87037037037036924E-002,
444 -8.24074074074073876E-002, -5.27777777777777637E-002, -1.38888888888888829E-002, 4.90740740740740269E-002, 1.38888888888889020E-002, -4.72222222222222071E-002, 0.18981481481481480, 5.27777777777777637E-002, -5.27777777777777637E-002,
445 -5.64814814814814728E-002, -1.38888888888889037E-002, 1.38888888888888951E-002, -7.87037037037036785E-002, -4.72222222222222002E-002, 4.72222222222221932E-002, 1.20370370370370128E-002, -1.38888888888888725E-002, 1.38888888888888812E-002,
446 4.90740740740740408E-002, 4.72222222222222002E-002, -1.38888888888888951E-002, -8.24074074074073876E-002, 1.38888888888888812E-002, 5.27777777777777637E-002, -5.27777777777777429E-002, -8.24074074074073876E-002, -1.38888888888888829E-002,
447 -1.38888888888889072E-002, -5.64814814814814728E-002, 1.38888888888888968E-002, 5.27777777777777637E-002, 0.18981481481481480, -5.27777777777777568E-002, 1.38888888888888916E-002, 4.90740740740740339E-002, -4.72222222222222210E-002,
448 -4.72222222222221932E-002, -7.87037037037036924E-002, 4.72222222222222002E-002, 1.38888888888888742E-002, -8.24074074074073876E-002, 5.27777777777777429E-002, 4.72222222222222002E-002, 4.90740740740740269E-002, -1.38888888888888951E-002,
449 -1.38888888888888846E-002, 1.20370370370370267E-002, 1.38888888888888916E-002, 1.38888888888888725E-002, 1.38888888888888725E-002, 1.20370370370370180E-002, -4.72222222222221932E-002, -1.38888888888888951E-002, 4.90740740740740131E-002,
450 -5.27777777777777637E-002, -5.27777777777777568E-002, 0.18981481481481480, -1.38888888888888968E-002, -4.72222222222221932E-002, 4.90740740740740339E-002, 4.72222222222221932E-002, 4.72222222222222071E-002, -7.87037037037036646E-002,
451 -1.38888888888888742E-002, 5.27777777777777499E-002, -8.24074074074073737E-002, 1.38888888888888933E-002, 1.38888888888889020E-002, -5.64814814814814589E-002, 5.27777777777777568E-002, -1.38888888888888794E-002, -8.24074074074073876E-002,
452 4.90740740740740339E-002, -1.38888888888889037E-002, 4.72222222222222002E-002, -8.24074074074073876E-002, 5.27777777777777637E-002, 1.38888888888888812E-002, -5.64814814814814728E-002, 1.38888888888888916E-002, -1.38888888888888968E-002,
453 0.18981481481481480, -5.27777777777777499E-002, 5.27777777777777707E-002, 1.20370370370370180E-002, 1.38888888888888812E-002, -1.38888888888888812E-002, -7.87037037037036785E-002, 4.72222222222222002E-002, -4.72222222222222071E-002,
454 -8.24074074074073876E-002, -1.38888888888888742E-002, -5.27777777777777568E-002, 4.90740740740740616E-002, -4.72222222222222002E-002, 1.38888888888888846E-002, 1.38888888888889124E-002, -5.64814814814814728E-002, 1.38888888888888985E-002,
455 5.27777777777777568E-002, -8.24074074074073876E-002, -1.38888888888888708E-002, -1.38888888888889037E-002, 4.90740740740740339E-002, -4.72222222222221932E-002, -5.27777777777777499E-002, 0.18981481481481480, -5.27777777777777568E-002,
456 -1.38888888888888673E-002, -8.24074074074073598E-002, 5.27777777777777429E-002, 4.72222222222222002E-002, -7.87037037037036785E-002, 4.72222222222222002E-002, 1.38888888888888708E-002, 1.20370370370370128E-002, 1.38888888888888760E-002,
457 -4.72222222222222002E-002, 4.90740740740740478E-002, -1.38888888888888951E-002, 4.72222222222222071E-002, -1.38888888888888985E-002, 4.90740740740740339E-002, -1.38888888888888812E-002, 1.38888888888888881E-002, 1.20370370370370267E-002,
458 1.38888888888888951E-002, -4.72222222222222210E-002, 4.90740740740740339E-002, 5.27777777777777707E-002, -5.27777777777777568E-002, 0.18981481481481477, 1.38888888888888829E-002, 5.27777777777777707E-002, -8.24074074074073598E-002,
459 -4.72222222222222140E-002, 4.72222222222222140E-002, -7.87037037037036646E-002, -5.27777777777777707E-002, -1.38888888888888829E-002, -8.24074074074073876E-002, -1.38888888888888881E-002, 1.38888888888888881E-002, -5.64814814814814589E-002,
460 4.90740740740740339E-002, 4.72222222222221932E-002, -1.38888888888888985E-002, -8.24074074074073876E-002, 1.38888888888888742E-002, 5.27777777777777568E-002, -7.87037037037036785E-002, -4.72222222222221932E-002, 4.72222222222221932E-002,
461 1.20370370370370180E-002, -1.38888888888888673E-002, 1.38888888888888829E-002, 0.18981481481481469, 5.27777777777777429E-002, -5.27777777777777429E-002, -5.64814814814814659E-002, -1.38888888888889055E-002, 1.38888888888889055E-002,
462 -8.24074074074074153E-002, -5.27777777777777429E-002, -1.38888888888888760E-002, 4.90740740740740408E-002, 1.38888888888888968E-002, -4.72222222222222071E-002, 4.72222222222221932E-002, 4.90740740740740478E-002, -1.38888888888888968E-002,
463 -1.38888888888888742E-002, 1.20370370370370232E-002, 1.38888888888888812E-002, -4.72222222222222002E-002, -7.87037037037036924E-002, 4.72222222222222071E-002, 1.38888888888888812E-002, -8.24074074074073598E-002, 5.27777777777777707E-002,
464 5.27777777777777429E-002, 0.18981481481481477, -5.27777777777777499E-002, 1.38888888888889107E-002, 4.90740740740740478E-002, -4.72222222222221932E-002, -5.27777777777777568E-002, -8.24074074074074153E-002, -1.38888888888888812E-002,
465 -1.38888888888888846E-002, -5.64814814814814659E-002, 1.38888888888888812E-002, 1.38888888888888968E-002, 1.38888888888888968E-002, -5.64814814814814520E-002, 5.27777777777777499E-002, -1.38888888888888812E-002, -8.24074074074073876E-002,
466 4.72222222222221932E-002, 4.72222222222222002E-002, -7.87037037037036646E-002, -1.38888888888888812E-002, 5.27777777777777429E-002, -8.24074074074073598E-002, -5.27777777777777429E-002, -5.27777777777777499E-002, 0.18981481481481474,
467 -1.38888888888888985E-002, -4.72222222222221863E-002, 4.90740740740740339E-002, 1.38888888888888829E-002, 1.38888888888888777E-002, 1.20370370370370249E-002, -4.72222222222222002E-002, -1.38888888888888933E-002, 4.90740740740740339E-002,
468 -8.24074074074073876E-002, -1.38888888888888673E-002, -5.27777777777777568E-002, 4.90740740740740269E-002, -4.72222222222221863E-002, 1.38888888888889124E-002, 1.20370370370370128E-002, 1.38888888888888742E-002, -1.38888888888888742E-002,
469 -7.87037037037036785E-002, 4.72222222222222002E-002, -4.72222222222222140E-002, -5.64814814814814659E-002, 1.38888888888889107E-002, -1.38888888888888985E-002, 0.18981481481481474, -5.27777777777777499E-002, 5.27777777777777499E-002,
470 4.90740740740740339E-002, -1.38888888888889055E-002, 4.72222222222221932E-002, -8.24074074074074153E-002, 5.27777777777777499E-002, 1.38888888888888829E-002, 1.38888888888888673E-002, 1.20370370370370058E-002, 1.38888888888888777E-002,
471 -4.72222222222221863E-002, 4.90740740740740339E-002, -1.38888888888889055E-002, -1.38888888888888725E-002, -8.24074074074073876E-002, 5.27777777777777499E-002, 4.72222222222222002E-002, -7.87037037037036785E-002, 4.72222222222222140E-002,
472 -1.38888888888889055E-002, 4.90740740740740478E-002, -4.72222222222221863E-002, -5.27777777777777499E-002, 0.18981481481481469, -5.27777777777777499E-002, 1.38888888888889072E-002, -5.64814814814814659E-002, 1.38888888888889003E-002,
473 5.27777777777777429E-002, -8.24074074074074153E-002, -1.38888888888888812E-002, -5.27777777777777429E-002, -1.38888888888888742E-002, -8.24074074074073876E-002, -1.38888888888889089E-002, 1.38888888888888933E-002, -5.64814814814814589E-002,
474 1.38888888888888812E-002, 5.27777777777777429E-002, -8.24074074074073737E-002, -4.72222222222222071E-002, 4.72222222222222002E-002, -7.87037037037036646E-002, 1.38888888888889055E-002, -4.72222222222221932E-002, 4.90740740740740339E-002,
475 5.27777777777777499E-002, -5.27777777777777499E-002, 0.18981481481481474, 4.72222222222222002E-002, -1.38888888888888985E-002, 4.90740740740740339E-002, -1.38888888888888846E-002, 1.38888888888888812E-002, 1.20370370370370284E-002,
476 -7.87037037037036785E-002, -4.72222222222221932E-002, -4.72222222222222002E-002, 1.20370370370370162E-002, -1.38888888888888812E-002, -1.38888888888888812E-002, 4.90740740740740408E-002, 4.72222222222222002E-002, 1.38888888888888933E-002,
477 -8.24074074074073876E-002, 1.38888888888888708E-002, -5.27777777777777707E-002, -8.24074074074074153E-002, -5.27777777777777568E-002, 1.38888888888888829E-002, 4.90740740740740339E-002, 1.38888888888889072E-002, 4.72222222222222002E-002,
478 0.18981481481481477, 5.27777777777777429E-002, 5.27777777777777568E-002, -5.64814814814814659E-002, -1.38888888888888846E-002, -1.38888888888888881E-002, -4.72222222222221932E-002, -7.87037037037036785E-002, -4.72222222222221932E-002,
479 1.38888888888888673E-002, -8.24074074074073876E-002, -5.27777777777777568E-002, 4.72222222222222002E-002, 4.90740740740740269E-002, 1.38888888888889020E-002, -1.38888888888888742E-002, 1.20370370370370128E-002, -1.38888888888888829E-002,
480 -5.27777777777777429E-002, -8.24074074074074153E-002, 1.38888888888888777E-002, -1.38888888888889055E-002, -5.64814814814814659E-002, -1.38888888888888985E-002, 5.27777777777777429E-002, 0.18981481481481469, 5.27777777777777429E-002,
481 1.38888888888888933E-002, 4.90740740740740339E-002, 4.72222222222222071E-002, -4.72222222222222071E-002, -4.72222222222222002E-002, -7.87037037037036646E-002, 1.38888888888888742E-002, -5.27777777777777568E-002, -8.24074074074073737E-002,
482 -1.38888888888888951E-002, -1.38888888888888951E-002, -5.64814814814814589E-002, -5.27777777777777568E-002, 1.38888888888888760E-002, -8.24074074074073876E-002, -1.38888888888888760E-002, -1.38888888888888812E-002, 1.20370370370370249E-002,
483 4.72222222222221932E-002, 1.38888888888889003E-002, 4.90740740740740339E-002, 5.27777777777777568E-002, 5.27777777777777429E-002, 0.18981481481481474, 1.38888888888888933E-002, 4.72222222222222071E-002, 4.90740740740740339E-002,
484 1.20370370370370180E-002, 1.38888888888888742E-002, 1.38888888888888794E-002, -7.87037037037036785E-002, 4.72222222222222071E-002, 4.72222222222222002E-002, -8.24074074074073876E-002, -1.38888888888888846E-002, 5.27777777777777568E-002,
485 4.90740740740740616E-002, -4.72222222222222002E-002, -1.38888888888888881E-002, 4.90740740740740408E-002, -1.38888888888888846E-002, -4.72222222222222002E-002, -8.24074074074074153E-002, 5.27777777777777429E-002, -1.38888888888888846E-002,
486 -5.64814814814814659E-002, 1.38888888888888933E-002, 1.38888888888888933E-002, 0.18981481481481477, -5.27777777777777568E-002, -5.27777777777777637E-002, -1.38888888888888742E-002, -8.24074074074073598E-002, -5.27777777777777568E-002,
487 4.72222222222222002E-002, -7.87037037037036924E-002, -4.72222222222222002E-002, 1.38888888888888812E-002, 1.20370370370370267E-002, -1.38888888888888794E-002, -4.72222222222222002E-002, 4.90740740740740478E-002, 1.38888888888888881E-002,
488 1.38888888888888968E-002, -5.64814814814814659E-002, -1.38888888888888933E-002, 5.27777777777777499E-002, -8.24074074074074153E-002, 1.38888888888888812E-002, -1.38888888888888846E-002, 4.90740740740740339E-002, 4.72222222222222071E-002,
489 -5.27777777777777568E-002, 0.18981481481481477, 5.27777777777777637E-002, -1.38888888888888829E-002, -5.27777777777777568E-002, -8.24074074074073598E-002, 4.72222222222222071E-002, -4.72222222222222140E-002, -7.87037037037036924E-002,
490 5.27777777777777637E-002, 1.38888888888888916E-002, -8.24074074074073876E-002, 1.38888888888888846E-002, -1.38888888888888951E-002, -5.64814814814814589E-002, -4.72222222222222071E-002, 1.38888888888888812E-002, 4.90740740740740339E-002,
491 1.38888888888888829E-002, -1.38888888888888812E-002, 1.20370370370370284E-002, -1.38888888888888881E-002, 4.72222222222222071E-002, 4.90740740740740339E-002, -5.27777777777777637E-002, 5.27777777777777637E-002, 0.18981481481481477,
492 };
493
494 PetscFunctionBeginUser;
495 PetscCall(PetscArraycpy(dd, DD, 576));
496 PetscFunctionReturn(PETSC_SUCCESS);
497 }
498
499 /*TEST
500
501 testset:
502 requires: !complex
503 args: -ne 11 -alpha 1.e-3 -ksp_type cg -pc_type gamg -pc_gamg_agg_nsmooths 1 -two_solves -ksp_converged_reason -use_mat_nearnullspace -mg_levels_ksp_max_it 1 -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_esteig 0,0.2,0,1.05 -mg_levels_sub_pc_type lu -pc_gamg_asm_use_agg -mg_levels_pc_asm_overlap 0 -pc_gamg_parallel_coarse_grid_solver -mg_coarse_pc_type jacobi -mg_coarse_ksp_type cg -pc_gamg_mat_coarsen_type hem -pc_gamg_mat_coarsen_max_it 5 -ksp_rtol 1e-4 -ksp_norm_type unpreconditioned -pc_gamg_threshold .001 -pc_gamg_mat_coarsen_strength_index 1,2
504 test:
505 suffix: 1
506 nsize: 1
507 filter: sed -e "s/Linear solve converged due to CONVERGED_RTOL iterations 15/Linear solve converged due to CONVERGED_RTOL iterations 14/g"
508 test:
509 suffix: 2
510 nsize: 8
511 filter: sed -e "s/Linear solve converged due to CONVERGED_RTOL iterations 1[3|4]/Linear solve converged due to CONVERGED_RTOL iterations 15/g"
512
513 testset:
514 nsize: 8
515 args: -ne 15 -alpha 1.e-3 -ksp_type cg -ksp_converged_reason -use_mat_nearnullspace -ksp_rtol 1e-4 -ksp_norm_type unpreconditioned -two_solves
516 test:
517 requires: hypre !complex !defined(PETSC_HAVE_HYPRE_DEVICE)
518 suffix: hypre
519 args: -pc_type hypre -pc_hypre_boomeramg_relax_type_all l1scaled-Jacobi
520 test:
521 suffix: gamg
522 args: -pc_type gamg -mg_levels_ksp_type richardson -mg_levels_pc_type jacobi -mg_levels_pc_jacobi_type rowl1 -mg_levels_pc_jacobi_rowl1_scale .5 -mg_levels_pc_jacobi_fixdiagonal
523 test:
524 suffix: baij
525 filter: grep -v variant
526 args: -pc_type jacobi -pc_jacobi_type rowl1 -ksp_type cg -mat_type baij -ksp_view -ksp_rtol 1e-1 -two_solves false
527
528 test:
529 suffix: latebs
530 filter: grep -v variant
531 nsize: 8
532 args: -test_late_bs 0 -ne 9 -alpha 1.e-3 -ksp_type cg -pc_type gamg -pc_gamg_agg_nsmooths 1 -pc_gamg_reuse_interpolation true -two_solves -ksp_converged_reason -use_mat_nearnullspace false -mg_levels_ksp_max_it 2 -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_esteig 0,0.2,0,1.05 -pc_gamg_esteig_ksp_max_it 10 -pc_gamg_threshold -0.01 -pc_gamg_coarse_eq_limit 200 -pc_gamg_process_eq_limit 30 -pc_gamg_repartition false -pc_mg_cycle_type v -pc_gamg_parallel_coarse_grid_solver -mg_coarse_pc_type jacobi -mg_coarse_ksp_type cg -ksp_monitor_short -ksp_view -pc_gamg_injection_index 1,2 -mg_fine_ksp_type richardson -mg_fine_pc_type jacobi -mg_fine_pc_jacobi_type rowl1 -mg_fine_pc_jacobi_rowl1_scale .25
533
534 test:
535 suffix: latebs-2
536 filter: grep -v variant
537 nsize: 8
538 args: -test_late_bs -ne 9 -alpha 1.e-3 -ksp_type cg -pc_type gamg -pc_gamg_agg_nsmooths 1 -pc_gamg_reuse_interpolation true -two_solves -ksp_converged_reason -use_mat_nearnullspace -mg_levels_ksp_max_it 2 -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_esteig 0,0.2,0,1.05 -pc_gamg_esteig_ksp_max_it 10 -pc_gamg_threshold -0.01 -pc_gamg_coarse_eq_limit 200 -pc_gamg_process_eq_limit 30 -pc_gamg_repartition false -pc_mg_cycle_type v -pc_gamg_parallel_coarse_grid_solver -mg_coarse_pc_type jacobi -mg_coarse_ksp_type cg -ksp_monitor_short -ksp_view
539
540 test:
541 suffix: ml
542 nsize: 8
543 args: -ne 9 -alpha 1.e-3 -ksp_type cg -pc_type ml -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_esteig 0,0.05,0,1.05 -mg_levels_pc_type sor -ksp_monitor_short
544 requires: ml
545
546 test:
547 suffix: nns
548 args: -ne 9 -alpha 1.e-3 -ksp_converged_reason -ksp_type cg -ksp_max_it 50 -pc_type gamg -pc_gamg_esteig_ksp_type cg -pc_gamg_esteig_ksp_max_it 10 -pc_gamg_type agg -pc_gamg_agg_nsmooths 1 -pc_gamg_coarse_eq_limit 1000 -mg_levels_ksp_type chebyshev -mg_levels_pc_type sor -pc_gamg_reuse_interpolation true -two_solves -use_mat_nearnullspace -pc_gamg_use_sa_esteig 0 -mg_levels_esteig_ksp_max_it 10
549
550 test:
551 suffix: nns_telescope
552 nsize: 2
553 args: -use_mat_nearnullspace -pc_type telescope -pc_telescope_reduction_factor 2 -telescope_pc_type gamg -telescope_pc_gamg_esteig_ksp_type cg -telescope_pc_gamg_esteig_ksp_max_it 10
554 output_file: output/empty.out
555
556 test:
557 suffix: nns_gdsw
558 filter: grep -v "variant HERMITIAN"
559 nsize: 8
560 args: -use_mat_nearnullspace -ksp_monitor_short -pc_type mg -pc_mg_levels 2 -pc_mg_adapt_interp_coarse_space gdsw -pc_mg_galerkin -mg_levels_pc_type bjacobi -ne 3 -ksp_view
561
562 test:
563 suffix: seqaijmkl
564 nsize: 8
565 requires: mkl_sparse
566 args: -ne 9 -alpha 1.e-3 -ksp_type cg -pc_type gamg -pc_gamg_agg_nsmooths 1 -pc_gamg_reuse_interpolation true -two_solves -ksp_converged_reason -use_mat_nearnullspace -mg_levels_ksp_max_it 2 -mg_levels_ksp_type chebyshev -mg_levels_pc_type jacobi -mg_levels_ksp_chebyshev_esteig 0,0.2,0,1.05 -pc_gamg_esteig_ksp_max_it 10 -pc_gamg_threshold 0.01 -pc_gamg_coarse_eq_limit 2000 -pc_gamg_process_eq_limit 200 -pc_gamg_repartition false -pc_mg_cycle_type v -mat_seqaij_type seqaijmkl
567
568 testset:
569 nsize: {{1 8}separate output}
570 args: -ne 7 -pc_type gamg -rap_mg_levels_pc_type pbjacobi -rap_ksp_monitor -rap_mg_coarse_ksp_type cg -rap_mg_coarse_pc_type pbjacobi -use_mat_nearnullspace -test_rap_bs -rap_ksp_view -rap_mat_view ::ascii_info -rap_ksp_converged_reason
571 filter: grep -v "variant HERMITIAN"
572 test:
573 suffix: rap_bs
574
575 test:
576 requires: kokkos_kernels
577 suffix: rap_bs_kokkos
578 args: -mat_type aijkokkos
579
580 test:
581 requires: cuda
582 suffix: rap_bs_cuda
583 args: -mat_type aijcusparse -rap_mg_coarse_pc_type jacobi -rap_mg_levels_pc_type jacobi -rap_mg_levels_ksp_type richardson -rap_mg_levels_pc_jacobi_type rowl1 -rap_mg_levels_pc_jacobi_rowl1_scale .5
584
585 test:
586 requires: hip
587 suffix: rap_bs_hip
588 args: -mat_type aijhipsparse -rap_mg_coarse_pc_type jacobi -rap_mg_levels_pc_type jacobi -rap_mg_levels_ksp_type richardson -rap_mg_levels_pc_jacobi_type rowl1 -rap_mg_levels_pc_jacobi_rowl1_scale .5
589
590 TEST*/
591