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
 60:       ubar  = 0.5*(u[i+1]+u[i]),
 61:       gradu = (u[i+1]-u[i])/hx,
 62:       g     = 1. + gradu*gradu,
 63:       w     = 1./(1.+ubar) + 1./g;
 64:     f[i] = hx*(PetscExpScalar(k[i]-1.0) + k[i] - 1./w);
 65:   }
 66:   return 0;
 67: }

 69: static PetscErrorCode FormFunction_All(SNES snes,Vec X,Vec F,void *ctx)
 70: {
 71:   User           user = (User)ctx;
 72:   DM             dau,dak;
 73:   DMDALocalInfo  infou,infok;
 74:   PetscScalar    *u,*k;
 75:   PetscScalar    *fu,*fk;
 76:   Vec            Uloc,Kloc,Fu,Fk;

 79:   DMCompositeGetEntries(user->pack,&dau,&dak);
 80:   DMDAGetLocalInfo(dau,&infou);
 81:   DMDAGetLocalInfo(dak,&infok);
 82:   DMCompositeGetLocalVectors(user->pack,&Uloc,&Kloc);
 83:   switch (user->ptype) {
 84:   case 0:
 85:     DMGlobalToLocalBegin(dau,X,INSERT_VALUES,Uloc);
 86:     DMGlobalToLocalEnd  (dau,X,INSERT_VALUES,Uloc);
 87:     DMDAVecGetArray(dau,Uloc,&u);
 88:     DMDAVecGetArray(dak,user->Kloc,&k);
 89:     DMDAVecGetArray(dau,F,&fu);
 90:     FormFunctionLocal_U(user,&infou,u,k,fu);
 91:     DMDAVecRestoreArray(dau,F,&fu);
 92:     DMDAVecRestoreArray(dau,Uloc,&u);
 93:     DMDAVecRestoreArray(dak,user->Kloc,&k);
 94:     break;
 95:   case 1:
 96:     DMGlobalToLocalBegin(dak,X,INSERT_VALUES,Kloc);
 97:     DMGlobalToLocalEnd  (dak,X,INSERT_VALUES,Kloc);
 98:     DMDAVecGetArray(dau,user->Uloc,&u);
 99:     DMDAVecGetArray(dak,Kloc,&k);
100:     DMDAVecGetArray(dak,F,&fk);
101:     FormFunctionLocal_K(user,&infok,u,k,fk);
102:     DMDAVecRestoreArray(dak,F,&fk);
103:     DMDAVecRestoreArray(dau,user->Uloc,&u);
104:     DMDAVecRestoreArray(dak,Kloc,&k);
105:     break;
106:   case 2:
107:     DMCompositeScatter(user->pack,X,Uloc,Kloc);
108:     DMDAVecGetArray(dau,Uloc,&u);
109:     DMDAVecGetArray(dak,Kloc,&k);
110:     DMCompositeGetAccess(user->pack,F,&Fu,&Fk);
111:     DMDAVecGetArray(dau,Fu,&fu);
112:     DMDAVecGetArray(dak,Fk,&fk);
113:     FormFunctionLocal_U(user,&infou,u,k,fu);
114:     FormFunctionLocal_K(user,&infok,u,k,fk);
115:     DMDAVecRestoreArray(dau,Fu,&fu);
116:     DMDAVecRestoreArray(dak,Fk,&fk);
117:     DMCompositeRestoreAccess(user->pack,F,&Fu,&Fk);
118:     DMDAVecRestoreArray(dau,Uloc,&u);
119:     DMDAVecRestoreArray(dak,Kloc,&k);
120:     break;
121:   }
122:   DMCompositeRestoreLocalVectors(user->pack,&Uloc,&Kloc);
123:   return 0;
124: }

126: static PetscErrorCode FormJacobianLocal_U(User user,DMDALocalInfo *info,const PetscScalar u[],const PetscScalar k[],Mat Buu)
127: {
128:   PetscReal      hx = 1./info->mx;
129:   PetscInt       i;

132:   for (i=info->xs; i<info->xs+info->xm; i++) {
133:     PetscInt    row = i-info->gxs,cols[] = {row-1,row,row+1};
134:     PetscScalar val = 1./hx;
135:     if (i == 0) {
136:       MatSetValuesLocal(Buu,1,&row,1,&row,&val,INSERT_VALUES);
137:     } else if (i == info->mx-1) {
138:       MatSetValuesLocal(Buu,1,&row,1,&row,&val,INSERT_VALUES);
139:     } else {
140:       PetscScalar vals[] = {-k[i-1]/hx,(k[i-1]+k[i])/hx,-k[i]/hx};
141:       MatSetValuesLocal(Buu,1,&row,3,cols,vals,INSERT_VALUES);
142:     }
143:   }
144:   return 0;
145: }

