Actual source code: ex12.c

petsc-main 2021-04-20
Report Typos and Errors
  1: static char help[] = "Poisson Problem in 2d and 3d with simplicial finite elements.\n\
  2: We solve the Poisson problem in a rectangular\n\
  3: domain, using a parallel unstructured mesh (DMPLEX) to discretize it.\n\
  4: This example supports discretized auxiliary fields (conductivity) as well as\n\
  5: multilevel nonlinear solvers.\n\n\n";

  7: /*
  8: A visualization of the adaptation can be accomplished using:

 10:   -dm_adapt_view hdf5:$PWD/adapt.h5 -sol_adapt_view hdf5:$PWD/adapt.h5::append -dm_adapt_pre_view hdf5:$PWD/orig.h5 -sol_adapt_pre_view hdf5:$PWD/orig.h5::append

 12: Information on refinement:

 14:    -info :~sys,vec,is,mat,ksp,snes,ts
 15: */

 17: #include <petscdmplex.h>
 18: #include <petscdmadaptor.h>
 19: #include <petscsnes.h>
 20: #include <petscds.h>
 21: #include <petscviewerhdf5.h>

 23: typedef enum {NEUMANN, DIRICHLET, NONE} BCType;
 24: typedef enum {RUN_FULL, RUN_EXACT, RUN_TEST, RUN_PERF} RunType;
 25: typedef enum {COEFF_NONE, COEFF_ANALYTIC, COEFF_FIELD, COEFF_NONLINEAR, COEFF_CIRCLE, COEFF_CROSS, COEFF_CHECKERBOARD_0, COEFF_CHECKERBOARD_1} CoeffType;

 27: typedef struct {
 28:   PetscInt       debug;             /* The debugging level */
 29:   RunType        runType;           /* Whether to run tests, or solve the full problem */
 30:   PetscBool      jacobianMF;        /* Whether to calculate the Jacobian action on the fly */
 31:   PetscLogEvent  createMeshEvent;
 32:   PetscBool      showInitial, showSolution, restart, quiet, nonzInit;
 33:   /* Domain and mesh definition */
 34:   PetscInt       dim;               /* The topological mesh dimension */
 35:   DMBoundaryType periodicity[3];    /* The domain periodicity */
 36:   PetscInt       cells[3];          /* The initial domain division */
 37:   char           filename[2048];    /* The optional mesh file */
 38:   PetscBool      interpolate;       /* Generate intermediate mesh elements */
 39:   PetscReal      refinementLimit;   /* The largest allowable cell volume */
 40:   PetscBool      viewHierarchy;     /* Whether to view the hierarchy */
 41:   PetscBool      simplex;           /* Simplicial mesh */
 42:   /* Problem definition */
 43:   BCType         bcType;
 44:   CoeffType      variableCoefficient;
 45:   PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx);
 46:   PetscBool      fieldBC;
 47:   void           (**exactFields)(PetscInt, PetscInt, PetscInt,
 48:                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
 49:                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
 50:                                  PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]);
 51:   PetscBool      bdIntegral;        /* Compute the integral of the solution on the boundary */
 52:   /* Reproducing tests from SISC 40(3), pp. A1473-A1493, 2018 */
 53:   PetscInt       div;               /* Number of divisions */
 54:   PetscInt       k;                 /* Parameter for checkerboard coefficient */
 55:   PetscInt      *kgrid;             /* Random parameter grid */
 56:   /* Solver */
 57:   PC             pcmg;              /* This is needed for error monitoring */
 58:   PetscBool      checkksp;          /* Whether to check the KSPSolve for runType == RUN_TEST */
 59: } AppCtx;

 61: static PetscErrorCode zero(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 62: {
 63:   u[0] = 0.0;
 64:   return 0;
 65: }

 67: static PetscErrorCode ecks(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 68: {
 69:   u[0] = x[0];
 70:   return 0;
 71: }

 73: /*
 74:   In 2D for Dirichlet conditions, we use exact solution:

 76:     u = x^2 + y^2
 77:     f = 4

 79:   so that

 81:     -\Delta u + f = -4 + 4 = 0

 83:   For Neumann conditions, we have

 85:     -\nabla u \cdot -\hat y |_{y=0} =  (2y)|_{y=0} =  0 (bottom)
 86:     -\nabla u \cdot  \hat y |_{y=1} = -(2y)|_{y=1} = -2 (top)
 87:     -\nabla u \cdot -\hat x |_{x=0} =  (2x)|_{x=0} =  0 (left)
 88:     -\nabla u \cdot  \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)

 90:   Which we can express as

 92:     \nabla u \cdot  \hat n|_\Gamma = {2 x, 2 y} \cdot \hat n = 2 (x + y)

 94:   The boundary integral of this solution is (assuming we are not orienting the edges)

 96:     \int^1_0 x^2 dx + \int^1_0 (1 + y^2) dy + \int^1_0 (x^2 + 1) dx + \int^1_0 y^2 dy = 1/3 + 4/3 + 4/3 + 1/3 = 3 1/3
 97: */
 98: static PetscErrorCode quadratic_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 99: {
100:   *u = x[0]*x[0] + x[1]*x[1];
101:   return 0;
102: }

104: static void quadratic_u_field_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
105:                                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
106:                                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
107:                                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
108: {
109:   uexact[0] = a[0];
110: }

112: static PetscErrorCode circle_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
113: {
114:   const PetscReal alpha   = 500.;
115:   const PetscReal radius2 = PetscSqr(0.15);
116:   const PetscReal r2      = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5);
117:   const PetscReal xi      = alpha*(radius2 - r2);

119:   *u = PetscTanhScalar(xi) + 1.0;
120:   return 0;
121: }

123: static PetscErrorCode cross_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
124: {
125:   const PetscReal alpha = 50*4;
126:   const PetscReal xy    = (x[0]-0.5)*(x[1]-0.5);

128:   *u = PetscSinReal(alpha*xy) * (alpha*PetscAbsReal(xy) < 2*PETSC_PI ? 1 : 0.01);
129:   return 0;
130: }

132: static void f0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
133:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
134:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
135:                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
136: {
137:   f0[0] = 4.0;
138: }

140: static void f0_circle_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
141:                         const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
142:                         const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
143:                         PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
144: {
145:   const PetscReal alpha   = 500.;
146:   const PetscReal radius2 = PetscSqr(0.15);
147:   const PetscReal r2      = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5);
148:   const PetscReal xi      = alpha*(radius2 - r2);

150:   f0[0] = (-4.0*alpha - 8.0*PetscSqr(alpha)*r2*PetscTanhReal(xi)) * PetscSqr(1.0/PetscCoshReal(xi));
151: }

153: static void f0_cross_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
154:                        const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
155:                        const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
156:                        PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
157: {
158:   const PetscReal alpha = 50*4;
159:   const PetscReal xy    = (x[0]-0.5)*(x[1]-0.5);

161:   f0[0] = PetscSinReal(alpha*xy) * (alpha*PetscAbsReal(xy) < 2*PETSC_PI ? 1 : 0.01);
162: }

164: static void f0_checkerboard_0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
165:                                 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
166:                                 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
167:                                 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
168: {
169:   f0[0] = -20.0*PetscExpReal(-(PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5)));
170: }

172: static void f0_bd_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
173:                     const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
174:                     const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
175:                     PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
176: {
177:   PetscInt d;
178:   for (d = 0, f0[0] = 0.0; d < dim; ++d) f0[0] += -n[d]*2.0*x[d];
179: }

181: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
182: static void f1_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
183:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
184:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
185:                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
186: {
187:   PetscInt d;
188:   for (d = 0; d < dim; ++d) f1[d] = u_x[d];
189: }

191: /* < \nabla v, \nabla u + {\nabla u}^T >
192:    This just gives \nabla u, give the perdiagonal for the transpose */
193: static void g3_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
194:                   const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
195:                   const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
196:                   PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
197: {
198:   PetscInt d;
199:   for (d = 0; d < dim; ++d) g3[d*dim+d] = 1.0;
200: }

202: /*
203:   In 2D for x periodicity and y Dirichlet conditions, we use exact solution:

205:     u = sin(2 pi x)
206:     f = -4 pi^2 sin(2 pi x)

208:   so that

210:     -\Delta u + f = 4 pi^2 sin(2 pi x) - 4 pi^2 sin(2 pi x) = 0
211: */
212: static PetscErrorCode xtrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
213: {
214:   *u = PetscSinReal(2.0*PETSC_PI*x[0]);
215:   return 0;
216: }

218: static void f0_xtrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
219:                        const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
220:                        const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
221:                        PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
222: {
223:   f0[0] = -4.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
224: }

226: /*
227:   In 2D for x-y periodicity, we use exact solution:

229:     u = sin(2 pi x) sin(2 pi y)
230:     f = -8 pi^2 sin(2 pi x)

232:   so that

234:     -\Delta u + f = 4 pi^2 sin(2 pi x) sin(2 pi y) + 4 pi^2 sin(2 pi x) sin(2 pi y) - 8 pi^2 sin(2 pi x) = 0
235: */
236: static PetscErrorCode xytrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
237: {
238:   *u = PetscSinReal(2.0*PETSC_PI*x[0])*PetscSinReal(2.0*PETSC_PI*x[1]);
239:   return 0;
240: }

242: static void f0_xytrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
243:                         const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
244:                         const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
245:                         PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
246: {
247:   f0[0] = -8.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
248: }

250: /*
251:   In 2D for Dirichlet conditions with a variable coefficient, we use exact solution:

253:     u  = x^2 + y^2
254:     f  = 6 (x + y)
255:     nu = (x + y)

257:   so that

259:     -\div \nu \grad u + f = -6 (x + y) + 6 (x + y) = 0
260: */
261: static PetscErrorCode nu_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
262: {
263:   *u = x[0] + x[1];
264:   return 0;
265: }

