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*/