147: static PetscErrorCode FormJacobianLocal_K(User user,DMDALocalInfo *info,const PetscScalar u[],const PetscScalar k[],Mat Bkk)
148: {
149:   PetscReal      hx = 1./info->mx;
150:   PetscInt       i;

153:   for (i=info->xs; i<info->xs+info->xm; i++) {
154:     PetscInt    row    = i-info->gxs;
155:     PetscScalar vals[] = {hx*(PetscExpScalar(k[i]-1.)+ (PetscScalar)1.)};
156:     MatSetValuesLocal(Bkk,1,&row,1,&row,vals,INSERT_VALUES);
157:   }
158:   return 0;
159: }

161: static PetscErrorCode FormJacobianLocal_UK(User user,DMDALocalInfo *info,DMDALocalInfo *infok,const PetscScalar u[],const PetscScalar k[],Mat Buk)
162: {
163:   PetscReal      hx = 1./info->mx;
164:   PetscInt       i;
165:   PetscInt       row,cols[2];
166:   PetscScalar    vals[2];

169:   if (!Buk) return 0; /* Not assembling this block */
170:   for (i=info->xs; i<info->xs+info->xm; i++) {
171:     if (i == 0 || i == info->mx-1) continue;
172:     row     = i-info->gxs;
173:     cols[0] = i-1-infok->gxs;  vals[0] = (u[i]-u[i-1])/hx;
174:     cols[1] = i-infok->gxs;    vals[1] = (u[i]-u[i+1])/hx;
175:     MatSetValuesLocal(Buk,1,&row,2,cols,vals,INSERT_VALUES);
176:   }
177:   return 0;
178: }

180: static PetscErrorCode FormJacobianLocal_KU(User user,DMDALocalInfo *info,DMDALocalInfo *infok,const PetscScalar u[],const PetscScalar k[],Mat Bku)
181: {
182:   PetscInt       i;
183:   PetscReal      hx = 1./(info->mx-1);

186:   if (!Bku) return 0; /* Not assembling this block */
187:   for (i=infok->xs; i<infok->xs+infok->xm; i++) {
188:     PetscInt    row = i-infok->gxs,cols[2];
189:     PetscScalar vals[2];
190:     const PetscScalar
191:       ubar     = 0.5*(u[i]+u[i+1]),
192:       ubar_L   = 0.5,
193:       ubar_R   = 0.5,
194:       gradu    = (u[i+1]-u[i])/hx,
195:       gradu_L  = -1./hx,
196:       gradu_R  = 1./hx,
197:       g        = 1. + PetscSqr(gradu),
198:       g_gradu  = 2.*gradu,
199:       w        = 1./(1.+ubar) + 1./g,
200:       w_ubar   = -1./PetscSqr(1.+ubar),
201:       w_gradu  = -g_gradu/PetscSqr(g),
202:       iw       = 1./w,
203:       iw_ubar  = -w_ubar * PetscSqr(iw),
204:       iw_gradu = -w_gradu * PetscSqr(iw);
205:     cols[0] = i-info->gxs;         vals[0] = -hx*(iw_ubar*ubar_L + iw_gradu*gradu_L);
206:     cols[1] = i+1-info->gxs;       vals[1] = -hx*(iw_ubar*ubar_R + iw_gradu*gradu_R);
207:     MatSetValuesLocal(Bku,1,&row,2,cols,vals,INSERT_VALUES);
208:   }
209:   return 0;
210: }