267: static PetscErrorCode checkerboardCoeff(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
268: {
269:   AppCtx  *user = (AppCtx *) ctx;
270:   PetscInt div  = user->div;
271:   PetscInt k    = user->k;
272:   PetscInt mask = 0, ind = 0, d;

275:   for (d = 0; d < dim; ++d) mask = (mask + (PetscInt) (x[d]*div)) % 2;
276:   if (user->kgrid) {
277:     for (d = 0; d < dim; ++d) {
278:       if (d > 0) ind *= dim;
279:       ind += (PetscInt) (x[d]*div);
280:     }
281:     k = user->kgrid[ind];
282:   }
283:   u[0] = mask ? 1.0 : PetscPowRealInt(10.0, -k);
284:   return(0);
285: }

287: void f0_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
288:                    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
289:                    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
290:                    PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
291: {
292:   f0[0] = 6.0*(x[0] + x[1]);
293: }

295: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
296: void f1_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
297:                    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
298:                    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
299:                    PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
300: {
301:   PetscInt d;
302:   for (d = 0; d < dim; ++d) f1[d] = (x[0] + x[1])*u_x[d];
303: }

305: void f1_field_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
306:                 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
307:                 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
308:                 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
309: {
310:   PetscInt d;
311:   for (d = 0; d < dim; ++d) f1[d] = a[0]*u_x[d];
312: }

314: /* < \nabla v, \nabla u + {\nabla u}^T >
315:    This just gives \nabla u, give the perdiagonal for the transpose */
316: void g3_analytic_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
317:                     const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
318:                     const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
319:                     PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
320: {
321:   PetscInt d;
322:   for (d = 0; d < dim; ++d) g3[d*dim+d] = x[0] + x[1];
323: }

325: void g3_field_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
326:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
327:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
328:                  PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
329: {
330:   PetscInt d;
331:   for (d = 0; d < dim; ++d) g3[d*dim+d] = a[0];
332: }

334: /*
335:   In 2D for Dirichlet conditions with a nonlinear coefficient (p-Laplacian with p = 4), we use exact solution:

337:     u  = x^2 + y^2
338:     f  = 16 (x^2 + y^2)
339:     nu = 1/2 |grad u|^2

341:   so that

343:     -\div \nu \grad u + f = -16 (x^2 + y^2) + 16 (x^2 + y^2) = 0
344: */
345: void f0_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
346:                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
347:                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
348:                              PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
349: {
350:   f0[0] = 16.0*(x[0]*x[0] + x[1]*x[1]);
351: }

353: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
354: void f1_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
355:                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
356:                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
357:                              PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
358: {
359:   PetscScalar nu = 0.0;
360:   PetscInt    d;
361:   for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
362:   for (d = 0; d < dim; ++d) f1[d] = 0.5*nu*u_x[d];
363: }

365: /*
366:   grad (u + eps w) - grad u = eps grad w

368:   1/2 |grad (u + eps w)|^2 grad (u + eps w) - 1/2 |grad u|^2 grad u
369: = 1/2 (|grad u|^2 + 2 eps <grad u,grad w>) (grad u + eps grad w) - 1/2 |grad u|^2 grad u
370: = 1/2 (eps |grad u|^2 grad w + 2 eps <grad u,grad w> grad u)
371: = eps (1/2 |grad u|^2 grad w + grad u <grad u,grad w>)
372: */
373: void g3_analytic_nonlinear_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
374:                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
375:                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
376:                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
377: {
378:   PetscScalar nu = 0.0;
379:   PetscInt    d, e;
380:   for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
381:   for (d = 0; d < dim; ++d) {
382:     g3[d*dim+d] = 0.5*nu;
383:     for (e = 0; e < dim; ++e) {
384:       g3[d*dim+e] += u_x[d]*u_x[e];
385:     }
386:   }
387: }

389: /*
390:   In 3D for Dirichlet conditions we use exact solution:

392:     u = 2/3 (x^2 + y^2 + z^2)
393:     f = 4

395:   so that

397:     -\Delta u + f = -2/3 * 6 + 4 = 0

399:   For Neumann conditions, we have

401:     -\nabla u \cdot -\hat z |_{z=0} =  (2z)|_{z=0} =  0 (bottom)
402:     -\nabla u \cdot  \hat z |_{z=1} = -(2z)|_{z=1} = -2 (top)
403:     -\nabla u \cdot -\hat y |_{y=0} =  (2y)|_{y=0} =  0 (front)
404:     -\nabla u \cdot  \hat y |_{y=1} = -(2y)|_{y=1} = -2 (back)
405:     -\nabla u \cdot -\hat x |_{x=0} =  (2x)|_{x=0} =  0 (left)
406:     -\nabla u \cdot  \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)

408:   Which we can express as

410:     \nabla u \cdot  \hat n|_\Gamma = {2 x, 2 y, 2z} \cdot \hat n = 2 (x + y + z)
411: */
412: static PetscErrorCode quadratic_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
413: {
414:   *u = 2.0*(x[0]*x[0] + x[1]*x[1] + x[2]*x[2])/3.0;
415:   return 0;
416: }

418: static void quadratic_u_field_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
419:                                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
420:                                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
421:                                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
422: {
423:   uexact[0] = a[0];
424: }

426: static void bd_integral_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
427:                            const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
428:                            const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
429:                            PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar *uint)
430: {
431:   uint[0] = u[0];
432: }

434: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
435: {
436:   const char    *bcTypes[3]  = {"neumann", "dirichlet", "none"};
437:   const char    *runTypes[4] = {"full", "exact", "test", "perf"};
438:   const char    *coeffTypes[8] = {"none", "analytic", "field", "nonlinear", "circle", "cross", "checkerboard_0", "checkerboard_1"};
439:   PetscInt       bd, bc, run, coeff, n;
440:   PetscBool      rand = PETSC_FALSE, flg;

444:   options->debug               = 0;
445:   options->runType             = RUN_FULL;
446:   options->dim                 = 2;
447:   options->periodicity[0]      = DM_BOUNDARY_NONE;
448:   options->periodicity[1]      = DM_BOUNDARY_NONE;
449:   options->periodicity[2]      = DM_BOUNDARY_NONE;
450:   options->cells[0]            = 2;
451:   options->cells[1]            = 2;
452:   options->cells[2]            = 2;
453:   options->filename[0]         = '\0';
454:   options->interpolate         = PETSC_TRUE;
455:   options->refinementLimit     = 0.0;
456:   options->bcType              = DIRICHLET;
457:   options->variableCoefficient = COEFF_NONE;
458:   options->fieldBC             = PETSC_FALSE;
459:   options->jacobianMF          = PETSC_FALSE;
460:   options->showInitial         = PETSC_FALSE;
461:   options->showSolution        = PETSC_FALSE;
462:   options->restart             = PETSC_FALSE;
463:   options->viewHierarchy       = PETSC_FALSE;
464:   options->simplex             = PETSC_TRUE;
465:   options->quiet               = PETSC_FALSE;
466:   options->nonzInit            = PETSC_FALSE;
467:   options->bdIntegral          = PETSC_FALSE;
468:   options->checkksp            = PETSC_FALSE;
469:   options->div                 = 4;
470:   options->k                   = 1;
471:   options->kgrid               = NULL;

473:   PetscOptionsBegin(comm, "", "Poisson Problem Options", "DMPLEX");
474:   PetscOptionsInt("-debug", "The debugging level", "ex12.c", options->debug, &options->debug, NULL);
475:   run  = options->runType;
476:   PetscOptionsEList("-run_type", "The run type", "ex12.c", runTypes, 4, runTypes[options->runType], &run, NULL);

478:   options->runType = (RunType) run;

480:   PetscOptionsInt("-dim", "The topological mesh dimension", "ex12.c", options->dim, &options->dim, NULL);
481:   bd = options->periodicity[0];
482:   PetscOptionsEList("-x_periodicity", "The x-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[0]], &bd, NULL);
483:   options->periodicity[0] = (DMBoundaryType) bd;
484:   bd = options->periodicity[1];
485:   PetscOptionsEList("-y_periodicity", "The y-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[1]], &bd, NULL);
486:   options->periodicity[1] = (DMBoundaryType) bd;
487:   bd = options->periodicity[2];
488:   PetscOptionsEList("-z_periodicity", "The z-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[2]], &bd, NULL);
489:   options->periodicity[2] = (DMBoundaryType) bd;
490:   n = 3;
491:   PetscOptionsIntArray("-cells", "The initial mesh division", "ex12.c", options->cells, &n, NULL);
492:   PetscOptionsString("-f", "Mesh filename to read", "ex12.c", options->filename, options->filename, sizeof(options->filename), &flg);
493:   PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex12.c", options->interpolate, &options->interpolate, NULL);
494:   PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex12.c", options->refinementLimit, &options->refinementLimit, NULL);
495:   bc   = options->bcType;
496:   PetscOptionsEList("-bc_type","Type of boundary condition","ex12.c",bcTypes,3,bcTypes[options->bcType],&bc,NULL);
497:   options->bcType = (BCType) bc;
498:   coeff = options->variableCoefficient;
499:   PetscOptionsEList("-variable_coefficient","Type of variable coefficent","ex12.c",coeffTypes,8,coeffTypes[options->variableCoefficient],&coeff,NULL);
500:   options->variableCoefficient = (CoeffType) coeff;

502:   PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL);
503:   PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL);
504:   PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL);
505:   PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL);
506:   PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL);
507:   PetscOptionsBool("-dm_view_hierarchy", "View the coarsened hierarchy", "ex12.c", options->viewHierarchy, &options->viewHierarchy, NULL);
508:   PetscOptionsBool("-simplex", "Simplicial (true) or tensor (false) mesh", "ex12.c", options->simplex, &options->simplex, NULL);
509:   PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL);
510:   PetscOptionsBool("-nonzero_initial_guess", "nonzero initial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL);
511:   PetscOptionsBool("-bd_integral", "Compute the integral of the solution on the boundary", "ex12.c", options->bdIntegral, &options->bdIntegral, NULL);
512:   if (options->runType == RUN_TEST) {
513:     PetscOptionsBool("-run_test_check_ksp", "Check solution of KSP", "ex12.c", options->checkksp, &options->checkksp, NULL);
514:   }
515:   PetscOptionsInt("-div", "The number of division for the checkerboard coefficient", "ex12.c", options->div, &options->div, NULL);
516:   PetscOptionsInt("-k", "The exponent for the checkerboard coefficient", "ex12.c", options->k, &options->k, NULL);
517:   PetscOptionsBool("-k_random", "Assign random k values to checkerboard", "ex12.c", rand, &rand, NULL);
518:   PetscOptionsEnd();
519:   PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);

521:   if (rand) {
522:     PetscRandom r;
523:     PetscReal   val;
524:     PetscInt    N = PetscPowInt(options->div, options->dim), i;

526:     PetscMalloc1(N, &options->kgrid);
527:     PetscRandomCreate(PETSC_COMM_SELF, &r);
528:     PetscRandomSetFromOptions(r);
529:     PetscRandomSetInterval(r, 0.0, options->k);
530:     PetscRandomSetSeed(r, 1973);
531:     PetscRandomSeed(r);
532:     for (i = 0; i < N; ++i) {
533:       PetscRandomGetValueReal(r, &val);
534:       options->kgrid[i] = 1 + (PetscInt) val;
535:     }
536:     PetscRandomDestroy(&r);
537:   }
538:   return(0);
539: }

