Actual source code: ex28.c
1: static const char help[] = "1D multiphysics prototype with analytic Jacobians to solve individual problems and a coupled problem.\n\n";
3: /* Solve a PDE coupled to an algebraic system in 1D
4: *
5: * PDE (U):
6: * -(k u_x)_x = 1 on (0,1), subject to u(0) = 0, u(1) = 1
7: * Algebraic (K):
8: * exp(k-1) + k = 1/(1/(1+u) + 1/(1+u_x^2))
9: *
10: * The discretization places k at staggered points, and a separate DMDA is used for each "physics".
11: *
12: * This example is a prototype for coupling in multi-physics problems, therefore residual evaluation and assembly for
13: * each problem (referred to as U and K) are written separately. This permits the same "physics" code to be used for
14: * solving each uncoupled problem as well as the coupled system. In particular, run with -problem_type 0 to solve only
15: * problem U (with K fixed), -problem_type 1 to solve only K (with U fixed), and -problem_type 2 to solve both at once.
16: *
17: * In all cases, a fully-assembled analytic Jacobian is available, so the systems can be solved with a direct solve or
18: * any other standard method. Additionally, by running with
19: *
20: * -pack_dm_mat_type nest
21: *
22: * The same code assembles a coupled matrix where each block is stored separately, which allows the use of PCFieldSplit
23: * without copying values to extract submatrices.
24: */
26: #include <petscsnes.h>
27: #include <petscdm.h>
28: #include <petscdmda.h>
29: #include <petscdmcomposite.h>
31: typedef struct _UserCtx *User;
32: struct _UserCtx {
33: PetscInt ptype;
34: DM pack;
35: Vec Uloc, Kloc;
36: };
38: static PetscErrorCode FormFunctionLocal_U(User user, DMDALocalInfo *info, const PetscScalar u[], const PetscScalar k[], PetscScalar f[])
39: {
40: PetscReal hx = 1. / info->mx;
41: PetscInt i;
44: for (i = info->xs; i < info->xs + info->xm; i++) {
45: if (i == 0) f[i] = 1. / hx * u[i];
46: else if (i == info->mx - 1) f[i] = 1. / hx * (u[i] - 1.0);
47: else f[i] = hx * ((k[i - 1] * (u[i] - u[i - 1]) - k[i] * (u[i + 1] - u[i])) / (hx * hx) - 1.0);
48: }
49: return 0;
50: }
52: static PetscErrorCode FormFunctionLocal_K(User user, DMDALocalInfo *info, const PetscScalar u[], const PetscScalar k[], PetscScalar f[])
53: {
54: PetscReal hx = 1. / info->mx;
55: PetscInt i;
58: for (i = info->xs; i < info->xs + info->xm; i++) {
59: const PetscScalar ubar = 0.5 * (u[i + 1] + u[i]), gradu = (u[i + 1] - u[i]) / hx, g = 1. + gradu * gradu, w = 1. / (1. + ubar) + 1. / g;
60: f[i] = hx * (PetscExpScalar(k[i] - 1.0) + k[i] - 1. / w);
61: }
62: return 0;
63: }
65: static PetscErrorCode FormFunction_All(SNES snes, Vec X, Vec F, void *ctx)
66: {
67: User user = (User)ctx;
68: DM dau, dak;
69: DMDALocalInfo infou, infok;
70: PetscScalar *u, *k;
71: PetscScalar *fu, *fk;
72: Vec Uloc, Kloc, Fu, Fk;
75: DMCompositeGetEntries(user->pack, &dau, &dak);
76: DMDAGetLocalInfo(dau, &infou);
77: DMDAGetLocalInfo(dak, &infok);
78: DMCompositeGetLocalVectors(user->pack, &Uloc, &Kloc);
79: switch (user->ptype) {
80: case 0:
81: DMGlobalToLocalBegin(dau, X, INSERT_VALUES, Uloc);
82: DMGlobalToLocalEnd(dau, X, INSERT_VALUES, Uloc);
83: DMDAVecGetArray(dau, Uloc, &u);
84: DMDAVecGetArray(dak, user->Kloc, &k);
85: DMDAVecGetArray(dau, F, &fu);
86: FormFunctionLocal_U(user, &infou, u, k, fu);
87: DMDAVecRestoreArray(dau, F, &fu);
88: DMDAVecRestoreArray(dau, Uloc, &u);
89: DMDAVecRestoreArray(dak, user->Kloc, &k);
90: break;
91: case 1:
92: DMGlobalToLocalBegin(dak, X, INSERT_VALUES, Kloc);
93: DMGlobalToLocalEnd(dak, X, INSERT_VALUES, Kloc);
94: DMDAVecGetArray(dau, user->Uloc, &u);
95: DMDAVecGetArray(dak, Kloc, &k);
96: DMDAVecGetArray(dak, F, &fk);
97: FormFunctionLocal_K(user, &infok, u, k, fk);
98: DMDAVecRestoreArray(dak, F, &fk);
99: DMDAVecRestoreArray(dau, user->Uloc, &u);
100: DMDAVecRestoreArray(dak, Kloc, &k);
101: break;
102: case 2:
103: DMCompositeScatter(user->pack, X, Uloc, Kloc);
104: DMDAVecGetArray(dau, Uloc, &u);
105: DMDAVecGetArray(dak, Kloc, &k);
106: DMCompositeGetAccess(user->pack, F, &Fu, &Fk);
107: DMDAVecGetArray(dau, Fu, &fu);
108: DMDAVecGetArray(dak, Fk, &fk);
109: FormFunctionLocal_U(user, &infou, u, k, fu);
110: FormFunctionLocal_K(user, &infok, u, k, fk);
111: DMDAVecRestoreArray(dau, Fu, &fu);
112: DMDAVecRestoreArray(dak, Fk, &fk);
113: DMCompositeRestoreAccess(user->pack, F, &Fu, &Fk);
114: DMDAVecRestoreArray(dau, Uloc, &u);
115: DMDAVecRestoreArray(dak, Kloc, &k);
116: break;
117: }
118: DMCompositeRestoreLocalVectors(user->pack, &Uloc, &Kloc);
119: return 0;
120: }
122: static PetscErrorCode FormJacobianLocal_U(User user, DMDALocalInfo *info, const PetscScalar u[], const PetscScalar k[], Mat Buu)
123: {
124: PetscReal hx = 1. / info->mx;
125: PetscInt i;
128: for (i = info->xs; i < info->xs + info->xm; i++) {
129: PetscInt row = i - info->gxs, cols[] = {row - 1, row, row + 1};
130: PetscScalar val = 1. / hx;
131: if (i == 0) {
132: MatSetValuesLocal(Buu, 1, &row, 1, &row, &val, INSERT_VALUES);
133: } else if (i == info->mx - 1) {
134: MatSetValuesLocal(Buu, 1, &row, 1, &row, &val, INSERT_VALUES);
135: } else {
136: PetscScalar vals[] = {-k[i - 1] / hx, (k[i - 1] + k[i]) / hx, -k[i] / hx};
137: MatSetValuesLocal(Buu, 1, &row, 3, cols, vals, INSERT_VALUES);
138: }
139: }
140: return 0;
141: }
143: static PetscErrorCode FormJacobianLocal_K(User user, DMDALocalInfo *info, const PetscScalar u[], const PetscScalar k[], Mat Bkk)
144: {
145: PetscReal hx = 1. / info->mx;
146: PetscInt i;
149: for (i = info->xs; i < info->xs + info->xm; i++) {
150: PetscInt row = i - info->gxs;
151: PetscScalar vals[] = {hx * (PetscExpScalar(k[i] - 1.) + (PetscScalar)1.)};
152: MatSetValuesLocal(Bkk, 1, &row, 1, &row, vals, INSERT_VALUES);
153: }
154: return 0;
155: }
157: static PetscErrorCode FormJacobianLocal_UK(User user, DMDALocalInfo *info, DMDALocalInfo *infok, const PetscScalar u[], const PetscScalar k[], Mat Buk)
158: {
159: PetscReal hx = 1. / info->mx;
160: PetscInt i;
161: PetscInt row, cols[2];
162: PetscScalar vals[2];
165: if (!Buk) return 0; /* Not assembling this block */
166: for (i = info->xs; i < info->xs + info->xm; i++) {
167: if (i == 0 || i == info->mx - 1) continue;
168: row = i - info->gxs;
169: cols[0] = i - 1 - infok->gxs;
170: vals[0] = (u[i] - u[i - 1]) / hx;
171: cols[1] = i - infok->gxs;
172: vals[1] = (u[i] - u[i + 1]) / hx;
173: MatSetValuesLocal(Buk, 1, &row, 2, cols, vals, INSERT_VALUES);
174: }
175: return 0;
176: }
178: static PetscErrorCode FormJacobianLocal_KU(User user, DMDALocalInfo *info, DMDALocalInfo *infok, const PetscScalar u[], const PetscScalar k[], Mat Bku)
179: {
180: PetscInt i;
181: PetscReal hx = 1. / (info->mx - 1);
184: if (!Bku) return 0; /* Not assembling this block */
185: for (i = infok->xs; i < infok->xs + infok->xm; i++) {
186: PetscInt row = i - infok->gxs, cols[2];
187: PetscScalar vals[2];
188: const PetscScalar ubar = 0.5 * (u[i] + u[i + 1]), ubar_L = 0.5, ubar_R = 0.5, gradu = (u[i + 1] - u[i]) / hx, gradu_L = -1. / hx, gradu_R = 1. / hx, g = 1. + PetscSqr(gradu), g_gradu = 2. * gradu, w = 1. / (1. + ubar) + 1. / g, w_ubar = -1. / PetscSqr(1. + ubar), w_gradu = -g_gradu / PetscSqr(g), iw = 1. / w, iw_ubar = -w_ubar * PetscSqr(iw), iw_gradu = -w_gradu * PetscSqr(iw);
189: cols[0] = i - info->gxs;
190: vals[0] = -hx * (iw_ubar * ubar_L + iw_gradu * gradu_L);
191: cols[1] = i + 1 - info->gxs;
192: vals[1] = -hx * (iw_ubar * ubar_R + iw_gradu * gradu_R);
193: MatSetValuesLocal(Bku, 1, &row, 2, cols, vals, INSERT_VALUES);
194: }
195: return 0;
196: }
198: static PetscErrorCode FormJacobian_All(SNES snes, Vec X, Mat J, Mat B, void *ctx)
199: {
200: User user = (User)ctx;
201: DM dau, dak;
202: DMDALocalInfo infou, infok;
203: PetscScalar *u, *k;
204: Vec Uloc, Kloc;
207: DMCompositeGetEntries(user->pack, &dau, &dak);
208: DMDAGetLocalInfo(dau, &infou);
209: DMDAGetLocalInfo(dak, &infok);
210: DMCompositeGetLocalVectors(user->pack, &Uloc, &Kloc);
211: switch (user->ptype) {
212: case 0:
213: DMGlobalToLocalBegin(dau, X, INSERT_VALUES, Uloc);
214: DMGlobalToLocalEnd(dau, X, INSERT_VALUES, Uloc);
215: DMDAVecGetArray(dau, Uloc, &u);
216: DMDAVecGetArray(dak, user->Kloc, &k);
217: FormJacobianLocal_U(user, &infou, u, k, B);
218: DMDAVecRestoreArray(dau, Uloc, &u);
219: DMDAVecRestoreArray(dak, user->Kloc, &k);
220: break;
221: case 1:
222: DMGlobalToLocalBegin(dak, X, INSERT_VALUES, Kloc);
223: DMGlobalToLocalEnd(dak, X, INSERT_VALUES, Kloc);
224: DMDAVecGetArray(dau, user->Uloc, &u);
225: DMDAVecGetArray(dak, Kloc, &k);
226: FormJacobianLocal_K(user, &infok, u, k, B);
227: DMDAVecRestoreArray(dau, user->Uloc, &u);
228: DMDAVecRestoreArray(dak, Kloc, &k);
229: break;
230: case 2: {
231: Mat Buu, Buk, Bku, Bkk;
232: PetscBool nest;
233: IS *is;
234: DMCompositeScatter(user->pack, X, Uloc, Kloc);
235: DMDAVecGetArray(dau, Uloc, &u);
236: DMDAVecGetArray(dak, Kloc, &k);
237: DMCompositeGetLocalISs(user->pack, &is);
238: MatGetLocalSubMatrix(B, is[0], is[0], &Buu);
239: MatGetLocalSubMatrix(B, is[0], is[1], &Buk);
240: MatGetLocalSubMatrix(B, is[1], is[0], &Bku);
241: MatGetLocalSubMatrix(B, is[1], is[1], &Bkk);
242: FormJacobianLocal_U(user, &infou, u, k, Buu);
243: PetscObjectTypeCompare((PetscObject)B, MATNEST, &nest);
244: if (!nest) {
245: /*
246: DMCreateMatrix_Composite() for a nested matrix does not generate off-block matrices that one can call MatSetValuesLocal() on, it just creates dummy
247: matrices with no entries; there cannot insert values into them. Commit b6480e041dd2293a65f96222772d68cdb4ed6306
248: changed Mat_Nest() from returning NULL pointers for these submatrices to dummy matrices because PCFIELDSPLIT could not
249: handle the returned null matrices.
250: */
251: FormJacobianLocal_UK(user, &infou, &infok, u, k, Buk);
252: FormJacobianLocal_KU(user, &infou, &infok, u, k, Bku);
253: }
254: FormJacobianLocal_K(user, &infok, u, k, Bkk);
255: MatRestoreLocalSubMatrix(B, is[0], is[0], &Buu);
256: MatRestoreLocalSubMatrix(B, is[0], is[1], &Buk);
257: MatRestoreLocalSubMatrix(B, is[1], is[0], &Bku);
258: MatRestoreLocalSubMatrix(B, is[1], is[1], &Bkk);
259: DMDAVecRestoreArray(dau, Uloc, &u);
260: DMDAVecRestoreArray(dak, Kloc, &k);
262: ISDestroy(&is[0]);
263: ISDestroy(&is[1]);
264: PetscFree(is);
265: } break;
266: }
267: DMCompositeRestoreLocalVectors(user->pack, &Uloc, &Kloc);
268: MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
269: MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
270: if (J != B) {
271: MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY);
272: MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY);
273: }
274: return 0;
275: }
277: static PetscErrorCode FormInitial_Coupled(User user, Vec X)
278: {
279: DM dau, dak;
280: DMDALocalInfo infou, infok;
281: Vec Xu, Xk;
282: PetscScalar *u, *k, hx;
283: PetscInt i;
286: DMCompositeGetEntries(user->pack, &dau, &dak);
287: DMCompositeGetAccess(user->pack, X, &Xu, &Xk);
288: DMDAVecGetArray(dau, Xu, &u);
289: DMDAVecGetArray(dak, Xk, &k);
290: DMDAGetLocalInfo(dau, &infou);
291: DMDAGetLocalInfo(dak, &infok);
292: hx = 1. / (infok.mx);
293: for (i = infou.xs; i < infou.xs + infou.xm; i++) u[i] = (PetscScalar)i * hx * (1. - (PetscScalar)i * hx);
294: for (i = infok.xs; i < infok.xs + infok.xm; i++) k[i] = 1.0 + 0.5 * PetscSinScalar(2 * PETSC_PI * i * hx);
295: DMDAVecRestoreArray(dau, Xu, &u);
296: DMDAVecRestoreArray(dak, Xk, &k);
297: DMCompositeRestoreAccess(user->pack, X, &Xu, &Xk);
298: DMCompositeScatter(user->pack, X, user->Uloc, user->Kloc);
299: return 0;
300: }
302: int main(int argc, char *argv[])
303: {
304: DM dau, dak, pack;
305: const PetscInt *lxu;
306: PetscInt *lxk, m, sizes;
307: User user;
308: SNES snes;
309: Vec X, F, Xu, Xk, Fu, Fk;
310: Mat B;
311: IS *isg;
312: PetscBool pass_dm;
315: PetscInitialize(&argc, &argv, 0, help);
316: DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 10, 1, 1, NULL, &dau);
317: DMSetOptionsPrefix(dau, "u_");
318: DMSetFromOptions(dau);
319: DMSetUp(dau);
320: DMDAGetOwnershipRanges(dau, &lxu, 0, 0);
321: DMDAGetInfo(dau, 0, &m, 0, 0, &sizes, 0, 0, 0, 0, 0, 0, 0, 0);
322: PetscMalloc1(sizes, &lxk);
323: PetscArraycpy(lxk, lxu, sizes);
324: lxk[0]--;
325: DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, m - 1, 1, 1, lxk, &dak);
326: DMSetOptionsPrefix(dak, "k_");
327: DMSetFromOptions(dak);
328: DMSetUp(dak);
329: PetscFree(lxk);
331: DMCompositeCreate(PETSC_COMM_WORLD, &pack);
332: DMSetOptionsPrefix(pack, "pack_");
333: DMCompositeAddDM(pack, dau);
334: DMCompositeAddDM(pack, dak);
335: DMDASetFieldName(dau, 0, "u");
336: DMDASetFieldName(dak, 0, "k");
337: DMSetFromOptions(pack);
339: DMCreateGlobalVector(pack, &X);
340: VecDuplicate(X, &F);
342: PetscNew(&user);
344: user->pack = pack;
346: DMCompositeGetGlobalISs(pack, &isg);
347: DMCompositeGetLocalVectors(pack, &user->Uloc, &user->Kloc);
348: DMCompositeScatter(pack, X, user->Uloc, user->Kloc);
350: PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "Coupled problem options", "SNES");
351: {
352: user->ptype = 0;
353: pass_dm = PETSC_TRUE;
355: PetscOptionsInt("-problem_type", "0: solve for u only, 1: solve for k only, 2: solve for both", 0, user->ptype, &user->ptype, NULL);
356: PetscOptionsBool("-pass_dm", "Pass the packed DM to SNES to use when determining splits and forward into splits", 0, pass_dm, &pass_dm, NULL);
357: }
358: PetscOptionsEnd();
360: FormInitial_Coupled(user, X);
362: SNESCreate(PETSC_COMM_WORLD, &snes);
363: switch (user->ptype) {
364: case 0:
365: DMCompositeGetAccess(pack, X, &Xu, 0);
366: DMCompositeGetAccess(pack, F, &Fu, 0);
367: DMCreateMatrix(dau, &B);
368: SNESSetFunction(snes, Fu, FormFunction_All, user);
369: SNESSetJacobian(snes, B, B, FormJacobian_All, user);
370: SNESSetFromOptions(snes);
371: SNESSetDM(snes, dau);
372: SNESSolve(snes, NULL, Xu);
373: DMCompositeRestoreAccess(pack, X, &Xu, 0);
374: DMCompositeRestoreAccess(pack, F, &Fu, 0);
375: break;
376: case 1:
377: DMCompositeGetAccess(pack, X, 0, &Xk);
378: DMCompositeGetAccess(pack, F, 0, &Fk);
379: DMCreateMatrix(dak, &B);
380: SNESSetFunction(snes, Fk, FormFunction_All, user);
381: SNESSetJacobian(snes, B, B, FormJacobian_All, user);
382: SNESSetFromOptions(snes);
383: SNESSetDM(snes, dak);
384: SNESSolve(snes, NULL, Xk);
385: DMCompositeRestoreAccess(pack, X, 0, &Xk);
386: DMCompositeRestoreAccess(pack, F, 0, &Fk);
387: break;
388: case 2:
389: DMCreateMatrix(pack, &B);
390: /* This example does not correctly allocate off-diagonal blocks. These options allows new nonzeros (slow). */
391: MatSetOption(B, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_FALSE);
392: MatSetOption(B, MAT_NEW_NONZERO_ALLOCATION_ERR, PETSC_FALSE);
393: SNESSetFunction(snes, F, FormFunction_All, user);
394: SNESSetJacobian(snes, B, B, FormJacobian_All, user);
395: SNESSetFromOptions(snes);
396: if (!pass_dm) { /* Manually provide index sets and names for the splits */
397: KSP ksp;
398: PC pc;
399: SNESGetKSP(snes, &ksp);
400: KSPGetPC(ksp, &pc);
401: PCFieldSplitSetIS(pc, "u", isg[0]);
402: PCFieldSplitSetIS(pc, "k", isg[1]);
403: } else {
404: /* The same names come from the options prefix for dau and dak. This option can support geometric multigrid inside
405: * of splits, but it requires using a DM (perhaps your own implementation). */
406: SNESSetDM(snes, pack);
407: }
408: SNESSolve(snes, NULL, X);
409: break;
410: }
411: VecViewFromOptions(X, NULL, "-view_sol");
413: if (0) {
414: PetscInt col = 0;
415: PetscBool mult_dup = PETSC_FALSE, view_dup = PETSC_FALSE;
416: Mat D;
417: Vec Y;
419: PetscOptionsGetInt(NULL, 0, "-col", &col, 0);
420: PetscOptionsGetBool(NULL, 0, "-mult_dup", &mult_dup, 0);
421: PetscOptionsGetBool(NULL, 0, "-view_dup", &view_dup, 0);
423: VecDuplicate(X, &Y);
424: /* MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY); */
425: /* MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY); */
426: MatConvert(B, MATAIJ, MAT_INITIAL_MATRIX, &D);
427: VecZeroEntries(X);
428: VecSetValue(X, col, 1.0, INSERT_VALUES);
429: VecAssemblyBegin(X);
430: VecAssemblyEnd(X);
431: MatMult(mult_dup ? D : B, X, Y);
432: MatView(view_dup ? D : B, PETSC_VIEWER_STDOUT_WORLD);
433: /* VecView(X,PETSC_VIEWER_STDOUT_WORLD); */
434: VecView(Y, PETSC_VIEWER_STDOUT_WORLD);
435: MatDestroy(&D);
436: VecDestroy(&Y);
437: }
439: DMCompositeRestoreLocalVectors(pack, &user->Uloc, &user->Kloc);
440: PetscFree(user);
442: ISDestroy(&isg[0]);
443: ISDestroy(&isg[1]);
444: PetscFree(isg);
445: VecDestroy(&X);
446: VecDestroy(&F);
447: MatDestroy(&B);
448: DMDestroy(&dau);
449: DMDestroy(&dak);
450: DMDestroy(&pack);
451: SNESDestroy(&snes);
452: PetscFinalize();
453: return 0;
454: }
456: /*TEST
458: test:
459: suffix: 0
460: nsize: 3
461: args: -u_da_grid_x 20 -snes_converged_reason -snes_monitor_short -problem_type 0
463: test:
464: suffix: 1
465: nsize: 3
466: args: -u_da_grid_x 20 -snes_converged_reason -snes_monitor_short -problem_type 1
468: test:
469: suffix: 2
470: nsize: 3
471: args: -u_da_grid_x 20 -snes_converged_reason -snes_monitor_short -problem_type 2
473: testset:
474: suffix: 3
475: nsize: 3
476: output_file: output/ex28_3.out
477: args: -u_da_grid_x 20 -snes_converged_reason -snes_monitor_short -ksp_monitor_short -problem_type 2 -snes_mf_operator -pc_type fieldsplit -pc_fieldsplit_dm_splits -pc_fieldsplit_type additive -fieldsplit_u_ksp_type gmres -fieldsplit_k_pc_type jacobi
479: test:
480: suffix: std
481: args: -pack_dm_mat_type {{aij nest}}
483: test:
484: suffix: cuda
485: requires: cuda
486: args: -pack_dm_mat_type aijcusparse -pack_dm_vec_type cuda
488: test:
489: suffix: hip
490: requires: hip
491: args: -pack_dm_mat_type aijhipsparse -pack_dm_vec_type hip
493: test:
494: suffix: kok
495: requires: kokkos_kernels
496: args: -pack_dm_mat_type aijkokkos -pack_dm_vec_type kokkos
498: test:
499: suffix: 4
500: nsize: 6
501: args: -u_da_grid_x 257 -snes_converged_reason -snes_monitor_short -ksp_monitor_short -problem_type 2 -snes_mf_operator -pack_dm_mat_type aij -pc_type fieldsplit -pc_fieldsplit_type multiplicative -fieldsplit_u_ksp_type gmres -fieldsplit_u_ksp_pc_side right -fieldsplit_u_pc_type mg -fieldsplit_u_pc_mg_levels 4 -fieldsplit_u_mg_levels_ksp_type richardson -fieldsplit_u_mg_levels_ksp_max_it 1 -fieldsplit_u_mg_levels_pc_type sor -fieldsplit_u_pc_mg_galerkin pmat -fieldsplit_u_ksp_converged_reason -fieldsplit_k_pc_type jacobi
502: requires: !single
504: test:
505: suffix: glvis_composite_da_1d
506: args: -u_da_grid_x 20 -problem_type 0 -snes_monitor_solution glvis:
508: test:
509: suffix: cuda
510: nsize: 1
511: requires: cuda
512: args: -u_da_grid_x 20 -snes_converged_reason -snes_monitor_short -problem_type 2 -pack_dm_mat_type aijcusparse -pack_dm_vec_type cuda
514: test:
515: suffix: viennacl
516: nsize: 1
517: requires: viennacl
518: args: -u_da_grid_x 20 -snes_converged_reason -snes_monitor_short -problem_type 2 -pack_dm_mat_type aijviennacl -pack_dm_vec_type viennacl
520: TEST*/