212: static PetscErrorCode FormJacobian_All(SNES snes,Vec X,Mat J,Mat B,void *ctx)
213: {
214:   User           user = (User)ctx;
215:   DM             dau,dak;
216:   DMDALocalInfo  infou,infok;
217:   PetscScalar    *u,*k;
218:   Vec            Uloc,Kloc;

221:   DMCompositeGetEntries(user->pack,&dau,&dak);
222:   DMDAGetLocalInfo(dau,&infou);
223:   DMDAGetLocalInfo(dak,&infok);
224:   DMCompositeGetLocalVectors(user->pack,&Uloc,&Kloc);
225:   switch (user->ptype) {
226:   case 0:
227:     DMGlobalToLocalBegin(dau,X,INSERT_VALUES,Uloc);
228:     DMGlobalToLocalEnd  (dau,X,INSERT_VALUES,Uloc);
229:     DMDAVecGetArray(dau,Uloc,&u);
230:     DMDAVecGetArray(dak,user->Kloc,&k);
231:     FormJacobianLocal_U(user,&infou,u,k,B);
232:     DMDAVecRestoreArray(dau,Uloc,&u);
233:     DMDAVecRestoreArray(dak,user->Kloc,&k);
234:     break;
235:   case 1:
236:     DMGlobalToLocalBegin(dak,X,INSERT_VALUES,Kloc);
237:     DMGlobalToLocalEnd  (dak,X,INSERT_VALUES,Kloc);
238:     DMDAVecGetArray(dau,user->Uloc,&u);
239:     DMDAVecGetArray(dak,Kloc,&k);
240:     FormJacobianLocal_K(user,&infok,u,k,B);
241:     DMDAVecRestoreArray(dau,user->Uloc,&u);
242:     DMDAVecRestoreArray(dak,Kloc,&k);
243:     break;
244:   case 2: {
245:     Mat       Buu,Buk,Bku,Bkk;
246:     PetscBool nest;
247:     IS        *is;
248:     DMCompositeScatter(user->pack,X,Uloc,Kloc);
249:     DMDAVecGetArray(dau,Uloc,&u);
250:     DMDAVecGetArray(dak,Kloc,&k);
251:     DMCompositeGetLocalISs(user->pack,&is);
252:     MatGetLocalSubMatrix(B,is[0],is[0],&Buu);
253:     MatGetLocalSubMatrix(B,is[0],is[1],&Buk);
254:     MatGetLocalSubMatrix(B,is[1],is[0],&Bku);
255:     MatGetLocalSubMatrix(B,is[1],is[1],&Bkk);
256:     FormJacobianLocal_U(user,&infou,u,k,Buu);
257:     PetscObjectTypeCompare((PetscObject)B,MATNEST,&nest);
258:     if (!nest) {
259:       /*
260:          DMCreateMatrix_Composite()  for a nested matrix does not generate off-block matrices that one can call MatSetValuesLocal() on, it just creates dummy
261:          matrices with no entries; there cannot insert values into them. Commit b6480e041dd2293a65f96222772d68cdb4ed6306
262:          changed Mat_Nest() from returning NULL pointers for these submatrices to dummy matrices because PCFIELDSPLIT could not
263:          handle the returned null matrices.
264:       */
265:       FormJacobianLocal_UK(user,&infou,&infok,u,k,Buk);
266:       FormJacobianLocal_KU(user,&infou,&infok,u,k,Bku);
267:     }
268:     FormJacobianLocal_K(user,&infok,u,k,Bkk);
269:     MatRestoreLocalSubMatrix(B,is[0],is[0],&Buu);
270:     MatRestoreLocalSubMatrix(B,is[0],is[1],&Buk);
271:     MatRestoreLocalSubMatrix(B,is[1],is[0],&Bku);
272:     MatRestoreLocalSubMatrix(B,is[1],is[1],&Bkk);
273:     DMDAVecRestoreArray(dau,Uloc,&u);
274:     DMDAVecRestoreArray(dak,Kloc,&k);

276:     ISDestroy(&is[0]);
277:     ISDestroy(&is[1]);
278:     PetscFree(is);
279:   } break;
280:   }
281:   DMCompositeRestoreLocalVectors(user->pack,&Uloc,&Kloc);
282:   MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
283:   MatAssemblyEnd  (B,MAT_FINAL_ASSEMBLY);
284:   if (J != B) {
285:     MatAssemblyBegin(J,MAT_FINAL_ASSEMBLY);
286:     MatAssemblyEnd  (J,MAT_FINAL_ASSEMBLY);
287:   }
288:   return 0;
289: }

291: static PetscErrorCode FormInitial_Coupled(User user,Vec X)
292: {
293:   DM             dau,dak;
294:   DMDALocalInfo  infou,infok;
295:   Vec            Xu,Xk;
296:   PetscScalar    *u,*k,hx;
297:   PetscInt       i;

300:   DMCompositeGetEntries(user->pack,&dau,&dak);
301:   DMCompositeGetAccess(user->pack,X,&Xu,&Xk);
302:   DMDAVecGetArray(dau,Xu,&u);
303:   DMDAVecGetArray(dak,Xk,&k);
304:   DMDAGetLocalInfo(dau,&infou);
305:   DMDAGetLocalInfo(dak,&infok);
306:   hx   = 1./(infok.mx);
307:   for (i=infou.xs; i<infou.xs+infou.xm; i++) u[i] = (PetscScalar)i*hx * (1.-(PetscScalar)i*hx);
308:   for (i=infok.xs; i<infok.xs+infok.xm; i++) k[i] = 1.0 + 0.5*PetscSinScalar(2*PETSC_PI*i*hx);
309:   DMDAVecRestoreArray(dau,Xu,&u);
310:   DMDAVecRestoreArray(dak,Xk,&k);
311:   DMCompositeRestoreAccess(user->pack,X,&Xu,&Xk);
312:   DMCompositeScatter(user->pack,X,user->Uloc,user->Kloc);
313:   return 0;
314: }