541: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
542: {
543:   DM             plex;
544:   DMLabel        label;

548:   DMCreateLabel(dm, name);
549:   DMGetLabel(dm, name, &label);
550:   DMConvert(dm, DMPLEX, &plex);
551:   DMPlexMarkBoundaryFaces(plex, 1, label);
552:   DMDestroy(&plex);
553:   return(0);
554: }

556: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
557: {
558:   PetscInt       dim             = user->dim;
559:   const char    *filename        = user->filename;
560:   PetscBool      interpolate     = user->interpolate;
561:   PetscReal      refinementLimit = user->refinementLimit;
562:   size_t         len;

566:   PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
567:   PetscStrlen(filename, &len);
568:   if (!len) {
569:     PetscInt d;

571:     if (user->periodicity[0] || user->periodicity[1] || user->periodicity[2]) for (d = 0; d < dim; ++d) user->cells[d] = PetscMax(user->cells[d], 3);
572:     DMPlexCreateBoxMesh(comm, dim, user->simplex, user->cells, NULL, NULL, user->periodicity, interpolate, dm);
573:     PetscObjectSetName((PetscObject) *dm, "Mesh");
574:   } else {
575:     DMPlexCreateFromFile(comm, filename, interpolate, dm);
576:     DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
577:   }
578:   {
579:     PetscPartitioner part;
580:     DM               refinedMesh     = NULL;
581:     DM               distributedMesh = NULL;

583:     /* Refine mesh using a volume constraint */
584:     if (refinementLimit > 0.0) {
585:       DMPlexSetRefinementLimit(*dm, refinementLimit);
586:       DMRefine(*dm, comm, &refinedMesh);
587:       if (refinedMesh) {
588:         const char *name;

590:         PetscObjectGetName((PetscObject) *dm,         &name);
591:         PetscObjectSetName((PetscObject) refinedMesh,  name);
592:         DMDestroy(dm);
593:         *dm  = refinedMesh;
594:       }
595:     }
596:     /* Distribute mesh over processes */
597:     DMPlexGetPartitioner(*dm,&part);
598:     PetscPartitionerSetFromOptions(part);
599:     DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
600:     if (distributedMesh) {
601:       DMDestroy(dm);
602:       *dm  = distributedMesh;
603:     }
604:   }
605:   if (interpolate) {
606:     if (user->bcType == NEUMANN) {
607:       DMLabel   label;

609:       DMCreateLabel(*dm, "boundary");
610:       DMGetLabel(*dm, "boundary", &label);
611:       DMPlexMarkBoundaryFaces(*dm, 1, label);
612:     } else if (user->bcType == DIRICHLET) {
613:       PetscBool hasLabel;

615:       DMHasLabel(*dm,"marker",&hasLabel);
616:       if (!hasLabel) {
617:         //DMSetUp(*dm);
618:         CreateBCLabel(*dm, "marker");
619:       }
620:     }
621:   }
622:   {
623:     char      convType[256];
624:     PetscBool flg;

626:     PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
627:     PetscOptionsFList("-dm_plex_convert_type","Convert DMPlex to another format","ex12",DMList,DMPLEX,convType,256,&flg);
628:     PetscOptionsEnd();
629:     if (flg) {
630:       DM dmConv;

632:       DMConvert(*dm,convType,&dmConv);
633:       if (dmConv) {
634:         DMDestroy(dm);
635:         *dm  = dmConv;
636:       }
637:     }
638:   }
639:   DMLocalizeCoordinates(*dm); /* needed for periodic */
640:   DMSetFromOptions(*dm);
641:   DMViewFromOptions(*dm, NULL, "-dm_view");
642:   if (user->viewHierarchy) {
643:     DM       cdm = *dm;
644:     PetscInt i   = 0;
645:     char     buf[256];

647:     while (cdm) {
648:       DMSetUp(cdm);
649:       DMGetCoarseDM(cdm, &cdm);
650:       ++i;
651:     }
652:     cdm = *dm;
653:     while (cdm) {
654:       PetscViewer       viewer;
655:       PetscBool   isHDF5, isVTK;

657:       --i;
658:       PetscViewerCreate(comm,&viewer);
659:       PetscViewerSetType(viewer,PETSCVIEWERHDF5);
660:       PetscViewerSetOptionsPrefix(viewer,"hierarchy_");
661:       PetscViewerSetFromOptions(viewer);
662:       PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&isHDF5);
663:       PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERVTK,&isVTK);
664:       if (isHDF5) {
665:         PetscSNPrintf(buf, 256, "ex12-%d.h5", i);
666:       } else if (isVTK) {
667:         PetscSNPrintf(buf, 256, "ex12-%d.vtu", i);
668:         PetscViewerPushFormat(viewer,PETSC_VIEWER_VTK_VTU);
669:       } else {
670:         PetscSNPrintf(buf, 256, "ex12-%d", i);
671:       }
672:       PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
673:       PetscViewerFileSetName(viewer,buf);
674:       DMView(cdm, viewer);
675:       PetscViewerDestroy(&viewer);
676:       DMGetCoarseDM(cdm, &cdm);
677:     }
678:   }
679:   PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
680:   return(0);
681: }

683: static PetscErrorCode SetupProblem(DM dm, AppCtx *user)
684: {
685:   PetscDS        ds;
686:   DMLabel        label;
687:   PetscWeakForm  wf;
688:   const PetscInt id = 1;
689:   PetscInt       bd;

693:   DMGetDS(dm, &ds);
694:   switch (user->variableCoefficient) {
695:   case COEFF_NONE:
696:     if (user->periodicity[0]) {
697:       if (user->periodicity[1]) {
698:         PetscDSSetResidual(ds, 0, f0_xytrig_u, f1_u);
699:         PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
700:       } else {
701:         PetscDSSetResidual(ds, 0, f0_xtrig_u,  f1_u);
702:         PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
703:       }
704:     } else {
705:       PetscDSSetResidual(ds, 0, f0_u, f1_u);
706:       PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
707:     }
708:     break;
709:   case COEFF_ANALYTIC:
710:     PetscDSSetResidual(ds, 0, f0_analytic_u, f1_analytic_u);
711:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_analytic_uu);
712:     break;
713:   case COEFF_FIELD:
714:     PetscDSSetResidual(ds, 0, f0_analytic_u, f1_field_u);
715:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_field_uu);
716:     break;
717:   case COEFF_NONLINEAR:
718:     PetscDSSetResidual(ds, 0, f0_analytic_nonlinear_u, f1_analytic_nonlinear_u);
719:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_analytic_nonlinear_uu);
720:     break;
721:   case COEFF_CIRCLE:
722:     PetscDSSetResidual(ds, 0, f0_circle_u, f1_u);
723:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
724:     break;
725:   case COEFF_CROSS:
726:     PetscDSSetResidual(ds, 0, f0_cross_u, f1_u);
727:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
728:     break;
729:   case COEFF_CHECKERBOARD_0:
730:     PetscDSSetResidual(ds, 0, f0_checkerboard_0_u, f1_field_u);
731:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_field_uu);
732:     break;
733:   default: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid variable coefficient type %d", user->variableCoefficient);
734:   }
735:   switch (user->dim) {
736:   case 2:
737:     switch (user->variableCoefficient) {
738:     case COEFF_CIRCLE:
739:       user->exactFuncs[0]  = circle_u_2d;break;
740:     case COEFF_CROSS:
741:       user->exactFuncs[0]  = cross_u_2d;break;
742:     case COEFF_CHECKERBOARD_0:
743:       user->exactFuncs[0]  = zero;break;
744:     default:
745:       if (user->periodicity[0]) {
746:         if (user->periodicity[1]) {
747:           user->exactFuncs[0] = xytrig_u_2d;
748:         } else {
749:           user->exactFuncs[0] = xtrig_u_2d;
750:         }
751:       } else {
752:         user->exactFuncs[0]  = quadratic_u_2d;
753:         user->exactFields[0] = quadratic_u_field_2d;
754:       }
755:     }
756:     if (user->bcType == NEUMANN) {
757:       DMGetLabel(dm, "boundary", &label);
758:       DMAddBoundary(dm, DM_BC_NATURAL, "wall", label, 1, &id, 0, 0, NULL, NULL, NULL, user, &bd);
759:       PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL);
760:       PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, f0_bd_u, 0, NULL);
761:     }
762:     break;
763:   case 3:
764:     user->exactFuncs[0]  = quadratic_u_3d;
765:     user->exactFields[0] = quadratic_u_field_3d;
766:     if (user->bcType == NEUMANN) {
767:       DMGetLabel(dm, "boundary", &label);
768:       DMAddBoundary(dm, DM_BC_NATURAL, "wall", label, 1, &id, 0, 0, NULL, NULL, NULL, user, &bd);
769:       PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL);
770:       PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, f0_bd_u, 0, NULL);
771:     }
772:     break;
773:   default:
774:     SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %d", user->dim);
775:   }
776:   /* Setup constants */
777:   switch (user->variableCoefficient) {
778:   case COEFF_CHECKERBOARD_0:
779:   {
780:     PetscScalar constants[2];

782:     constants[0] = user->div;
783:     constants[1] = user->k;
784:     PetscDSSetConstants(ds, 2, constants);
785:   }
786:   break;
787:   default: break;
788:   }
789:   PetscDSSetExactSolution(ds, 0, user->exactFuncs[0], user);
790:   /* Setup Boundary Conditions */
791:   if (user->bcType == DIRICHLET) {
792:     DMGetLabel(dm, "marker", &label);
793:     if (!label) {
794:       /* Right now, p4est cannot create labels immediately */
795:       PetscDSAddBoundaryByName(ds, user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL, "wall", "marker", 1, &id, 0, 0, NULL, user->fieldBC ? (void (*)(void)) user->exactFields[0] : (void (*)(void)) user->exactFuncs[0], NULL, user, NULL);
796:     } else {
797:       DMAddBoundary(dm, user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL, "wall", label, 1, &id, 0, 0, NULL, user->fieldBC ? (void (*)(void)) user->exactFields[0] : (void (*)(void)) user->exactFuncs[0], NULL, user, NULL);
798:     }
799:   }
800:   return(0);
801: }

803: static PetscErrorCode SetupMaterial(DM dm, DM dmAux, AppCtx *user)
804: {
805:   PetscErrorCode (*matFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {nu_2d};
806:   void            *ctx[1];
807:   Vec              nu;
808:   PetscErrorCode   ierr;

811:   ctx[0] = user;
812:   if (user->variableCoefficient == COEFF_CHECKERBOARD_0) {matFuncs[0] = checkerboardCoeff;}
813:   DMCreateLocalVector(dmAux, &nu);
814:   PetscObjectSetName((PetscObject) nu, "Coefficient");
815:   DMProjectFunctionLocal(dmAux, 0.0, matFuncs, ctx, INSERT_ALL_VALUES, nu);
816:   DMSetAuxiliaryVec(dm, NULL, 0, nu);
817:   VecDestroy(&nu);
818:   return(0);
819: }

821: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
822: {
823:   PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
824:   Vec            uexact;
825:   PetscInt       dim;

829:   DMGetDimension(dm, &dim);
830:   if (dim == 2) bcFuncs[0] = quadratic_u_2d;
831:   else          bcFuncs[0] = quadratic_u_3d;
832:   DMCreateLocalVector(dmAux, &uexact);
833:   DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
834:   DMSetAuxiliaryVec(dm, NULL, 0, uexact);
835:   VecDestroy(&uexact);
836:   return(0);
837: }

839: static PetscErrorCode SetupAuxDM(DM dm, PetscFE feAux, AppCtx *user)
840: {
841:   DM             dmAux, coordDM;

845:   /* MUST call DMGetCoordinateDM() in order to get p4est setup if present */
846:   DMGetCoordinateDM(dm, &coordDM);
847:   if (!feAux) return(0);
848:   DMClone(dm, &dmAux);
849:   DMSetCoordinateDM(dmAux, coordDM);
850:   DMSetField(dmAux, 0, NULL, (PetscObject) feAux);
851:   DMCreateDS(dmAux);
852:   if (user->fieldBC) {SetupBC(dm, dmAux, user);}
853:   else               {SetupMaterial(dm, dmAux, user);}
854:   DMDestroy(&dmAux);
855:   return(0);
856: }

858: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
859: {
860:   DM             cdm = dm;
861:   const PetscInt dim = user->dim;
862:   PetscFE        fe, feAux = NULL;
863:   PetscBool      simplex   = user->simplex;
864:   MPI_Comm       comm;

868:   /* Create finite element for each field and auxiliary field */
869:   PetscObjectGetComm((PetscObject) dm, &comm);
870:   PetscFECreateDefault(comm, dim, 1, simplex, NULL, -1, &fe);
871:   PetscObjectSetName((PetscObject) fe, "potential");
872:   if (user->variableCoefficient == COEFF_FIELD || user->variableCoefficient == COEFF_CHECKERBOARD_0) {
873:     PetscFECreateDefault(comm, dim, 1, simplex, "mat_", -1, &feAux);
874:     PetscObjectSetName((PetscObject) feAux, "coefficient");
875:     PetscFECopyQuadrature(fe, feAux);
876:   } else if (user->fieldBC) {
877:     PetscFECreateDefault(comm, dim, 1, simplex, "bc_", -1, &feAux);
878:     PetscFECopyQuadrature(fe, feAux);
879:   }
880:   /* Set discretization and boundary conditions for each mesh */
881:   DMSetField(dm, 0, NULL, (PetscObject) fe);
882:   DMCreateDS(dm);
883:   SetupProblem(dm, user);
884:   while (cdm) {
885:     SetupAuxDM(cdm, feAux, user);
886:     if (user->bcType == DIRICHLET && user->interpolate) {
887:       PetscBool hasLabel;

889:       DMHasLabel(cdm, "marker", &hasLabel);
890:       if (!hasLabel) {CreateBCLabel(cdm, "marker");}
891:     }
892:     DMCopyDisc(dm, cdm);
893:     DMGetCoarseDM(cdm, &cdm);
894:   }
895:   PetscFEDestroy(&fe);
896:   PetscFEDestroy(&feAux);
897:   return(0);
898: }

900: #include "petsc/private/petscimpl.h"

902: /*
903:   MonitorError - Outputs the error at each iteration of an iterative solver.

905:   Collective on KSP

907:   Input Parameters:
908: + ksp   - the KSP
909: . its   - iteration number
910: . rnorm - 2-norm, preconditioned residual value (may be estimated).
911: - ctx   - monitor context

913:   Level: intermediate

915: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual(), KSPMonitorResidual()
916: */
917: static PetscErrorCode MonitorError(KSP ksp, PetscInt its, PetscReal rnorm, void *ctx)
918: {
919:   AppCtx        *user = (AppCtx *) ctx;
920:   DM             dm;
921:   Vec            du = NULL, r;
922:   PetscInt       level = 0;
923:   PetscBool      hasLevel;
924: #if defined(PETSC_HAVE_HDF5)
925:   PetscViewer    viewer;
926:   char           buf[256];
927: #endif

931:   KSPGetDM(ksp, &dm);
932:   /* Calculate solution */
933:   {
934:     PC        pc = user->pcmg; /* The MG PC */
935:     DM        fdm = NULL,  cdm = NULL;
936:     KSP       fksp, cksp;
937:     Vec       fu,   cu = NULL;
938:     PetscInt  levels, l;

940:     KSPBuildSolution(ksp, NULL, &du);
941:     PetscObjectComposedDataGetInt((PetscObject) ksp, PetscMGLevelId, level, hasLevel);
942:     PCMGGetLevels(pc, &levels);
943:     PCMGGetSmoother(pc, levels-1, &fksp);
944:     KSPBuildSolution(fksp, NULL, &fu);
945:     for (l = levels-1; l > level; --l) {
946:       Mat R;
947:       Vec s;

949:       PCMGGetSmoother(pc, l-1, &cksp);
950:       KSPGetDM(cksp, &cdm);
951:       DMGetGlobalVector(cdm, &cu);
952:       PCMGGetRestriction(pc, l, &R);
953:       PCMGGetRScale(pc, l, &s);
954:       MatRestrict(R, fu, cu);
955:       VecPointwiseMult(cu, cu, s);
956:       if (l < levels-1) {DMRestoreGlobalVector(fdm, &fu);}
957:       fdm  = cdm;
958:       fu   = cu;
959:     }
960:     if (levels-1 > level) {
961:       VecAXPY(du, 1.0, cu);
962:       DMRestoreGlobalVector(cdm, &cu);
963:     }
964:   }
965:   /* Calculate error */
966:   DMGetGlobalVector(dm, &r);
967:   DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
968:   VecAXPY(r,-1.0,du);
969:   PetscObjectSetName((PetscObject) r, "solution error");
970:   /* View error */
971: #if defined(PETSC_HAVE_HDF5)
972:   PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
973:   PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
974:   VecView(r, viewer);
975:   PetscViewerDestroy(&viewer);
976: #endif
977:   DMRestoreGlobalVector(dm, &r);
978:   return(0);
979: }

981: /*@C
982:   SNESMonitorError - Outputs the error at each iteration of an iterative solver.

984:   Collective on SNES

986:   Input Parameters:
987: + snes  - the SNES
988: . its   - iteration number
989: . rnorm - 2-norm of residual
990: - ctx   - user context

992:   Level: intermediate

994: .seealso: SNESMonitorDefault(), SNESMonitorSet(), SNESMonitorSolution()
995: @*/
996: static PetscErrorCode SNESMonitorError(SNES snes, PetscInt its, PetscReal rnorm, void *ctx)
997: {
998:   AppCtx        *user = (AppCtx *) ctx;
999:   DM             dm;
1000:   Vec            u, r;
1001:   PetscInt       level = -1;
1002:   PetscBool      hasLevel;
1003: #if defined(PETSC_HAVE_HDF5)
1004:   PetscViewer    viewer;
1005: #endif
1006:   char           buf[256];

1010:   SNESGetDM(snes, &dm);
1011:   /* Calculate error */
1012:   SNESGetSolution(snes, &u);
1013:   DMGetGlobalVector(dm, &r);
1014:   PetscObjectSetName((PetscObject) r, "solution error");
1015:   DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
1016:   VecAXPY(r, -1.0, u);
1017:   /* View error */
1018:   PetscObjectComposedDataGetInt((PetscObject) snes, PetscMGLevelId, level, hasLevel);
1019:   PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
1020: #if defined(PETSC_HAVE_HDF5)
1021:   PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
1022:   VecView(r, viewer);
1023:   PetscViewerDestroy(&viewer);
1024:   /* Cleanup */
1025:   DMRestoreGlobalVector(dm, &r);
1026:   return(0);
1027: #else
1028:   SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"You need to configure with --download-hdf5");
1029: #endif
1030: }

1032: int main(int argc, char **argv)
1033: {
1034:   DM             dm;          /* Problem specification */
1035:   SNES           snes;        /* nonlinear solver */
1036:   Vec            u;           /* solution vector */
1037:   Mat            A,J;         /* Jacobian matrix */
1038:   MatNullSpace   nullSpace;   /* May be necessary for Neumann conditions */
1039:   AppCtx         user;        /* user-defined work context */
1040:   JacActionCtx   userJ;       /* context for Jacobian MF action */
1041:   PetscReal      error = 0.0; /* L_2 error in the solution */
1042:   PetscBool      isFAS;

1045:   PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
1046:   ProcessOptions(PETSC_COMM_WORLD, &user);
1047:   SNESCreate(PETSC_COMM_WORLD, &snes);
1048:   CreateMesh(PETSC_COMM_WORLD, &user, &dm);
1049:   SNESSetDM(snes, dm);
1050:   DMSetApplicationContext(dm, &user);

1052:   PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
1053:   SetupDiscretization(dm, &user);

1055:   DMCreateGlobalVector(dm, &u);
1056:   PetscObjectSetName((PetscObject) u, "potential");

1058:   DMCreateMatrix(dm, &J);
1059:   if (user.jacobianMF) {
1060:     PetscInt M, m, N, n;

1062:     MatGetSize(J, &M, &N);
1063:     MatGetLocalSize(J, &m, &n);
1064:     MatCreate(PETSC_COMM_WORLD, &A);
1065:     MatSetSizes(A, m, n, M, N);
1066:     MatSetType(A, MATSHELL);
1067:     MatSetUp(A);
1068: #if 0
1069:     MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction);
1070: #endif

1072:     userJ.dm   = dm;
1073:     userJ.J    = J;
1074:     userJ.user = &user;

1076:     DMCreateLocalVector(dm, &userJ.u);
1077:     if (user.fieldBC) {DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);}
1078:     else              {DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);}
1079:     MatShellSetContext(A, &userJ);
1080:   } else {
1081:     A = J;
1082:   }