316: int main(int argc, char *argv[])
317: {
319:   DM             dau,dak,pack;
320:   const PetscInt *lxu;
321:   PetscInt       *lxk,m,sizes;
322:   User           user;
323:   SNES           snes;
324:   Vec            X,F,Xu,Xk,Fu,Fk;
325:   Mat            B;
326:   IS             *isg;
327:   PetscBool      pass_dm;

329:   PetscInitialize(&argc,&argv,0,help);
330:   DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,10,1,1,NULL,&dau);
331:   DMSetOptionsPrefix(dau,"u_");
332:   DMSetFromOptions(dau);
333:   DMSetUp(dau);
334:   DMDAGetOwnershipRanges(dau,&lxu,0,0);
335:   DMDAGetInfo(dau,0, &m,0,0, &sizes,0,0, 0,0,0,0,0,0);
336:   PetscMalloc1(sizes,&lxk);
337:   PetscArraycpy(lxk,lxu,sizes);
338:   lxk[0]--;
339:   DMDACreate1d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,m-1,1,1,lxk,&dak);
340:   DMSetOptionsPrefix(dak,"k_");
341:   DMSetFromOptions(dak);
342:   DMSetUp(dak);
343:   PetscFree(lxk);

345:   DMCompositeCreate(PETSC_COMM_WORLD,&pack);
346:   DMSetOptionsPrefix(pack,"pack_");
347:   DMCompositeAddDM(pack,dau);
348:   DMCompositeAddDM(pack,dak);
349:   DMDASetFieldName(dau,0,"u");
350:   DMDASetFieldName(dak,0,"k");
351:   DMSetFromOptions(pack);

353:   DMCreateGlobalVector(pack,&X);
354:   VecDuplicate(X,&F);

356:   PetscNew(&user);

358:   user->pack = pack;

360:   DMCompositeGetGlobalISs(pack,&isg);
361:   DMCompositeGetLocalVectors(pack,&user->Uloc,&user->Kloc);
362:   DMCompositeScatter(pack,X,user->Uloc,user->Kloc);

364:   PetscOptionsBegin(PETSC_COMM_WORLD,NULL,"Coupled problem options","SNES");
365:   {
366:     user->ptype = 0; pass_dm = PETSC_TRUE;

368:     PetscOptionsInt("-problem_type","0: solve for u only, 1: solve for k only, 2: solve for both",0,user->ptype,&user->ptype,NULL);
369:     PetscOptionsBool("-pass_dm","Pass the packed DM to SNES to use when determining splits and forward into splits",0,pass_dm,&pass_dm,NULL);
370:   }
371:   PetscOptionsEnd();

373:   FormInitial_Coupled(user,X);