1084:   nullSpace = NULL;
1085:   if (user.bcType != DIRICHLET) {
1086:     MatNullSpaceCreate(PetscObjectComm((PetscObject) dm), PETSC_TRUE, 0, NULL, &nullSpace);
1087:     MatSetNullSpace(A, nullSpace);
1088:   }

1090:   DMPlexSetSNESLocalFEM(dm,&user,&user,&user);
1091:   SNESSetJacobian(snes, A, J, NULL, NULL);

1093:   SNESSetFromOptions(snes);

1095:   if (user.fieldBC) {DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u);}
1096:   else              {DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);}
1097:   if (user.restart) {
1098: #if defined(PETSC_HAVE_HDF5)
1099:     PetscViewer viewer;

1101:     PetscViewerCreate(PETSC_COMM_WORLD, &viewer);
1102:     PetscViewerSetType(viewer, PETSCVIEWERHDF5);
1103:     PetscViewerFileSetMode(viewer, FILE_MODE_READ);
1104:     PetscViewerFileSetName(viewer, user.filename);
1105:     PetscViewerHDF5PushGroup(viewer, "/fields");
1106:     VecLoad(u, viewer);
1107:     PetscViewerHDF5PopGroup(viewer);
1108:     PetscViewerDestroy(&viewer);
1109: #endif
1110:   }
1111:   if (user.showInitial) {
1112:     Vec lv;
1113:     DMGetLocalVector(dm, &lv);
1114:     DMGlobalToLocalBegin(dm, u, INSERT_VALUES, lv);
1115:     DMGlobalToLocalEnd(dm, u, INSERT_VALUES, lv);
1116:     DMPrintLocalVec(dm, "Local function", 1.0e-10, lv);
1117:     DMRestoreLocalVector(dm, &lv);
1118:   }
1119:   if (user.viewHierarchy) {
1120:     SNES      lsnes;
1121:     KSP       ksp;
1122:     PC        pc;
1123:     PetscInt  numLevels, l;
1124:     PetscBool isMG;

1126:     PetscObjectTypeCompare((PetscObject) snes, SNESFAS, &isFAS);
1127:     if (isFAS) {
1128:       SNESFASGetLevels(snes, &numLevels);
1129:       for (l = 0; l < numLevels; ++l) {
1130:         SNESFASGetCycleSNES(snes, l, &lsnes);
1131:         SNESMonitorSet(lsnes, SNESMonitorError, &user, NULL);
1132:       }
1133:     } else {
1134:       SNESGetKSP(snes, &ksp);
1135:       KSPGetPC(ksp, &pc);
1136:       PetscObjectTypeCompare((PetscObject) pc, PCMG, &isMG);
1137:       if (isMG) {
1138:         user.pcmg = pc;
1139:         PCMGGetLevels(pc, &numLevels);
1140:         for (l = 0; l < numLevels; ++l) {
1141:           PCMGGetSmootherDown(pc, l, &ksp);
1142:           KSPMonitorSet(ksp, MonitorError, &user, NULL);
1143:         }
1144:       }
1145:     }
1146:   }
1147:   if (user.runType == RUN_FULL || user.runType == RUN_EXACT) {
1148:     PetscErrorCode (*initialGuess[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {zero};

1150:     if (user.nonzInit) initialGuess[0] = ecks;
1151:     if (user.runType == RUN_FULL) {
1152:       DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
1153:     }
1154:     if (user.debug) {
1155:       PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1156:       VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1157:     }
1158:     VecViewFromOptions(u, NULL, "-guess_vec_view");
1159:     SNESSolve(snes, NULL, u);
1160:     SNESGetSolution(snes, &u);
1161:     SNESGetDM(snes, &dm);

1163:     if (user.showSolution) {
1164:       PetscPrintf(PETSC_COMM_WORLD, "Solution\n");
1165:       VecChop(u, 3.0e-9);
1166:       VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1167:     }
1168:   } else if (user.runType == RUN_PERF) {
1169:     Vec       r;
1170:     PetscReal res = 0.0;

1172:     SNESGetFunction(snes, &r, NULL, NULL);
1173:     SNESComputeFunction(snes, u, r);
1174:     PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1175:     VecChop(r, 1.0e-10);
1176:     VecNorm(r, NORM_2, &res);
1177:     PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res);
1178:   } else {
1179:     Vec       r;
1180:     PetscReal res = 0.0, tol = 1.0e-11;

1182:     /* Check discretization error */
1183:     SNESGetFunction(snes, &r, NULL, NULL);
1184:     PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1185:     if (!user.quiet) {VecView(u, PETSC_VIEWER_STDOUT_WORLD);}
1186:     DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
1187:     if (error < tol) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < %2.1e\n", (double)tol);}
1188:     else             {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", (double)error);}
1189:     /* Check residual */
1190:     SNESComputeFunction(snes, u, r);
1191:     PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1192:     VecChop(r, 1.0e-10);
1193:     if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1194:     VecNorm(r, NORM_2, &res);
1195:     PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res);
1196:     /* Check Jacobian */
1197:     {
1198:       Vec b;

1200:       SNESComputeJacobian(snes, u, A, A);
1201:       VecDuplicate(u, &b);
1202:       VecSet(r, 0.0);
1203:       SNESComputeFunction(snes, r, b);
1204:       MatMult(A, u, r);
1205:       VecAXPY(r, 1.0, b);
1206:       PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n");
1207:       VecChop(r, 1.0e-10);
1208:       if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1209:       VecNorm(r, NORM_2, &res);
1210:       PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", (double)res);
1211:       /* check solver */
1212:       if (user.checkksp) {
1213:         KSP ksp;

1215:         if (nullSpace) {
1216:           MatNullSpaceRemove(nullSpace, u);
1217:         }
1218:         SNESComputeJacobian(snes, u, A, J);
1219:         MatMult(A, u, b);
1220:         SNESGetKSP(snes, &ksp);
1221:         KSPSetOperators(ksp, A, J);
1222:         KSPSolve(ksp, b, r);
1223:         VecAXPY(r, -1.0, u);
1224:         VecNorm(r, NORM_2, &res);
1225:         PetscPrintf(PETSC_COMM_WORLD, "KSP Error: %g\n", (double)res);
1226:       }
1227:       VecDestroy(&b);
1228:     }
1229:   }
1230:   VecViewFromOptions(u, NULL, "-vec_view");
1231:   {
1232:     Vec nu;

1234:     DMGetAuxiliaryVec(dm, NULL, 0, &nu);
1235:     if (nu) {VecViewFromOptions(nu, NULL, "-coeff_view");}
1236:   }

1238:   if (user.bdIntegral) {
1239:     DMLabel   label;
1240:     PetscInt  id = 1;
1241:     PetscScalar bdInt = 0.0;
1242:     PetscReal   exact = 3.3333333333;

1244:     DMGetLabel(dm, "marker", &label);
1245:     DMPlexComputeBdIntegral(dm, u, label, 1, &id, bd_integral_2d, &bdInt, NULL);
1246:     PetscPrintf(PETSC_COMM_WORLD, "Solution boundary integral: %.4g\n", (double) PetscAbsScalar(bdInt));
1247:     if (PetscAbsReal(PetscAbsScalar(bdInt) - exact) > PETSC_SQRT_MACHINE_EPSILON) SETERRQ2(PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Invalid boundary integral %g != %g", (double) PetscAbsScalar(bdInt), (double)exact);
1248:   }

1250:   MatNullSpaceDestroy(&nullSpace);
1251:   if (user.jacobianMF) {VecDestroy(&userJ.u);}
1252:   if (A != J) {MatDestroy(&A);}
1253:   MatDestroy(&J);
1254:   VecDestroy(&u);
1255:   SNESDestroy(&snes);
1256:   DMDestroy(&dm);
1257:   PetscFree2(user.exactFuncs, user.exactFields);
1258:   PetscFree(user.kgrid);
1259:   PetscFinalize();
1260:   return ierr;
1261: }

1263: /*TEST
1264:   # 2D serial P1 test 0-4
1265:   test:
1266:     suffix: 2d_p1_0
1267:     requires: triangle
1268:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1270:   test:
1271:     suffix: 2d_p1_1
1272:     requires: triangle
1273:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1275:   test:
1276:     suffix: 2d_p1_2
1277:     requires: triangle
1278:     args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1280:   test:
1281:     suffix: 2d_p1_neumann_0
1282:     requires: triangle
1283:     args: -run_type test -refinement_limit 0.0    -bc_type neumann   -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail

1285:   test:
1286:     suffix: 2d_p1_neumann_1
1287:     requires: triangle
1288:     args: -run_type test -refinement_limit 0.0625 -bc_type neumann   -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1290:   # 2D serial P2 test 5-8
1291:   test:
1292:     suffix: 2d_p2_0
1293:     requires: triangle
1294:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1296:   test:
1297:     suffix: 2d_p2_1
1298:     requires: triangle
1299:     args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1301:   test:
1302:     suffix: 2d_p2_neumann_0
1303:     requires: triangle
1304:     args: -run_type test -refinement_limit 0.0    -bc_type neumann   -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail

1306:   test:
1307:     suffix: 2d_p2_neumann_1
1308:     requires: triangle
1309:     args: -run_type test -refinement_limit 0.0625 -bc_type neumann   -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail

1311:   test:
1312:     suffix: bd_int_0
1313:     requires: triangle
1314:     args: -run_type test -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet

1316:   test:
1317:     suffix: bd_int_1
1318:     requires: triangle
1319:     args: -run_type test -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet

1321:   # 3D serial P1 test 9-12
1322:   test:
1323:     suffix: 3d_p1_0
1324:     requires: ctetgen
1325:     args: -run_type test -dim 3 -refinement_limit 0.0    -bc_type dirichlet -interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1

1327:   test:
1328:     suffix: 3d_p1_1
1329:     requires: ctetgen
1330:     args: -run_type test -dim 3 -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1

1332:   test:
1333:     suffix: 3d_p1_2
1334:     requires: ctetgen
1335:     args: -run_type test -dim 3 -refinement_limit 0.0125 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1

1337:   test:
1338:     suffix: 3d_p1_neumann_0
1339:     requires: ctetgen
1340:     args: -run_type test -dim 3 -bc_type neumann   -interpolate 1 -petscspace_degree 1 -snes_fd -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1

1342:   # Analytic variable coefficient 13-20
1343:   test:
1344:     suffix: 13
1345:     requires: triangle
1346:     args: -run_type test -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1347:   test:
1348:     suffix: 14
1349:     requires: triangle
1350:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1351:   test:
1352:     suffix: 15
1353:     requires: triangle
1354:     args: -run_type test -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1355:   test:
1356:     suffix: 16
1357:     requires: triangle
1358:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1359:   test:
1360:     suffix: 17
1361:     requires: ctetgen
1362:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1364:   test:
1365:     suffix: 18
1366:     requires: ctetgen
1367:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1369:   test:
1370:     suffix: 19
1371:     requires: ctetgen
1372:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1374:   test:
1375:     suffix: 20
1376:     requires: ctetgen
1377:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1379:   # P1 variable coefficient 21-28
1380:   test:
1381:     suffix: 21
1382:     requires: triangle
1383:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1385:   test:
1386:     suffix: 22
1387:     requires: triangle
1388:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1390:   test:
1391:     suffix: 23
1392:     requires: triangle
1393:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1395:   test:
1396:     suffix: 24
1397:     requires: triangle
1398:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1400:   test:
1401:     suffix: 25
1402:     requires: ctetgen
1403:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1405:   test:
1406:     suffix: 26
1407:     requires: ctetgen
1408:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1410:   test:
1411:     suffix: 27
1412:     requires: ctetgen
1413:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1415:   test:
1416:     suffix: 28
1417:     requires: ctetgen
1418:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1420:   # P0 variable coefficient 29-36
1421:   test:
1422:     suffix: 29
1423:     requires: triangle
1424:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1426:   test:
1427:     suffix: 30
1428:     requires: triangle
1429:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1431:   test:
1432:     suffix: 31
1433:     requires: triangle
1434:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1436:   test:
1437:     requires: triangle
1438:     suffix: 32
1439:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1441:   test:
1442:     requires: ctetgen
1443:     suffix: 33
1444:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1446:   test:
1447:     suffix: 34
1448:     requires: ctetgen
1449:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1451:   test:
1452:     suffix: 35
1453:     requires: ctetgen
1454:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1456:   test:
1457:     suffix: 36
1458:     requires: ctetgen
1459:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1461:   # Full solve 39-44
1462:   test:
1463:     suffix: 39
1464:     requires: triangle !single
1465:     args: -run_type full -refinement_limit 0.015625 -interpolate 1 -petscspace_degree 2 -pc_type gamg -ksp_rtol 1.0e-10 -ksp_monitor_short -ksp_converged_reason -snes_monitor_short -snes_converged_reason ::ascii_info_detail
1466:   test:
1467:     suffix: 40
1468:     requires: triangle !single
1469:     args: -run_type full -refinement_limit 0.015625 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 2 -pc_type svd -ksp_rtol 1.0e-10 -snes_monitor_short -snes_converged_reason ::ascii_info_detail
1470:   test:
1471:     suffix: 41
1472:     requires: triangle !single
1473:     args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 2 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -dm_refine_hierarchy 1 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short
1474:   test:
1475:     suffix: 42
1476:     requires: triangle !single
1477:     args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -dm_refine_hierarchy 2 -dm_plex_print_fem 0 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short -fas_levels_2_snes_type newtonls -fas_levels_2_pc_type svd -fas_levels_2_ksp_rtol 1.0e-10 -fas_levels_2_snes_atol 1.0e-11 -fas_levels_2_snes_monitor_short
1478:   test:
1479:     suffix: 43
1480:     requires: triangle !single
1481:     nsize: 2
1482:     args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 2 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -dm_refine_hierarchy 1 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short

1484:   test:
1485:     suffix: 44
1486:     requires: triangle !single
1487:     nsize: 2
1488:     args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short  -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -dm_refine_hierarchy 2 -dm_plex_print_fem 0 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short -fas_levels_2_snes_type newtonls -fas_levels_2_pc_type svd -fas_levels_2_ksp_rtol 1.0e-10 -fas_levels_2_snes_atol 1.0e-11 -fas_levels_2_snes_monitor_short

1490:   # These tests use a loose tolerance just to exercise the PtAP operations for MATIS and multiple PCBDDC setup calls inside PCMG
1491:   testset:
1492:     requires: triangle !single
1493:     nsize: 3
1494:     args: -interpolate -run_type full -petscspace_degree 1 -dm_mat_type is -pc_type mg -pc_mg_levels 2 -mg_coarse_pc_type bddc -pc_mg_galerkin pmat -ksp_rtol 1.0e-2 -snes_converged_reason -dm_refine_hierarchy 2 -snes_max_it 4
1495:     test:
1496:       suffix: gmg_bddc
1497:       filter: sed -e "s/CONVERGED_FNORM_RELATIVE iterations 3/CONVERGED_FNORM_RELATIVE iterations 4/g"
1498:       args: -mg_levels_pc_type jacobi
1499:     test:
1500:       filter: sed -e "s/iterations [0-4]/iterations 4/g"
1501:       suffix: gmg_bddc_lev
1502:       args: -mg_levels_pc_type bddc

1504:   # Restarting
1505:   testset:
1506:     suffix: restart
1507:     requires: hdf5 triangle !complex
1508:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 1
1509:     test:
1510:       args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1511:     test:
1512:       args: -f sol.h5 -restart

1514:   # Periodicity
1515:   test:
1516:     suffix: periodic_0
1517:     requires: triangle
1518:     args: -run_type full -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -snes_converged_reason ::ascii_info_detail

1520:   test:
1521:     requires: !complex
1522:     suffix: periodic_1
1523:     args: -quiet -run_type test -simplex 0 -x_periodicity periodic -y_periodicity periodic -vec_view vtk:test.vtu:vtk_vtu -interpolate 1 -petscspace_degree 1 -dm_refine 1

1525:   # 2D serial P1 test with field bc
1526:   test:
1527:     suffix: field_bc_2d_p1_0
1528:     requires: triangle
1529:     args: -run_type test              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1531:   test:
1532:     suffix: field_bc_2d_p1_1
1533:     requires: triangle
1534:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1536:   test:
1537:     suffix: field_bc_2d_p1_neumann_0
1538:     requires: triangle
1539:     args: -run_type test              -interpolate 1 -bc_type neumann   -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1541:   test:
1542:     suffix: field_bc_2d_p1_neumann_1
1543:     requires: triangle
1544:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1546:   # 3D serial P1 test with field bc
1547:   test:
1548:     suffix: field_bc_3d_p1_0
1549:     requires: ctetgen
1550:     args: -run_type test -dim 3              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1552:   test:
1553:     suffix: field_bc_3d_p1_1
1554:     requires: ctetgen
1555:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1557:   test:
1558:     suffix: field_bc_3d_p1_neumann_0
1559:     requires: ctetgen
1560:     args: -run_type test -dim 3              -interpolate 1 -bc_type neumann   -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1562:   test:
1563:     suffix: field_bc_3d_p1_neumann_1
1564:     requires: ctetgen
1565:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1567:   # 2D serial P2 test with field bc
1568:   test:
1569:     suffix: field_bc_2d_p2_0
1570:     requires: triangle
1571:     args: -run_type test              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1573:   test:
1574:     suffix: field_bc_2d_p2_1
1575:     requires: triangle
1576:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1578:   test:
1579:     suffix: field_bc_2d_p2_neumann_0
1580:     requires: triangle
1581:     args: -run_type test              -interpolate 1 -bc_type neumann   -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1583:   test:
1584:     suffix: field_bc_2d_p2_neumann_1
1585:     requires: triangle
1586:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1588:   # 3D serial P2 test with field bc
1589:   test:
1590:     suffix: field_bc_3d_p2_0
1591:     requires: ctetgen
1592:     args: -run_type test -dim 3              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1594:   test:
1595:     suffix: field_bc_3d_p2_1
1596:     requires: ctetgen
1597:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1599:   test:
1600:     suffix: field_bc_3d_p2_neumann_0
1601:     requires: ctetgen
1602:     args: -run_type test -dim 3              -interpolate 1 -bc_type neumann   -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1604:   test:
1605:     suffix: field_bc_3d_p2_neumann_1
1606:     requires: ctetgen
1607:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1609:   # Full solve simplex: Convergence
1610:   test:
1611:     suffix: 3d_p1_conv
1612:     requires: ctetgen
1613:     args: -run_type full -dim 3 -cells 1,1,1 -dm_refine 1 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 \
1614:       -snes_convergence_estimate -convest_num_refine 1 -pc_type lu

1616:   # Full solve simplex: PCBDDC
1617:   test:
1618:     suffix: tri_bddc
1619:     requires: triangle !single
1620:     nsize: 5
1621:     args: -run_type full -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -dm_mat_type is -pc_type bddc -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0

1623:   # Full solve simplex: PCBDDC
1624:   test:
1625:     suffix: tri_parmetis_bddc
1626:     requires: triangle !single parmetis
1627:     nsize: 4
1628:     args: -run_type full -petscpartitioner_type parmetis -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -dm_mat_type is -pc_type bddc -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0

1630:   testset:
1631:     args: -run_type full -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -dm_mat_type is -pc_type bddc -ksp_type gmres -snes_monitor_short -ksp_monitor_short -snes_view -simplex 0 -petscspace_poly_tensor -pc_bddc_corner_selection -cells 3,3 -ksp_rtol 1.e-9 -pc_bddc_use_edges 0
1632:     nsize: 5
1633:     output_file: output/ex12_quad_bddc.out
1634:     filter: sed -e "s/aijcusparse/aij/g" -e "s/aijviennacl/aij/g" -e "s/factorization: cusparse/factorization: petsc/g"
1635:     test:
1636:       requires: !single
1637:       suffix: quad_bddc
1638:     test:
1639:       requires: !single cuda
1640:       suffix: quad_bddc_cuda
1641:       args: -matis_localmat_type aijcusparse -pc_bddc_dirichlet_pc_factor_mat_solver_type cusparse -pc_bddc_neumann_pc_factor_mat_solver_type cusparse
1642:     test:
1643:       requires: !single viennacl
1644:       suffix: quad_bddc_viennacl
1645:       args: -matis_localmat_type aijviennacl

1647:   # Full solve simplex: ASM
1648:   test:
1649:     suffix: tri_q2q1_asm_lu
1650:     requires: triangle !single
1651:     args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type asm -pc_asm_type restrict -pc_asm_blocks 4 -sub_pc_type lu -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0

1653:   test:
1654:     suffix: tri_q2q1_msm_lu
1655:     requires: triangle !single
1656:     args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type asm -pc_asm_type restrict -pc_asm_local_type multiplicative -pc_asm_blocks 4 -sub_pc_type lu -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0

1658:   test:
1659:     suffix: tri_q2q1_asm_sor
1660:     requires: triangle !single
1661:     args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type asm -pc_asm_type restrict -pc_asm_blocks 4 -sub_pc_type sor -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0

1663:   test:
1664:     suffix: tri_q2q1_msm_sor
1665:     requires: triangle !single
1666:     args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type asm -pc_asm_type restrict -pc_asm_local_type multiplicative -pc_asm_blocks 4 -sub_pc_type sor -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0

1668:   # Full solve simplex: FAS
1669:   test:
1670:     suffix: fas_newton_0
1671:     requires: triangle !single
1672:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 2 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -dm_refine_hierarchy 1 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short

1674:   test:
1675:     suffix: fas_newton_1
1676:     requires: triangle !single
1677:     args: -run_type full -dm_refine_hierarchy 3 -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type lu -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_snes_linesearch_type basic -fas_levels_ksp_rtol 1.0e-10 -fas_levels_snes_monitor_short
1678:     filter: sed -e "s/total number of linear solver iterations=14/total number of linear solver iterations=15/g"

1680:   test:
1681:     suffix: fas_ngs_0
1682:     requires: triangle !single
1683:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 2 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -dm_refine_hierarchy 1 -snes_view -fas_levels_1_snes_type ngs -fas_levels_1_snes_monitor_short

1685:   test:
1686:     suffix: fas_newton_coarse_0
1687:     requires: pragmatic triangle
1688:     TODO: broken
1689:     args: -run_type full -dm_refine 2 -dm_plex_hash_location -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -dm_coarsen_hierarchy 1 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short

1691:   test:
1692:     suffix: mg_newton_coarse_0
1693:     requires: triangle pragmatic
1694:     TODO: broken
1695:     args: -run_type full -dm_refine 3 -interpolate 1 -petscspace_degree 1 -snes_monitor_short -ksp_monitor_true_residual -snes_converged_reason ::ascii_info_detail -dm_coarsen_hierarchy 3 -dm_plex_hash_location -snes_view -dm_view -ksp_type richardson -pc_type mg  -pc_mg_levels 4 -snes_atol 1.0e-8 -ksp_atol 1.0e-8 -snes_rtol 0.0 -ksp_rtol 0.0 -ksp_norm_type unpreconditioned -mg_levels_ksp_type gmres -mg_levels_pc_type ilu -mg_levels_ksp_max_it 10

1697:   test:
1698:     suffix: mg_newton_coarse_1
1699:     requires: triangle pragmatic
1700:     TODO: broken
1701:     args: -run_type full -dm_refine 5 -interpolate 1 -petscspace_degree 1 -dm_coarsen_hierarchy 5 -dm_plex_hash_location -dm_plex_separate_marker -dm_plex_coarsen_bd_label marker -dm_plex_remesh_bd -ksp_type richardson -ksp_rtol 1.0e-12 -pc_type mg -pc_mg_levels 3 -mg_levels_ksp_max_it 2 -snes_converged_reason ::ascii_info_detail -snes_monitor -ksp_monitor_true_residual -mg_levels_ksp_monitor_true_residual -dm_view -ksp_view

1703:   test:
1704:     suffix: mg_newton_coarse_2
1705:     requires: triangle pragmatic
1706:     TODO: broken
1707:     args: -run_type full -dm_refine 5 -interpolate 1 -petscspace_degree 1 -dm_coarsen_hierarchy 5 -dm_plex_hash_location -dm_plex_separate_marker -dm_plex_remesh_bd -ksp_type richardson -ksp_rtol 1.0e-12 -pc_type mg -pc_mg_levels 3 -mg_levels_ksp_max_it 2 -snes_converged_reason ::ascii_info_detail -snes_monitor -ksp_monitor_true_residual -mg_levels_ksp_monitor_true_residual -dm_view -ksp_view

1709:   # Full solve tensor
1710:   test:
1711:     suffix: tensor_plex_2d
1712:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_refine_hierarchy 2 -cells 2,2

1714:   test:
1715:     suffix: tensor_p4est_2d
1716:     requires: p4est
1717:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_forest_initial_refinement 2 -dm_forest_minimum_refinement 0 -dm_plex_convert_type p4est -cells 2,2

1719:   test:
1720:     suffix: tensor_plex_3d
1721:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dim 3 -dm_refine_hierarchy 1 -cells 2,2,2

1723:   test:
1724:     suffix: tensor_p4est_3d
1725:     requires: p4est
1726:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_forest_initial_refinement 1 -dm_forest_minimum_refinement 0 -dim 3 -dm_plex_convert_type p8est -cells 2,2,2

1728:   test:
1729:     suffix: p4est_test_q2_conformal_serial
1730:     requires: p4est
1731:     args: -run_type test -interpolate 1 -petscspace_degree 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -cells 2,2

1733:   test:
1734:     suffix: p4est_test_q2_conformal_parallel
1735:     requires: p4est
1736:     nsize: 7
1737:     args: -run_type test -interpolate 1 -petscspace_degree 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type simple -cells 2,2

1739:   test:
1740:     suffix: p4est_test_q2_conformal_parallel_parmetis
1741:     requires: parmetis p4est
1742:     nsize: 4
1743:     args: -run_type test -interpolate 1 -petscspace_degree 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type parmetis -cells 2,2

1745:   test:
1746:     suffix: p4est_test_q2_nonconformal_serial
1747:     requires: p4est
1748:     filter: grep -v "CG or CGNE: variant"
1749:     args: -run_type test -interpolate 1 -petscspace_degree 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2

1751:   test:
1752:     suffix: p4est_test_q2_nonconformal_parallel
1753:     requires: p4est
1754:     filter: grep -v "CG or CGNE: variant"
1755:     nsize: 7
1756:     args: -run_type test -interpolate 1 -petscspace_degree 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -cells 2,2

1758:   test:
1759:     suffix: p4est_test_q2_nonconformal_parallel_parmetis
1760:     requires: parmetis p4est
1761:     nsize: 4
1762:     args: -run_type test -interpolate 1 -petscspace_degree 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type parmetis -cells 2,2

1764:   test:
1765:     suffix: p4est_exact_q2_conformal_serial
1766:     requires: p4est !single !complex !__float128
1767:     args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type none -fas_levels_ksp_type preonly -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -cells 2,2

1769:   test:
1770:     suffix: p4est_exact_q2_conformal_parallel
1771:     requires: p4est !single !complex !__float128
1772:     nsize: 4
1773:     args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type none -fas_levels_ksp_type preonly -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -cells 2,2

1775:   test:
1776:     suffix: p4est_exact_q2_conformal_parallel_parmetis
1777:     requires: parmetis p4est !single
1778:     nsize: 4
1779:     args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type none -fas_levels_ksp_type preonly -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type parmetis  -cells 2,2

1781:   test:
1782:     suffix: p4est_exact_q2_nonconformal_serial
1783:     requires: p4est
1784:     args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type none -fas_levels_ksp_type preonly -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2

1786:   test:
1787:     suffix: p4est_exact_q2_nonconformal_parallel
1788:     requires: p4est
1789:     nsize: 7
1790:     args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type none -fas_levels_ksp_type preonly -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -cells 2,2

1792:   test:
1793:     suffix: p4est_exact_q2_nonconformal_parallel_parmetis
1794:     requires: parmetis p4est
1795:     nsize: 4
1796:     args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type none -fas_levels_ksp_type preonly -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type parmetis -cells 2,2

1798:   test:
1799:     suffix: p4est_full_q2_nonconformal_serial
1800:     requires: p4est !single
1801:     filter: grep -v "variant HERMITIAN"
1802:     args: -run_type full -interpolate 1 -petscspace_degree 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type jacobi -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type jacobi -fas_levels_ksp_type cg -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2

1804:   test:
1805:     suffix: p4est_full_q2_nonconformal_parallel
1806:     requires: p4est !single
1807:     filter: grep -v "variant HERMITIAN"
1808:     nsize: 7
1809:     args: -run_type full -interpolate 1 -petscspace_degree 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type jacobi -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type jacobi -fas_levels_ksp_type cg -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -cells 2,2

1811:   test:
1812:     suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1813:     requires: p4est !single
1814:     filter: grep -v "variant HERMITIAN"
1815:     nsize: 7
1816:     args: -run_type full -interpolate 1 -petscspace_degree 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -dm_mat_type is -fas_coarse_pc_type bddc -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type bddc -fas_levels_ksp_type cg -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -cells 2,2

1818:   test:
1819:     suffix: p4est_full_q2_nonconformal_parallel_bddc
1820:     requires: p4est !single
1821:     filter: grep -v "variant HERMITIAN"
1822:     nsize: 7
1823:     args: -run_type full -interpolate 1 -petscspace_degree 2 -snes_max_it 20 -snes_type newtonls -dm_mat_type is -pc_type bddc -ksp_type cg -snes_monitor_short -snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -cells 2,2

1825:   test:
1826:     TODO: broken
1827:     suffix: p4est_fas_q2_conformal_serial
1828:     requires: p4est !complex !__float128
1829:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type gmres -fas_coarse_pc_type svd -fas_coarse_ksp_type gmres -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type svd -fas_levels_ksp_type gmres -fas_levels_snes_monitor_short -simplex 0 -dm_refine_hierarchy 3 -cells 2,2

1831:   test:
1832:     TODO: broken
1833:     suffix: p4est_fas_q2_nonconformal_serial
1834:     requires: p4est
1835:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type gmres -fas_coarse_pc_type jacobi -fas_coarse_ksp_type gmres -fas_coarse_ksp_monitor_true_residual -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type jacobi -fas_levels_ksp_type gmres -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2

1837:   test:
1838:     suffix: fas_newton_0_p4est
1839:     requires: p4est !single !__float128
1840:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 2 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2

1842:   # Full solve simplicial AMR
1843:   test:
1844:     suffix: tri_p1_adapt_0
1845:     requires: pragmatic
1846:     TODO: broken
1847:     args: -run_type exact -dim 2 -dm_refine 5 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -variable_coefficient circle -snes_converged_reason ::ascii_info_detail -pc_type lu -adaptor_refinement_factor 1.0 -dm_view -dm_adapt_view -snes_adapt_initial 1

1849:   test:
1850:     suffix: tri_p1_adapt_1
1851:     requires: pragmatic
1852:     TODO: broken
1853:     args: -run_type exact -dim 2 -dm_refine 5 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -variable_coefficient circle -snes_converged_reason ::ascii_info_detail -pc_type lu -adaptor_refinement_factor 1.0 -dm_view -dm_adapt_iter_view -dm_adapt_view -snes_adapt_sequence 2

1855:   test:
1856:     suffix: tri_p1_adapt_analytic_0
1857:     requires: pragmatic
1858:     TODO: broken
1859:     args: -run_type exact -dim 2 -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -adaptor_monitor -dm_view -dm_adapt_iter_view

1861:   # Full solve tensor AMR
1862:   test:
1863:     suffix: quad_q1_adapt_0
1864:     requires: p4est
1865:     args: -run_type exact -dim 2 -simplex 0 -dm_plex_convert_type p4est -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -variable_coefficient circle -snes_converged_reason ::ascii_info_detail -pc_type lu -dm_forest_initial_refinement 4   -snes_adapt_initial 1 -dm_view
1866:     filter: grep -v DM_

1868:   test:
1869:     suffix: amr_0
1870:     nsize: 5
1871:     args: -run_type test -petscpartitioner_type simple -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_refine 1 -cells 2,2

1873:   test:
1874:     suffix: amr_1
1875:     requires: p4est !complex
1876:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_plex_convert_type p4est -dm_p4est_refine_pattern center -dm_forest_maximum_refinement 5 -dm_view vtk:amr.vtu:vtk_vtu -vec_view vtk:amr.vtu:vtk_vtu:append -cells 2,2

1878:   test:
1879:     suffix: p4est_solve_bddc
1880:     requires: p4est !complex
1881:     args: -run_type full -variable_coefficient nonlinear -nonzero_initial_guess 1 -interpolate 1 -petscspace_degree 2 -snes_max_it 20 -snes_type newtonls -dm_mat_type is -pc_type bddc -ksp_type cg -snes_monitor_short -ksp_monitor -snes_linesearch_type bt -snes_converged_reason -snes_view -simplex 0 -petscspace_poly_tensor -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -pc_bddc_detect_disconnected
1882:     nsize: 4

1884:   test:
1885:     suffix: p4est_solve_fas
1886:     requires: p4est
1887:     args: -run_type full -variable_coefficient nonlinear -nonzero_initial_guess 1 -interpolate 1 -petscspace_degree 2 -snes_max_it 10 -snes_type fas -snes_linesearch_type bt -snes_fas_levels 3 -fas_coarse_snes_type newtonls -fas_coarse_snes_linesearch_type basic -fas_coarse_ksp_type cg -fas_coarse_pc_type jacobi -fas_coarse_snes_monitor_short -fas_levels_snes_max_it 4 -fas_levels_snes_type newtonls -fas_levels_snes_linesearch_type bt -fas_levels_ksp_type cg -fas_levels_pc_type jacobi -fas_levels_snes_monitor_short -fas_levels_cycle_snes_linesearch_type bt -snes_monitor_short -snes_converged_reason -snes_view -simplex 0 -petscspace_poly_tensor -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash
1888:     nsize: 4
1889:     TODO: identical machine two runs produce slightly different solver trackers

1891:   test:
1892:     suffix: p4est_convergence_test_1
1893:     requires: p4est
1894:     args:  -quiet -run_type test -interpolate 1 -petscspace_degree 1 -simplex 0 -petscspace_poly_tensor -dm_plex_convert_type p4est -dm_forest_minimum_refinement 2 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash
1895:     nsize: 4

1897:   test:
1898:     suffix: p4est_convergence_test_2
1899:     requires: p4est
1900:     args: -quiet -run_type test -interpolate 1 -petscspace_degree 1 -simplex 0 -petscspace_poly_tensor -dm_plex_convert_type p4est -dm_forest_minimum_refinement 3 -dm_forest_initial_refinement 3 -dm_forest_maximum_refinement 5 -dm_p4est_refine_pattern hash

1902:   test:
1903:     suffix: p4est_convergence_test_3
1904:     requires: p4est
1905:     args: -quiet -run_type test -interpolate 1 -petscspace_degree 1 -simplex 0 -petscspace_poly_tensor -dm_plex_convert_type p4est -dm_forest_minimum_refinement 4 -dm_forest_initial_refinement 4 -dm_forest_maximum_refinement 6 -dm_p4est_refine_pattern hash

1907:   test:
1908:     suffix: p4est_convergence_test_4
1909:     requires: p4est
1910:     args: -quiet -run_type test -interpolate 1 -petscspace_degree 1 -simplex 0 -petscspace_poly_tensor -dm_plex_convert_type p4est -dm_forest_minimum_refinement 5 -dm_forest_initial_refinement 5 -dm_forest_maximum_refinement 7 -dm_p4est_refine_pattern hash
1911:     timeoutfactor: 5

1913:   # Serial tests with GLVis visualization
1914:   test:
1915:     suffix: glvis_2d_tet_p1
1916:     args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -f ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_plex_gmsh_periodic 0
1917:   test:
1918:     suffix: glvis_2d_tet_p2
1919:     args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -f ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_plex_gmsh_periodic 0
1920:   test:
1921:     suffix: glvis_2d_hex_p1
1922:     args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -simplex 0 -dm_refine 1
1923:   test:
1924:     suffix: glvis_2d_hex_p2
1925:     args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -simplex 0 -dm_refine 1
1926:   test:
1927:     suffix: glvis_2d_hex_p2_p4est
1928:     requires: p4est
1929:     args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 1 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2 -viewer_glvis_dm_plex_enable_ncmesh
1930:   test:
1931:     suffix: glvis_2d_tet_p0
1932:     args: -run_type exact  -interpolate 1 -guess_vec_view glvis: -nonzero_initial_guess 1 -f ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -petscspace_degree 0
1933:   test:
1934:     suffix: glvis_2d_hex_p0
1935:     args: -run_type exact  -interpolate 1 -guess_vec_view glvis: -nonzero_initial_guess 1 -cells 5,7  -simplex 0 -petscspace_degree 0

1937:   # PCHPDDM tests
1938:   testset:
1939:     nsize: 4
1940:     requires: hpddm slepc !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1941:     args: -run_type test -run_test_check_ksp -quiet -petscspace_degree 1 -interpolate 1 -petscpartitioner_type simple -bc_type none -simplex 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 2 -pc_hpddm_coarse_p 1 -pc_hpddm_coarse_pc_type svd -ksp_rtol 1.e-10 -pc_hpddm_levels_1_st_pc_factor_shift_type INBLOCKS -ksp_converged_reason
1942:     test:
1943:       suffix: quad_singular_hpddm
1944:       args: -cells 6,7
1945:     test:
1946:       requires: p4est
1947:       suffix: p4est_singular_2d_hpddm
1948:       args: -dm_plex_convert_type p4est -dm_forest_minimum_refinement 1 -dm_forest_initial_refinement 3 -dm_forest_maximum_refinement 3
1949:     test:
1950:       requires: p4est
1951:       suffix: p4est_nc_singular_2d_hpddm
1952:       args: -dm_plex_convert_type p4est -dm_forest_minimum_refinement 1 -dm_forest_initial_refinement 1 -dm_forest_maximum_refinement 3 -dm_p4est_refine_pattern hash
1953:   testset:
1954:     nsize: 4
1955:     requires: hpddm slepc triangle !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1956:     args: -run_type full -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 4 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-1
1957:     test:
1958:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1959:       suffix: tri_hpddm_reuse_baij
1960:     test:
1961:       requires: !complex
1962:       suffix: tri_hpddm_reuse
1963:   testset:
1964:     nsize: 4
1965:     requires: hpddm slepc !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1966:     args: -run_type full -petscpartitioner_type simple -cells 7,5 -dm_refine 2 -simplex 0 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 4 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-1
1967:     test:
1968:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1969:       suffix: quad_hpddm_reuse_baij
1970:     test:
1971:       requires: !complex
1972:       suffix: quad_hpddm_reuse
1973:   testset:
1974:     nsize: 4
1975:     requires: hpddm slepc !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1976:     args: -run_type full -petscpartitioner_type simple -cells 7,5 -dm_refine 2 -simplex 0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_threshold 0.1 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-1
1977:     test:
1978:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1979:       suffix: quad_hpddm_reuse_threshold_baij
1980:     test:
1981:       requires: !complex
1982:       suffix: quad_hpddm_reuse_threshold
1983:   testset:
1984:     nsize: 4
1985:     requires: hpddm slepc parmetis !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1986:     filter: sed -e "s/linear solver iterations=17/linear solver iterations=16/g"
1987:     args: -run_type full -petscpartitioner_type parmetis -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -snes_converged_reason ::ascii_info_detail -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type icc -pc_hpddm_levels_1_eps_nev 20 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-10 -f ${PETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -pc_hpddm_levels_1_sub_pc_factor_levels 3 -variable_coefficient circle -dm_plex_gmsh_periodic 0
1988:     test:
1989:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1990:       filter: grep -v "      total: nonzeros=" | grep -v "      rows=" | sed -e "s/total number of linear solver iterations=1[5-7]/total number of linear solver iterations=16/g"
1991:       suffix: tri_parmetis_hpddm_baij
1992:     test:
1993:       filter: grep -v "      total: nonzeros=" | grep -v "      rows=" | sed -e "s/total number of linear solver iterations=1[5-7]/total number of linear solver iterations=16/g"
1994:       requires: !complex
1995:       suffix: tri_parmetis_hpddm

1997:   # 2D serial P1 tests for adaptive MG
1998:   test:
1999:     suffix: 2d_p1_adaptmg_0
2000:     requires: triangle bamg
2001:     args: -dm_refine_hierarchy 3 -cells 4,4 -interpolate -bc_type dirichlet -petscspace_degree 1 \
2002:           -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
2003:           -snes_max_it 1 -ksp_converged_reason \
2004:           -ksp_rtol 1e-8 -pc_type mg
2005:   # -ksp_monitor_true_residual -ksp_converged_reason -mg_levels_ksp_monitor_true_residual -pc_mg_mesp_monitor -dm_adapt_interp_view_fine draw -dm_adapt_interp_view_coarse draw -draw_pause 1
2006:   test:
2007:     suffix: 2d_p1_adaptmg_1
2008:     requires: triangle bamg
2009:     args: -dm_refine_hierarchy 3 -cells 4,4 -interpolate -bc_type dirichlet -petscspace_degree 1 \
2010:           -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
2011:           -snes_max_it 1 -ksp_converged_reason \
2012:           -ksp_rtol 1e-8 -pc_type mg -pc_mg_galerkin -pc_mg_adapt_interp -pc_mg_adapt_interp_coarse_space eigenvector -pc_mg_adapt_interp_n 1 \
2013:             -pc_mg_mesp_ksp_type richardson -pc_mg_mesp_ksp_richardson_self_scale -pc_mg_mesp_ksp_max_it 100 -pc_mg_mesp_pc_type none

2015: TEST*/