375:   SNESCreate(PETSC_COMM_WORLD,&snes);
376:   switch (user->ptype) {
377:   case 0:
378:     DMCompositeGetAccess(pack,X,&Xu,0);
379:     DMCompositeGetAccess(pack,F,&Fu,0);
380:     DMCreateMatrix(dau,&B);
381:     SNESSetFunction(snes,Fu,FormFunction_All,user);
382:     SNESSetJacobian(snes,B,B,FormJacobian_All,user);
383:     SNESSetFromOptions(snes);
384:     SNESSetDM(snes,dau);
385:     SNESSolve(snes,NULL,Xu);
386:     DMCompositeRestoreAccess(pack,X,&Xu,0);
387:     DMCompositeRestoreAccess(pack,F,&Fu,0);
388:     break;
389:   case 1:
390:     DMCompositeGetAccess(pack,X,0,&Xk);
391:     DMCompositeGetAccess(pack,F,0,&Fk);
392:     DMCreateMatrix(dak,&B);
393:     SNESSetFunction(snes,Fk,FormFunction_All,user);
394:     SNESSetJacobian(snes,B,B,FormJacobian_All,user);
395:     SNESSetFromOptions(snes);
396:     SNESSetDM(snes,dak);
397:     SNESSolve(snes,NULL,Xk);
398:     DMCompositeRestoreAccess(pack,X,0,&Xk);
399:     DMCompositeRestoreAccess(pack,F,0,&Fk);
400:     break;
401:   case 2:
402:     DMCreateMatrix(pack,&B);
403:     /* This example does not correctly allocate off-diagonal blocks. These options allows new nonzeros (slow). */
404:     MatSetOption(B,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE);
405:     MatSetOption(B,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
406:     SNESSetFunction(snes,F,FormFunction_All,user);
407:     SNESSetJacobian(snes,B,B,FormJacobian_All,user);
408:     SNESSetFromOptions(snes);
409:     if (!pass_dm) {             /* Manually provide index sets and names for the splits */
410:       KSP ksp;
411:       PC  pc;
412:       SNESGetKSP(snes,&ksp);
413:       KSPGetPC(ksp,&pc);
414:       PCFieldSplitSetIS(pc,"u",isg[0]);
415:       PCFieldSplitSetIS(pc,"k",isg[1]);
416:     } else {
417:       /* The same names come from the options prefix for dau and dak. This option can support geometric multigrid inside
418:        * of splits, but it requires using a DM (perhaps your own implementation). */
419:       SNESSetDM(snes,pack);
420:     }
421:     SNESSolve(snes,NULL,X);
422:     break;
423:   }
424:   VecViewFromOptions(X,NULL,"-view_sol");

426:   if (0) {
427:     PetscInt  col      = 0;
428:     PetscBool mult_dup = PETSC_FALSE,view_dup = PETSC_FALSE;
429:     Mat       D;
430:     Vec       Y;

432:     PetscOptionsGetInt(NULL,0,"-col",&col,0);
433:     PetscOptionsGetBool(NULL,0,"-mult_dup",&mult_dup,0);
434:     PetscOptionsGetBool(NULL,0,"-view_dup",&view_dup,0);

436:     VecDuplicate(X,&Y);
437:     /* MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY); */
438:     /* MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY); */
439:     MatConvert(B,MATAIJ,MAT_INITIAL_MATRIX,&D);
440:     VecZeroEntries(X);
441:     VecSetValue(X,col,1.0,INSERT_VALUES);
442:     VecAssemblyBegin(X);
443:     VecAssemblyEnd(X);
444:     MatMult(mult_dup ? D : B,X,Y);
445:     MatView(view_dup ? D : B,PETSC_VIEWER_STDOUT_WORLD);
446:     /* VecView(X,PETSC_VIEWER_STDOUT_WORLD); */
447:     VecView(Y,PETSC_VIEWER_STDOUT_WORLD);
448:     MatDestroy(&D);
449:     VecDestroy(&Y);
450:   }

452:   DMCompositeRestoreLocalVectors(pack,&user->Uloc,&user->Kloc);
453:   PetscFree(user);

455:   ISDestroy(&isg[0]);
456:   ISDestroy(&isg[1]);
457:   PetscFree(isg);
458:   VecDestroy(&X);
459:   VecDestroy(&F);
460:   MatDestroy(&B);
461:   DMDestroy(&dau);
462:   DMDestroy(&dak);
463:   DMDestroy(&pack);
464:   SNESDestroy(&snes);
465:   PetscFinalize();
466:   return 0;
467: }

469: /*TEST

471:    test:
472:       suffix: 0
473:       nsize: 3
474:       args: -u_da_grid_x 20 -snes_converged_reason -snes_monitor_short -problem_type 0

476:    test:
477:       suffix: 1
478:       nsize: 3
479:       args: -u_da_grid_x 20 -snes_converged_reason -snes_monitor_short -problem_type 1

481:    test:
482:       suffix: 2
483:       nsize: 3
484:       args: -u_da_grid_x 20 -snes_converged_reason -snes_monitor_short -problem_type 2

486:    test:
487:       suffix: 3
488:       nsize: 3
489:       args: -u_da_grid_x 20 -snes_converged_reason -snes_monitor_short -ksp_monitor_short -problem_type 2 -snes_mf_operator -pack_dm_mat_type {{aij nest}} -pc_type fieldsplit -pc_fieldsplit_dm_splits -pc_fieldsplit_type additive -fieldsplit_u_ksp_type gmres -fieldsplit_k_pc_type jacobi

491:    test:
492:       suffix: 4
493:       nsize: 6
494:       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
495:       requires: !single

497:    test:
498:       suffix: glvis_composite_da_1d
499:       args: -u_da_grid_x 20 -problem_type 0 -snes_monitor_solution glvis:

501:    test:
502:       suffix: cuda
503:       nsize: 1
504:       requires: cuda
505:       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

507:    test:
508:       suffix: viennacl
509:       nsize: 1
510:       requires: viennacl
511:       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

513: TEST*/