Actual source code: ex12.c

petsc-master 2017-11-16
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:  #include <petscdmplex.h>
  8:  #include <petscsnes.h>
  9:  #include <petscds.h>
 10: #include <petscviewerhdf5.h>

 12: typedef enum {NEUMANN, DIRICHLET, NONE} BCType;
 13: typedef enum {RUN_FULL, RUN_EXACT, RUN_TEST, RUN_PERF} RunType;
 14: typedef enum {COEFF_NONE, COEFF_ANALYTIC, COEFF_FIELD, COEFF_NONLINEAR} CoeffType;

 16: typedef struct {
 17:   PetscInt       debug;             /* The debugging level */
 18:   RunType        runType;           /* Whether to run tests, or solve the full problem */
 19:   PetscBool      jacobianMF;        /* Whether to calculate the Jacobian action on the fly */
 20:   PetscLogEvent  createMeshEvent;
 21:   PetscBool      showInitial, showSolution, restart, check, quiet, nonzInit;
 22:   /* Domain and mesh definition */
 23:   PetscInt       dim;               /* The topological mesh dimension */
 24:   DMBoundaryType periodicity[3];    /* The domain periodicity */
 25:   PetscInt       cells[3];          /* The initial domain division */
 26:   char           filename[2048];    /* The optional mesh file */
 27:   PetscBool      interpolate;       /* Generate intermediate mesh elements */
 28:   PetscReal      refinementLimit;   /* The largest allowable cell volume */
 29:   PetscBool      viewHierarchy;     /* Whether to view the hierarchy */
 30:   PetscBool      simplex;           /* Simplicial mesh */
 31:   /* Problem definition */
 32:   BCType         bcType;
 33:   CoeffType      variableCoefficient;
 34:   PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx);
 35:   PetscBool      fieldBC;
 36:   void           (**exactFields)(PetscInt, PetscInt, PetscInt,
 37:                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
 38:                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
 39:                                  PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]);
 40:   /* Solver */
 41:   PC            pcmg;              /* This is needed for error monitoring */
 42: } AppCtx;

 44: static PetscErrorCode zero(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 45: {
 46:   u[0] = 0.0;
 47:   return 0;
 48: }

 50: static PetscErrorCode ecks(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 51: {
 52:   u[0] = x[0];
 53:   return 0;
 54: }

 56: /*
 57:   In 2D for Dirichlet conditions, we use exact solution:

 59:     u = x^2 + y^2
 60:     f = 4

 62:   so that

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

 66:   For Neumann conditions, we have

 68:     -\nabla u \cdot -\hat y |_{y=0} =  (2y)|_{y=0} =  0 (bottom)
 69:     -\nabla u \cdot  \hat y |_{y=1} = -(2y)|_{y=1} = -2 (top)
 70:     -\nabla u \cdot -\hat x |_{x=0} =  (2x)|_{x=0} =  0 (left)
 71:     -\nabla u \cdot  \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)

 73:   Which we can express as

 75:     \nabla u \cdot  \hat n|_\Gamma = {2 x, 2 y} \cdot \hat n = 2 (x + y)
 76: */
 77: static PetscErrorCode quadratic_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 78: {
 79:   *u = x[0]*x[0] + x[1]*x[1];
 80:   return 0;
 81: }

 83: static void quadratic_u_field_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
 84:                                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
 85:                                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
 86:                                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
 87: {
 88:   uexact[0] = a[0];
 89: }

 91: static void f0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
 92:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
 93:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
 94:                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
 95: {
 96:   f0[0] = 4.0;
 97: }

 99: static void f0_bd_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
100:                     const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
101:                     const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
102:                     PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
103: {
104:   PetscInt d;
105:   for (d = 0, f0[0] = 0.0; d < dim; ++d) f0[0] += -n[d]*2.0*x[d];
106: }

108: static void f1_bd_zero(PetscInt dim, PetscInt Nf, PetscInt NfAux,
109:                        const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
110:                        const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
111:                        PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
112: {
113:   PetscInt comp;
114:   for (comp = 0; comp < dim; ++comp) f1[comp] = 0.0;
115: }

117: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
118: static void f1_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
119:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
120:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
121:                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
122: {
123:   PetscInt d;
124:   for (d = 0; d < dim; ++d) f1[d] = u_x[d];
125: }

127: /* < \nabla v, \nabla u + {\nabla u}^T >
128:    This just gives \nabla u, give the perdiagonal for the transpose */
129: static void g3_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
130:                   const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
131:                   const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
132:                   PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
133: {
134:   PetscInt d;
135:   for (d = 0; d < dim; ++d) g3[d*dim+d] = 1.0;
136: }

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

141:     u = sin(2 pi x)
142:     f = -4 pi^2 sin(2 pi x)

144:   so that

146:     -\Delta u + f = 4 pi^2 sin(2 pi x) - 4 pi^2 sin(2 pi x) = 0
147: */
148: static PetscErrorCode xtrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
149: {
150:   *u = PetscSinReal(2.0*PETSC_PI*x[0]);
151:   return 0;
152: }

154: static void f0_xtrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
155:                        const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
156:                        const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
157:                        PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
158: {
159:   f0[0] = -4.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
160: }

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

165:     u = sin(2 pi x) sin(2 pi y)
166:     f = -8 pi^2 sin(2 pi x)

168:   so that

170:     -\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
171: */
172: static PetscErrorCode xytrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
173: {
174:   *u = PetscSinReal(2.0*PETSC_PI*x[0])*PetscSinReal(2.0*PETSC_PI*x[1]);
175:   return 0;
176: }

178: static void f0_xytrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
179:                         const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
180:                         const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
181:                         PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
182: {
183:   f0[0] = -8.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
184: }

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

189:     u  = x^2 + y^2
190:     f  = 6 (x + y)
191:     nu = (x + y)

193:   so that

195:     -\div \nu \grad u + f = -6 (x + y) + 6 (x + y) = 0
196: */
197: static PetscErrorCode nu_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
198: {
199:   *u = x[0] + x[1];
200:   return 0;
201: }

203: void f0_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
204:                    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
205:                    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
206:                    PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
207: {
208:   f0[0] = 6.0*(x[0] + x[1]);
209: }

211: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
212: void f1_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
213:                    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
214:                    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
215:                    PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
216: {
217:   PetscInt d;
218:   for (d = 0; d < dim; ++d) f1[d] = (x[0] + x[1])*u_x[d];
219: }

221: void f1_field_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
222:                 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
223:                 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
224:                 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
225: {
226:   PetscInt d;
227:   for (d = 0; d < dim; ++d) f1[d] = a[0]*u_x[d];
228: }

230: /* < \nabla v, \nabla u + {\nabla u}^T >
231:    This just gives \nabla u, give the perdiagonal for the transpose */
232: void g3_analytic_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
233:                     const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
234:                     const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
235:                     PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
236: {
237:   PetscInt d;
238:   for (d = 0; d < dim; ++d) g3[d*dim+d] = x[0] + x[1];
239: }

241: void g3_field_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
242:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
243:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
244:                  PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
245: {
246:   PetscInt d;
247:   for (d = 0; d < dim; ++d) g3[d*dim+d] = a[0];
248: }

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

253:     u  = x^2 + y^2
254:     f  = 16 (x^2 + y^2)
255:     nu = 1/2 |grad u|^2

257:   so that

259:     -\div \nu \grad u + f = -16 (x^2 + y^2) + 16 (x^2 + y^2) = 0
260: */
261: void f0_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
262:                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
263:                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
264:                              PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
265: {
266:   f0[0] = 16.0*(x[0]*x[0] + x[1]*x[1]);
267: }

269: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
270: void f1_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
271:                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
272:                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
273:                              PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
274: {
275:   PetscScalar nu = 0.0;
276:   PetscInt    d;
277:   for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
278:   for (d = 0; d < dim; ++d) f1[d] = 0.5*nu*u_x[d];
279: }

281: /*
282:   grad (u + eps w) - grad u = eps grad w

284:   1/2 |grad (u + eps w)|^2 grad (u + eps w) - 1/2 |grad u|^2 grad u
285: = 1/2 (|grad u|^2 + 2 eps <grad u,grad w>) (grad u + eps grad w) - 1/2 |grad u|^2 grad u
286: = 1/2 (eps |grad u|^2 grad w + 2 eps <grad u,grad w> grad u)
287: = eps (1/2 |grad u|^2 grad w + grad u <grad u,grad w>)
288: */
289: void g3_analytic_nonlinear_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
290:                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
291:                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
292:                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
293: {
294:   PetscScalar nu = 0.0;
295:   PetscInt    d, e;
296:   for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
297:   for (d = 0; d < dim; ++d) {
298:     g3[d*dim+d] = 0.5*nu;
299:     for (e = 0; e < dim; ++e) {
300:       g3[d*dim+e] += u_x[d]*u_x[e];
301:     }
302:   }
303: }

305: /*
306:   In 3D for Dirichlet conditions we use exact solution:

308:     u = 2/3 (x^2 + y^2 + z^2)
309:     f = 4

311:   so that

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

315:   For Neumann conditions, we have

317:     -\nabla u \cdot -\hat z |_{z=0} =  (2z)|_{z=0} =  0 (bottom)
318:     -\nabla u \cdot  \hat z |_{z=1} = -(2z)|_{z=1} = -2 (top)
319:     -\nabla u \cdot -\hat y |_{y=0} =  (2y)|_{y=0} =  0 (front)
320:     -\nabla u \cdot  \hat y |_{y=1} = -(2y)|_{y=1} = -2 (back)
321:     -\nabla u \cdot -\hat x |_{x=0} =  (2x)|_{x=0} =  0 (left)
322:     -\nabla u \cdot  \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)

324:   Which we can express as

326:     \nabla u \cdot  \hat n|_\Gamma = {2 x, 2 y, 2z} \cdot \hat n = 2 (x + y + z)
327: */
328: static PetscErrorCode quadratic_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
329: {
330:   *u = 2.0*(x[0]*x[0] + x[1]*x[1] + x[2]*x[2])/3.0;
331:   return 0;
332: }

334: static void quadratic_u_field_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
335:                                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
336:                                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
337:                                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
338: {
339:   uexact[0] = a[0];
340: }

342: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
343: {
344:   const char    *bcTypes[3]  = {"neumann", "dirichlet", "none"};
345:   const char    *runTypes[4] = {"full", "exact", "test", "perf"};
346:   const char    *coeffTypes[4] = {"none", "analytic", "field", "nonlinear"};
347:   PetscInt       bd, bc, run, coeff, n;
348:   PetscBool      flg;

352:   options->debug               = 0;
353:   options->runType             = RUN_FULL;
354:   options->dim                 = 2;
355:   options->periodicity[0]      = DM_BOUNDARY_NONE;
356:   options->periodicity[1]      = DM_BOUNDARY_NONE;
357:   options->periodicity[2]      = DM_BOUNDARY_NONE;
358:   options->cells[0]            = 2;
359:   options->cells[1]            = 2;
360:   options->cells[2]            = 2;
361:   options->filename[0]         = '\0';
362:   options->interpolate         = PETSC_FALSE;
363:   options->refinementLimit     = 0.0;
364:   options->bcType              = DIRICHLET;
365:   options->variableCoefficient = COEFF_NONE;
366:   options->fieldBC             = PETSC_FALSE;
367:   options->jacobianMF          = PETSC_FALSE;
368:   options->showInitial         = PETSC_FALSE;
369:   options->showSolution        = PETSC_FALSE;
370:   options->restart             = PETSC_FALSE;
371:   options->check               = PETSC_FALSE;
372:   options->viewHierarchy       = PETSC_FALSE;
373:   options->simplex             = PETSC_TRUE;
374:   options->quiet               = PETSC_FALSE;
375:   options->nonzInit            = PETSC_FALSE;

377:   PetscOptionsBegin(comm, "", "Poisson Problem Options", "DMPLEX");
378:   PetscOptionsInt("-debug", "The debugging level", "ex12.c", options->debug, &options->debug, NULL);
379:   run  = options->runType;
380:   PetscOptionsEList("-run_type", "The run type", "ex12.c", runTypes, 3, runTypes[options->runType], &run, NULL);

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

384:   PetscOptionsInt("-dim", "The topological mesh dimension", "ex12.c", options->dim, &options->dim, NULL);
385:   bd = options->periodicity[0];
386:   PetscOptionsEList("-x_periodicity", "The x-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[0]], &bd, NULL);
387:   options->periodicity[0] = (DMBoundaryType) bd;
388:   bd = options->periodicity[1];
389:   PetscOptionsEList("-y_periodicity", "The y-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[1]], &bd, NULL);
390:   options->periodicity[1] = (DMBoundaryType) bd;
391:   bd = options->periodicity[2];
392:   PetscOptionsEList("-z_periodicity", "The z-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[2]], &bd, NULL);
393:   options->periodicity[2] = (DMBoundaryType) bd;
394:   n = 3;
395:   PetscOptionsIntArray("-cells", "The initial mesh division", "ex12.c", options->cells, &n, NULL);
396:   PetscOptionsString("-f", "Mesh filename to read", "ex12.c", options->filename, options->filename, sizeof(options->filename), &flg);
397:   PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex12.c", options->interpolate, &options->interpolate, NULL);
398:   PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex12.c", options->refinementLimit, &options->refinementLimit, NULL);
399:   bc   = options->bcType;
400:   PetscOptionsEList("-bc_type","Type of boundary condition","ex12.c",bcTypes,3,bcTypes[options->bcType],&bc,NULL);
401:   options->bcType = (BCType) bc;
402:   coeff = options->variableCoefficient;
403:   PetscOptionsEList("-variable_coefficient","Type of variable coefficent","ex12.c",coeffTypes,4,coeffTypes[options->variableCoefficient],&coeff,NULL);
404:   options->variableCoefficient = (CoeffType) coeff;

406:   PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL);
407:   PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL);
408:   PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL);
409:   PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL);
410:   PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL);
411:   PetscOptionsBool("-check", "Compare with default integration routines", "ex12.c", options->check, &options->check, NULL);
412:   PetscOptionsBool("-dm_view_hierarchy", "View the coarsened hierarchy", "ex12.c", options->viewHierarchy, &options->viewHierarchy, NULL);
413:   PetscOptionsBool("-simplex", "Simplicial (true) or tensor (false) mesh", "ex12.c", options->simplex, &options->simplex, NULL);
414:   PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL);
415:   PetscOptionsBool("-nonzero_initial_guess", "nonzero intial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL);
416:   PetscOptionsEnd();
417:   PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
418:   return(0);
419: }

421: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
422: {
423:   DMLabel        label;

427:   DMCreateLabel(dm, name);
428:   DMGetLabel(dm, name, &label);
429:   DMPlexMarkBoundaryFaces(dm, label);
430:   DMPlexLabelComplete(dm, label);
431:   return(0);
432: }

434: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
435: {
436:   PetscInt       dim             = user->dim;
437:   const char    *filename        = user->filename;
438:   PetscBool      interpolate     = user->interpolate;
439:   PetscReal      refinementLimit = user->refinementLimit;
440:   size_t         len;

444:   PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
445:   PetscStrlen(filename, &len);
446:   if (!len) {
447:     PetscInt d;

449:     if (user->periodicity[0] || user->periodicity[1] || user->periodicity[2]) for (d = 0; d < dim; ++d) user->cells[d] = PetscMax(user->cells[d], 3);
450:     DMPlexCreateBoxMesh(comm, dim, user->simplex, user->cells, NULL, NULL, user->periodicity, interpolate, dm);
451:     PetscObjectSetName((PetscObject) *dm, "Mesh");
452:   } else {
453:     DMPlexCreateFromFile(comm, filename, interpolate, dm);
454:     DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
455:   }
456:   {
457:     PetscPartitioner part;
458:     DM               refinedMesh     = NULL;
459:     DM               distributedMesh = NULL;

461:     /* Refine mesh using a volume constraint */
462:     if (refinementLimit > 0.0) {
463:       DMPlexSetRefinementLimit(*dm, refinementLimit);
464:       DMRefine(*dm, comm, &refinedMesh);
465:       if (refinedMesh) {
466:         const char *name;

468:         PetscObjectGetName((PetscObject) *dm,         &name);
469:         PetscObjectSetName((PetscObject) refinedMesh,  name);
470:         DMDestroy(dm);
471:         *dm  = refinedMesh;
472:       }
473:     }
474:     /* Distribute mesh over processes */
475:     DMPlexGetPartitioner(*dm,&part);
476:     PetscPartitionerSetFromOptions(part);
477:     DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
478:     if (distributedMesh) {
479:       DMDestroy(dm);
480:       *dm  = distributedMesh;
481:     }
482:   }
483:   if (user->bcType == NEUMANN) {
484:     DMLabel   label;

486:     DMCreateLabel(*dm, "boundary");
487:     DMGetLabel(*dm, "boundary", &label);
488:     DMPlexMarkBoundaryFaces(*dm, label);
489:   } else if (user->bcType == DIRICHLET) {
490:     PetscBool hasLabel;

492:     DMHasLabel(*dm,"marker",&hasLabel);
493:     if (!hasLabel) {CreateBCLabel(*dm, "marker");}
494:   }
495:   {
496:     char      convType[256];
497:     PetscBool flg;

499:     PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
500:     PetscOptionsFList("-dm_plex_convert_type","Convert DMPlex to another format","ex12",DMList,DMPLEX,convType,256,&flg);
501:     PetscOptionsEnd();
502:     if (flg) {
503:       DM dmConv;

505:       DMConvert(*dm,convType,&dmConv);
506:       if (dmConv) {
507:         DMDestroy(dm);
508:         *dm  = dmConv;
509:       }
510:     }
511:   }
512:   DMLocalizeCoordinates(*dm); /* needed for periodic */
513:   DMSetFromOptions(*dm);
514:   DMViewFromOptions(*dm, NULL, "-dm_view");
515:   if (user->viewHierarchy) {
516:     DM       cdm = *dm;
517:     PetscInt i   = 0;
518:     char     buf[256];

520:     while (cdm) {
521:       DMSetUp(cdm);
522:       DMGetCoarseDM(cdm, &cdm);
523:       ++i;
524:     }
525:     cdm = *dm;
526:     while (cdm) {
527:       PetscViewer       viewer;
528:       PetscBool   isHDF5, isVTK;

530:       --i;
531:       PetscViewerCreate(comm,&viewer);
532:       PetscViewerSetType(viewer,PETSCVIEWERHDF5);
533:       PetscViewerSetOptionsPrefix(viewer,"hierarchy_");
534:       PetscViewerSetFromOptions(viewer);
535:       PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&isHDF5);
536:       PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERVTK,&isVTK);
537:       if (isHDF5) {
538:         PetscSNPrintf(buf, 256, "ex12-%d.h5", i);
539:       }
540:       else if (isVTK) {
541:         PetscSNPrintf(buf, 256, "ex12-%d.vtu", i);
542:         PetscViewerPushFormat(viewer,PETSC_VIEWER_VTK_VTU);
543:       }
544:       else {
545:         PetscSNPrintf(buf, 256, "ex12-%d", i);
546:       }
547:       PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
548:       PetscViewerFileSetName(viewer,buf);
549:       DMView(cdm, viewer);
550:       PetscViewerDestroy(&viewer);
551:       DMGetCoarseDM(cdm, &cdm);
552:     }
553:   }
554:   PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
555:   return(0);
556: }

558: static PetscErrorCode SetupProblem(PetscDS prob, AppCtx *user)
559: {
560:   const PetscInt id = 1;

564:   switch (user->variableCoefficient) {
565:   case COEFF_NONE:
566:     if (user->periodicity[0]) {
567:       if (user->periodicity[1]) {
568:         PetscDSSetResidual(prob, 0, f0_xytrig_u, f1_u);
569:         PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
570:       } else {
571:         PetscDSSetResidual(prob, 0, f0_xtrig_u,  f1_u);
572:         PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
573:       }
574:     } else {
575:       PetscDSSetResidual(prob, 0, f0_u, f1_u);
576:       PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
577:     }
578:     break;
579:   case COEFF_ANALYTIC:
580:     PetscDSSetResidual(prob, 0, f0_analytic_u, f1_analytic_u);
581:     PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_analytic_uu);
582:     break;
583:   case COEFF_FIELD:
584:     PetscDSSetResidual(prob, 0, f0_analytic_u, f1_field_u);
585:     PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_field_uu);
586:     break;
587:   case COEFF_NONLINEAR:
588:     PetscDSSetResidual(prob, 0, f0_analytic_nonlinear_u, f1_analytic_nonlinear_u);
589:     PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_analytic_nonlinear_uu);
590:     break;
591:   default: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid variable coefficient type %d", user->variableCoefficient);
592:   }
593:   switch (user->dim) {
594:   case 2:
595:     if (user->periodicity[0]) {
596:       if (user->periodicity[1]) {
597:         user->exactFuncs[0] = xytrig_u_2d;
598:       } else {
599:         user->exactFuncs[0] = xtrig_u_2d;
600:       }
601:     } else {
602:       user->exactFuncs[0]  = quadratic_u_2d;
603:       user->exactFields[0] = quadratic_u_field_2d;
604:     }
605:     if (user->bcType == NEUMANN) {PetscDSSetBdResidual(prob, 0, f0_bd_u, f1_bd_zero);}
606:     break;
607:   case 3:
608:     user->exactFuncs[0]  = quadratic_u_3d;
609:     user->exactFields[0] = quadratic_u_field_3d;
610:     if (user->bcType == NEUMANN) {PetscDSSetBdResidual(prob, 0, f0_bd_u, f1_bd_zero);}
611:     break;
612:   default:
613:     SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %d", user->dim);
614:   }
615:   PetscDSAddBoundary(prob, user->bcType == DIRICHLET ? (user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL) : DM_BC_NATURAL,
616:                             "wall", user->bcType == DIRICHLET ? "marker" : "boundary", 0, 0, NULL,
617:                             user->fieldBC ? (void (*)()) user->exactFields[0] : (void (*)()) user->exactFuncs[0], 1, &id, user);
618:   return(0);
619: }

621: static PetscErrorCode SetupMaterial(DM dm, DM dmAux, AppCtx *user)
622: {
623:   PetscErrorCode (*matFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {nu_2d};
624:   Vec            nu;

628:   DMCreateLocalVector(dmAux, &nu);
629:   DMProjectFunctionLocal(dmAux, 0.0, matFuncs, NULL, INSERT_ALL_VALUES, nu);
630:   PetscObjectCompose((PetscObject) dm, "A", (PetscObject) nu);
631:   VecDestroy(&nu);
632:   return(0);
633: }

635: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
636: {
637:   PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
638:   Vec            uexact;
639:   PetscInt       dim;

643:   DMGetDimension(dm, &dim);
644:   if (dim == 2) bcFuncs[0] = quadratic_u_2d;
645:   else          bcFuncs[0] = quadratic_u_3d;
646:   DMCreateLocalVector(dmAux, &uexact);
647:   DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
648:   PetscObjectCompose((PetscObject) dm, "A", (PetscObject) uexact);
649:   VecDestroy(&uexact);
650:   return(0);
651: }

653: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
654: {
655:   DM             cdm   = dm;
656:   const PetscInt dim   = user->dim;
657:   PetscFE        feAux = NULL;
658:   PetscFE        feCh  = NULL;
659:   PetscFE        fe;
660:   PetscDS        prob, probAux = NULL, probCh = NULL;
661:   PetscBool      simplex = user->simplex;

665:   /* Create finite element */
666:   PetscFECreateDefault(dm, dim, 1, simplex, NULL, -1, &fe);
667:   PetscObjectSetName((PetscObject) fe, "potential");
668:   if (user->variableCoefficient == COEFF_FIELD) {
669:     PetscQuadrature q;

671:     PetscFECreateDefault(dm, dim, 1, simplex, "mat_", -1, &feAux);
672:     PetscFEGetQuadrature(fe, &q);
673:     PetscFESetQuadrature(feAux, q);
674:     PetscDSCreate(PetscObjectComm((PetscObject)dm),&probAux);
675:     PetscDSSetDiscretization(probAux, 0, (PetscObject) feAux);
676:   } else if (user->fieldBC) {
677:     PetscQuadrature q;

679:     PetscFECreateDefault(dm, dim, 1, simplex, "bc_", -1, &feAux);
680:     PetscFEGetQuadrature(fe, &q);
681:     PetscFESetQuadrature(feAux, q);
682:     PetscDSCreate(PetscObjectComm((PetscObject)dm),&probAux);
683:     PetscDSSetDiscretization(probAux, 0, (PetscObject) feAux);
684:   }
685:   if (user->check) {
686:     PetscFECreateDefault(dm, dim, 1, simplex, "ch_", -1, &feCh);
687:     PetscDSCreate(PetscObjectComm((PetscObject)dm),&probCh);
688:     PetscDSSetDiscretization(probCh, 0, (PetscObject) feCh);
689:   }
690:   /* Set discretization and boundary conditions for each mesh */
691:   DMGetDS(dm, &prob);
692:   PetscDSSetDiscretization(prob, 0, (PetscObject) fe);
693:   SetupProblem(prob, user);
694:   while (cdm) {
695:     DM coordDM;

697:     DMSetDS(cdm,prob);
698:     DMGetCoordinateDM(cdm,&coordDM);
699:     if (feAux) {
700:       DM      dmAux;

702:       DMClone(cdm, &dmAux);
703:       DMSetCoordinateDM(dmAux, coordDM);
704:       DMSetDS(dmAux, probAux);
705:       PetscObjectCompose((PetscObject) dm, "dmAux", (PetscObject) dmAux);
706:       if (user->fieldBC) {SetupBC(cdm, dmAux, user);}
707:       else               {SetupMaterial(cdm, dmAux, user);}
708:       DMDestroy(&dmAux);
709:     }
710:     if (feCh) {
711:       DM      dmCh;

713:       DMClone(cdm, &dmCh);
714:       DMSetCoordinateDM(dmCh, coordDM);
715:       DMSetDS(dmCh, probCh);
716:       PetscObjectCompose((PetscObject) dm, "dmCh", (PetscObject) dmCh);
717:       DMDestroy(&dmCh);
718:     }
719:     if (user->bcType == DIRICHLET) {
720:       PetscBool hasLabel;

722:       DMHasLabel(cdm, "marker", &hasLabel);
723:       if (!hasLabel) {CreateBCLabel(cdm, "marker");}
724:     }
725:     DMGetCoarseDM(cdm, &cdm);
726:   }
727:   PetscFEDestroy(&fe);
728:   PetscFEDestroy(&feAux);
729:   PetscFEDestroy(&feCh);
730:   PetscDSDestroy(&probAux);
731:   PetscDSDestroy(&probCh);
732:   return(0);
733: }

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

737: /*@C
738:   KSPMonitorError - Outputs the error at each iteration of an iterative solver.

740:   Collective on KSP

742:   Input Parameters:
743: + ksp   - the KSP
744: . its   - iteration number
745: . rnorm - 2-norm, preconditioned residual value (may be estimated).
746: - ctx   - monitor context

748:   Level: intermediate

750: .keywords: KSP, default, monitor, residual
751: .seealso: KSPMonitorSet(), KSPMonitorTrueResidualNorm(), KSPMonitorDefault()
752: @*/
753: static PetscErrorCode KSPMonitorError(KSP ksp, PetscInt its, PetscReal rnorm, void *ctx)
754: {
755:   AppCtx        *user = (AppCtx *) ctx;
756:   DM             dm;
757:   Vec            du, r;
758:   PetscInt       level = 0;
759:   PetscBool      hasLevel;
760: #if defined(PETSC_HAVE_HDF5)
761:   PetscViewer    viewer;
762: #endif
763:   char           buf[256];

767:   KSPGetDM(ksp, &dm);
768:   /* Calculate solution */
769:   {
770:     PC        pc = user->pcmg; /* The MG PC */
771:     DM        fdm = NULL,  cdm;
772:     KSP       fksp, cksp;
773:     Vec       fu,   cu;
774:     PetscInt  levels, l;

776:     KSPBuildSolution(ksp, NULL, &du);
777:     PetscObjectComposedDataGetInt((PetscObject) ksp, PetscMGLevelId, level, hasLevel);
778:     PCMGGetLevels(pc, &levels);
779:     PCMGGetSmoother(pc, levels-1, &fksp);
780:     KSPBuildSolution(fksp, NULL, &fu);
781:     for (l = levels-1; l > level; --l) {
782:       Mat R;
783:       Vec s;

785:       PCMGGetSmoother(pc, l-1, &cksp);
786:       KSPGetDM(cksp, &cdm);
787:       DMGetGlobalVector(cdm, &cu);
788:       PCMGGetRestriction(pc, l, &R);
789:       PCMGGetRScale(pc, l, &s);
790:       MatRestrict(R, fu, cu);
791:       VecPointwiseMult(cu, cu, s);
792:       if (l < levels-1) {DMRestoreGlobalVector(fdm, &fu);}
793:       fdm  = cdm;
794:       fu   = cu;
795:     }
796:     if (levels-1 > level) {
797:       VecAXPY(du, 1.0, cu);
798:       DMRestoreGlobalVector(cdm, &cu);
799:     }
800:   }
801:   /* Calculate error */
802:   DMGetGlobalVector(dm, &r);
803:   DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
804:   VecAXPY(r,-1.0,du);
805:   PetscObjectSetName((PetscObject) r, "solution error");
806:   /* View error */
807:   PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
808: #if defined(PETSC_HAVE_HDF5)
809:   PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
810:   VecView(r, viewer);
811:   PetscViewerDestroy(&viewer);
812: #else
813:   SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"You need to configure with --download-hdf5");
814: #endif
815:   /* Cleanup */
816:   DMRestoreGlobalVector(dm, &r);
817:   return(0);
818: }

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

823:   Collective on SNES

825:   Input Parameters:
826: + snes  - the SNES
827: . its   - iteration number
828: . rnorm - 2-norm of residual
829: - ctx   - user context

831:   Level: intermediate

833: .keywords: SNES, nonlinear, default, monitor, norm
834: .seealso: SNESMonitorDefault(), SNESMonitorSet(), SNESMonitorSolution()
835: @*/
836: static PetscErrorCode SNESMonitorError(SNES snes, PetscInt its, PetscReal rnorm, void *ctx)
837: {
838:   AppCtx        *user = (AppCtx *) ctx;
839:   DM             dm;
840:   Vec            u, r;
841:   PetscInt       level = -1;
842:   PetscBool      hasLevel;
843: #if defined(PETSC_HAVE_HDF5)
844:   PetscViewer    viewer;
845: #endif
846:   char           buf[256];

850:   SNESGetDM(snes, &dm);
851:   /* Calculate error */
852:   SNESGetSolution(snes, &u);
853:   DMGetGlobalVector(dm, &r);
854:   PetscObjectSetName((PetscObject) r, "solution error");
855:   DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
856:   VecAXPY(r, -1.0, u);
857:   /* View error */
858:   PetscObjectComposedDataGetInt((PetscObject) snes, PetscMGLevelId, level, hasLevel);
859:   PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
860: #if defined(PETSC_HAVE_HDF5)
861:   PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
862:   VecView(r, viewer);
863:   PetscViewerDestroy(&viewer);
864: #else
865:   SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"You need to configure with --download-hdf5");
866: #endif
867:   /* Cleanup */
868:   DMRestoreGlobalVector(dm, &r);
869:   return(0);
870: }

872: int main(int argc, char **argv)
873: {
874:   DM             dm;          /* Problem specification */
875:   SNES           snes;        /* nonlinear solver */
876:   Vec            u,r;         /* solution, residual vectors */
877:   Mat            A,J;         /* Jacobian matrix */
878:   MatNullSpace   nullSpace;   /* May be necessary for Neumann conditions */
879:   AppCtx         user;        /* user-defined work context */
880:   JacActionCtx   userJ;       /* context for Jacobian MF action */
881:   PetscInt       its;         /* iterations for convergence */
882:   PetscReal      error = 0.0; /* L_2 error in the solution */
883:   PetscBool      isFAS;

886:   PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
887:   ProcessOptions(PETSC_COMM_WORLD, &user);
888:   SNESCreate(PETSC_COMM_WORLD, &snes);
889:   CreateMesh(PETSC_COMM_WORLD, &user, &dm);
890:   SNESSetDM(snes, dm);
891:   DMSetApplicationContext(dm, &user);

893:   PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
894:   SetupDiscretization(dm, &user);

896:   DMCreateGlobalVector(dm, &u);
897:   PetscObjectSetName((PetscObject) u, "potential");
898:   VecDuplicate(u, &r);

900:   DMCreateMatrix(dm, &J);
901:   if (user.jacobianMF) {
902:     PetscInt M, m, N, n;

904:     MatGetSize(J, &M, &N);
905:     MatGetLocalSize(J, &m, &n);
906:     MatCreate(PETSC_COMM_WORLD, &A);
907:     MatSetSizes(A, m, n, M, N);
908:     MatSetType(A, MATSHELL);
909:     MatSetUp(A);
910: #if 0
911:     MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction);
912: #endif

914:     userJ.dm   = dm;
915:     userJ.J    = J;
916:     userJ.user = &user;

918:     DMCreateLocalVector(dm, &userJ.u);
919:     if (user.fieldBC) {DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);}
920:     else              {DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);}
921:     MatShellSetContext(A, &userJ);
922:   } else {
923:     A = J;
924:   }
925:   if (user.bcType == NEUMANN) {
926:     MatNullSpaceCreate(PetscObjectComm((PetscObject) dm), PETSC_TRUE, 0, NULL, &nullSpace);
927:     MatSetNullSpace(A, nullSpace);
928:   }

930:   DMPlexSetSNESLocalFEM(dm,&user,&user,&user);
931:   SNESSetJacobian(snes, A, J, NULL, NULL);

933:   SNESSetFromOptions(snes);

935:   if (user.fieldBC) {DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u);}
936:   else              {DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);}
937:   if (user.restart) {
938: #if defined(PETSC_HAVE_HDF5)
939:     PetscViewer viewer;

941:     PetscViewerCreate(PETSC_COMM_WORLD, &viewer);
942:     PetscViewerSetType(viewer, PETSCVIEWERHDF5);
943:     PetscViewerFileSetMode(viewer, FILE_MODE_READ);
944:     PetscViewerFileSetName(viewer, user.filename);
945:     PetscViewerHDF5PushGroup(viewer, "/fields");
946:     VecLoad(u, viewer);
947:     PetscViewerHDF5PopGroup(viewer);
948:     PetscViewerDestroy(&viewer);
949: #endif
950:   }
951:   if (user.showInitial) {
952:     Vec lv;
953:     DMGetLocalVector(dm, &lv);
954:     DMGlobalToLocalBegin(dm, u, INSERT_VALUES, lv);
955:     DMGlobalToLocalEnd(dm, u, INSERT_VALUES, lv);
956:     DMPrintLocalVec(dm, "Local function", 1.0e-10, lv);
957:     DMRestoreLocalVector(dm, &lv);
958:   }
959:   if (user.viewHierarchy) {
960:     SNES      lsnes;
961:     KSP       ksp;
962:     PC        pc;
963:     PetscInt  numLevels, l;
964:     PetscBool isMG;

966:     PetscObjectTypeCompare((PetscObject) snes, SNESFAS, &isFAS);
967:     if (isFAS) {
968:       SNESFASGetLevels(snes, &numLevels);
969:       for (l = 0; l < numLevels; ++l) {
970:         SNESFASGetCycleSNES(snes, l, &lsnes);
971:         SNESMonitorSet(lsnes, SNESMonitorError, &user, NULL);
972:       }
973:     } else {
974:       SNESGetKSP(snes, &ksp);
975:       KSPGetPC(ksp, &pc);
976:       PetscObjectTypeCompare((PetscObject) pc, PCMG, &isMG);
977:       if (isMG) {
978:         user.pcmg = pc;
979:         PCMGGetLevels(pc, &numLevels);
980:         for (l = 0; l < numLevels; ++l) {
981:           PCMGGetSmootherDown(pc, l, &ksp);
982:           KSPMonitorSet(ksp, KSPMonitorError, &user, NULL);
983:         }
984:       }
985:     }
986:   }
987:   if (user.runType == RUN_FULL || user.runType == RUN_EXACT) {
988:     PetscErrorCode (*initialGuess[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {zero};

990:     if (user.nonzInit) initialGuess[0] = ecks;
991:     if (user.runType == RUN_FULL) {
992:       DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
993:     }
994:     if (user.debug) {
995:       PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
996:       VecView(u, PETSC_VIEWER_STDOUT_WORLD);
997:     }
998:     SNESSolve(snes, NULL, u);
999:     SNESGetIterationNumber(snes, &its);
1000:     PetscPrintf(PETSC_COMM_WORLD, "Number of SNES iterations = %D\n", its);
1001:     DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
1002:     if (error < 1.0e-11) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < 1.0e-11\n");}
1003:     else                 {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", error);}
1004:     if (user.showSolution) {
1005:       PetscPrintf(PETSC_COMM_WORLD, "Solution\n");
1006:       VecChop(u, 3.0e-9);
1007:       VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1008:     }
1009:   } else if (user.runType == RUN_PERF) {
1010:     PetscReal res = 0.0;

1012:     SNESComputeFunction(snes, u, r);
1013:     PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1014:     VecChop(r, 1.0e-10);
1015:     VecNorm(r, NORM_2, &res);
1016:     PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1017:   } else {
1018:     PetscReal res = 0.0;

1020:     /* Check discretization error */
1021:     PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1022:     if (!user.quiet) {VecView(u, PETSC_VIEWER_STDOUT_WORLD);}
1023:     DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
1024:     if (error < 1.0e-11) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < 1.0e-11\n");}
1025:     else                 {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", error);}
1026:     /* Check residual */
1027:     SNESComputeFunction(snes, u, r);
1028:     PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1029:     VecChop(r, 1.0e-10);
1030:     if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1031:     VecNorm(r, NORM_2, &res);
1032:     PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1033:     /* Check Jacobian */
1034:     {
1035:       Vec          b;

1037:       SNESComputeJacobian(snes, u, A, A);
1038:       VecDuplicate(u, &b);
1039:       VecSet(r, 0.0);
1040:       SNESComputeFunction(snes, r, b);
1041:       MatMult(A, u, r);
1042:       VecAXPY(r, 1.0, b);
1043:       VecDestroy(&b);
1044:       PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n");
1045:       VecChop(r, 1.0e-10);
1046:       if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1047:       VecNorm(r, NORM_2, &res);
1048:       PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", res);
1049:     }
1050:   }
1051:   VecViewFromOptions(u, NULL, "-vec_view");

1053:   if (user.bcType == NEUMANN) {MatNullSpaceDestroy(&nullSpace);}
1054:   if (user.jacobianMF) {VecDestroy(&userJ.u);}
1055:   if (A != J) {MatDestroy(&A);}
1056:   MatDestroy(&J);
1057:   VecDestroy(&u);
1058:   VecDestroy(&r);
1059:   SNESDestroy(&snes);
1060:   DMDestroy(&dm);
1061:   PetscFree2(user.exactFuncs, user.exactFields);
1062:   PetscFinalize();
1063:   return ierr;
1064: }

1066: /*TEST
1067:   build:
1068:     requires: hdf5 triangle
1069:   # 2D serial P1 test 0-4
1070:   test:
1071:     suffix: 0
1072:     requires: hdf5 triangle
1073:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 0 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1074:   test:
1075:     suffix: 1
1076:     requires: hdf5 triangle
1077:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1078:   test:
1079:     suffix: 2
1080:     requires: hdf5 triangle
1081:     args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1082:   test:
1083:     suffix: 3
1084:     requires: hdf5 triangle
1085:     args: -run_type test -refinement_limit 0.0    -bc_type neumann   -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
1086:   test:
1087:     suffix: 4
1088:     requires: hdf5 triangle
1089:     args: -run_type test -refinement_limit 0.0625 -bc_type neumann   -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1090:   # 2D serial P2 test 5-8
1091:   test:
1092:     suffix: 5
1093:     requires: hdf5 triangle
1094:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1095:   test:
1096:     suffix: 6
1097:     requires: hdf5 triangle
1098:     args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1099:   test:
1100:     suffix: 7
1101:     requires: hdf5 triangle
1102:     args: -run_type test -refinement_limit 0.0    -bc_type neumann   -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
1103:   test:
1104:     suffix: 8
1105:     requires: hdf5 triangle
1106:     args: -run_type test -refinement_limit 0.0625 -bc_type neumann   -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
1107:   # 3D serial P1 test 9-12
1108:   test:
1109:     suffix: 9
1110:     requires: hdf5 ctetgen
1111:     args: -run_type test -dim 3 -refinement_limit 0.0    -bc_type dirichlet -interpolate 0 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1
1112:   test:
1113:     suffix: 10
1114:     requires: hdf5 ctetgen
1115:     args: -run_type test -dim 3 -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1
1116:   test:
1117:     suffix: 11
1118:     requires: hdf5 ctetgen
1119:     args: -run_type test -dim 3 -refinement_limit 0.0125 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1
1120:   test:
1121:     suffix: 12
1122:     requires: hdf5 ctetgen
1123:     args: -run_type test -dim 3 -refinement_limit 0.0    -bc_type neumann   -interpolate 1 -petscspace_order 1 -snes_fd -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1
1124:   # Analytic variable coefficient 13-20
1125:   test:
1126:     suffix: 13
1127:     requires: hdf5 triangle
1128:     args: -run_type test -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1129:   test:
1130:     suffix: 14
1131:     requires: hdf5 triangle
1132:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1133:   test:
1134:     suffix: 15
1135:     requires: hdf5 triangle
1136:     args: -run_type test -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1137:   test:
1138:     suffix: 16
1139:     requires: hdf5 triangle
1140:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1141:   test:
1142:     suffix: 17
1143:     requires: hdf5 ctetgen
1144:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1145:   test:
1146:     suffix: 18
1147:     requires: hdf5 ctetgen
1148:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1149:   test:
1150:     suffix: 19
1151:     requires: hdf5 ctetgen
1152:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1153:   test:
1154:     suffix: 20
1155:     requires: hdf5 ctetgen
1156:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1157:   # P1 variable coefficient 21-28
1158:   test:
1159:     suffix: 21
1160:     requires: hdf5 triangle
1161:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1162:   test:
1163:     suffix: 22
1164:     requires: hdf5 triangle
1165:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1166:   test:
1167:     suffix: 23
1168:     requires: hdf5 triangle
1169:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1170:   test:
1171:     suffix: 24
1172:     requires: hdf5 triangle
1173:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1174:   test:
1175:     suffix: 25
1176:     requires: hdf5 ctetgen
1177:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1178:   test:
1179:     suffix: 26
1180:     requires: hdf5 ctetgen
1181:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1182:   test:
1183:     suffix: 27
1184:     requires: hdf5 ctetgen
1185:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1186:   test:
1187:     suffix: 28
1188:     requires: hdf5 ctetgen
1189:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1190:   # P0 variable coefficient 29-36
1191:   test:
1192:     suffix: 29
1193:     requires: hdf5 triangle
1194:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1195:   test:
1196:     suffix: 30
1197:     requires: hdf5 triangle
1198:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1199:   test:
1200:     suffix: 31
1201:     requires: hdf5 triangle
1202:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1203:   test:
1204:     requires: hdf5 triangle
1205:     suffix: 32
1206:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1207:   test:
1208:     requires: hdf5 ctetgen
1209:     suffix: 33
1210:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1211:   test:
1212:     suffix: 34
1213:     requires: hdf5 ctetgen
1214:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1215:   test:
1216:     suffix: 35
1217:     requires: hdf5 ctetgen
1218:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1219:   test:
1220:     suffix: 36
1221:     requires: hdf5 ctetgen
1222:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1223:   # Full solve 39-44
1224:   test:
1225:     suffix: 39
1226:     requires: hdf5 triangle
1227:     args: -run_type full -refinement_limit 0.015625 -interpolate 1 -petscspace_order 2 -pc_type gamg -ksp_rtol 1.0e-10 -ksp_monitor_short -ksp_converged_reason -snes_monitor_short -snes_converged_reason
1228:   test:
1229:     suffix: 40
1230:     requires: hdf5 triangle
1231:     args: -run_type full -refinement_limit 0.015625 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 2 -pc_type svd -ksp_rtol 1.0e-10 -snes_monitor_short -snes_converged_reason
1232:   test:
1233:     suffix: 41
1234:     requires: hdf5 triangle
1235:     args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 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 -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1236:   test:
1237:     suffix: 42
1238:     requires: hdf5 triangle
1239:     args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 3 -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 -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1240:   test:
1241:     suffix: 43
1242:     requires: hdf5 triangle
1243:     nsize: 2
1244:     args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 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 -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1245:   test:
1246:     suffix: 44
1247:     requires: hdf5 triangle
1248:     nsize: 2
1249:     args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 3 -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 -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1250:   # Restarting
1251:   test:
1252:     suffix: restart
1253:     requires: hdf5 triangle !complex
1254:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_order 1
1255:     test:
1256:       args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1257:     test:
1258:       args: -f sol.h5 -restart
1259:   # Periodicity
1260:   test:
1261:     suffix: periodic_0
1262:     requires: hdf5 triangle
1263:     args: -run_type full -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_order 1
1264:   # 2D serial P1 test with field bc
1265:   test:
1266:     suffix: field_bc_p1_0
1267:     requires: hdf5 triangle
1268:     args: -run_type test              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1269:   test:
1270:     suffix: field_bc_p1_1
1271:     requires: hdf5 triangle
1272:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1273:   test:
1274:     suffix: field_bc_p1_2
1275:     requires: hdf5 triangle
1276:     args: -run_type test              -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1277:   test:
1278:     suffix: field_bc_p1_3
1279:     requires: hdf5 triangle
1280:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1281:   # 3D serial P1 test with field bc
1282:   test:
1283:     suffix: field_bc_p1_4
1284:     requires: hdf5 ctetgen
1285:     args: -run_type test -dim 3              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1286:   test:
1287:     suffix: field_bc_p1_5
1288:     requires: hdf5 ctetgen
1289:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1290:   test:
1291:     suffix: field_bc_p1_6
1292:     requires: hdf5 ctetgen
1293:     args: -run_type test -dim 3              -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1294:   test:
1295:     suffix: field_bc_p1_7
1296:     requires: hdf5 ctetgen
1297:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1298:   # 2D serial P2 test with field bc
1299:   test:
1300:     suffix: field_bc_p2_0
1301:     requires: hdf5 triangle
1302:     args: -run_type test              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1303:   test:
1304:     suffix: field_bc_p2_1
1305:     requires: hdf5 triangle
1306:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1307:   test:
1308:     suffix: field_bc_p2_2
1309:     requires: hdf5 triangle
1310:     args: -run_type test              -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1311:   test:
1312:     suffix: field_bc_p2_3
1313:     requires: hdf5 triangle
1314:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1315:   # 3D serial P2 test with field bc
1316:   test:
1317:     suffix: field_bc_p2_4
1318:     requires: hdf5 ctetgen
1319:     args: -run_type test -dim 3              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1320:   test:
1321:     suffix: field_bc_p2_5
1322:     requires: hdf5 ctetgen
1323:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1324:   test:
1325:     suffix: field_bc_p2_6
1326:     requires: hdf5 ctetgen
1327:     args: -run_type test -dim 3              -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1328:   test:
1329:     suffix: field_bc_p2_7
1330:     requires: hdf5 ctetgen
1331:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1332:   # Full solve simplex: Convergence
1333:   test:
1334:     suffix: tet_conv_p1_r0
1335:     requires: hdf5 ctetgen
1336:     args: -run_type full -dim 3 -dm_refine 0 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -dm_view -snes_converged_reason -pc_type lu -cells 1,1,1
1337:   test:
1338:     suffix: tet_conv_p1_r2
1339:     requires: hdf5 ctetgen
1340:     args: -run_type full -dim 3 -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -dm_view -snes_converged_reason -pc_type lu -cells 1,1,1
1341:   test:
1342:     suffix: tet_conv_p1_r3
1343:     requires: hdf5 ctetgen
1344:     args: -run_type full -dim 3 -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -dm_view -snes_converged_reason -pc_type lu -cells 1,1,1
1345:   test:
1346:     suffix: tet_conv_p2_r0
1347:     requires: hdf5 ctetgen
1348:     args: -run_type full -dim 3 -dm_refine 0 -bc_type dirichlet -interpolate 1 -petscspace_order 2 -dm_view -snes_converged_reason -pc_type lu -cells 1,1,1
1349:   test:
1350:     suffix: tet_conv_p2_r2
1351:     requires: hdf5 ctetgen
1352:     args: -run_type full -dim 3 -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_order 2 -dm_view -snes_converged_reason -pc_type lu -cells 1,1,1
1353:   # Full solve simplex: BDDC
1354:   test:
1355:     suffix: tri_bddc
1356:     requires: hdf5 triangle
1357:     nsize: 5
1358:     args: -run_type full -petscpartitioner_type simple -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 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 -ksp_converged_reason -snes_view -show_solution 0
1359:   # Full solve simplex: ASM
1360:   test:
1361:     suffix: tri_q2q1_asm_lu
1362:     requires: hdf5 triangle
1363:     args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 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 -ksp_converged_reason -snes_view -show_solution 0
1364:   test:
1365:     suffix: tri_q2q1_msm_lu
1366:     requires: hdf5 triangle
1367:     args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 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 -ksp_converged_reason -snes_view -show_solution 0
1368:   test:
1369:     suffix: tri_q2q1_asm_sor
1370:     requires: hdf5 triangle
1371:     args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 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 -ksp_converged_reason -snes_view -show_solution 0
1372:   test:
1373:     suffix: tri_q2q1_msm_sor
1374:     requires: hdf5 triangle
1375:     args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 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 -ksp_converged_reason -snes_view -show_solution 0
1376:   # Full solve simplex: FAS
1377:   test:
1378:     suffix: fas_newton_0
1379:     requires: hdf5 triangle
1380:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 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 -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1381:   test:
1382:     suffix: fas_newton_1
1383:     requires: hdf5 triangle
1384:     args: -run_type full -dm_refine_hierarchy 3 -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 3 -ksp_rtol 1.0e-10 -fas_coarse_pc_type lu -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1385:   test:
1386:     suffix: fas_ngs_0
1387:     requires: hdf5 triangle
1388:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 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 -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -dm_refine_hierarchy 1 -snes_view -fas_levels_1_snes_type ngs -fas_levels_1_snes_monitor_short
1389:   test:
1390:     suffix: fas_newton_coarse_0
1391:     requires: hdf5 pragmatic triangle
1392:     TODO: broken
1393:     args: -run_type full -dm_refine 2 -dm_plex_hash_location -variable_coefficient nonlinear -interpolate 1 -petscspace_order 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 -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1394:   test:
1395:     suffix: mg_newton_coarse_0
1396:     requires: hdf5 triangle pragmatic
1397:     args: -run_type full -dm_refine 3 -interpolate 1 -petscspace_order 1 -snes_monitor_short -ksp_monitor_true_residual -snes_linesearch_type basic -snes_converged_reason -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
1398:   test:
1399:     suffix: mg_newton_coarse_1
1400:     requires: hdf5 triangle pragmatic
1401:     args: -run_type full -dm_refine 5 -interpolate 1 -petscspace_order 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_monitor -ksp_monitor_true_residual -mg_levels_ksp_monitor_true_residual -dm_view -ksp_view
1402:   test:
1403:     suffix: mg_newton_coarse_2
1404:     requires: hdf5 triangle pragmatic
1405:     args: -run_type full -dm_refine 5 -interpolate 1 -petscspace_order 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_monitor -ksp_monitor_true_residual -mg_levels_ksp_monitor_true_residual -dm_view -ksp_view
1406:   # Full solve tensor
1407:   test:
1408:     suffix: tensor_plex_2d
1409:     requires: hdf5
1410:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 1 -dm_refine_hierarchy 2 -cells 2,2
1411:   test:
1412:     suffix: tensor_p4est_2d
1413:     requires: hdf5 p4est
1414:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 1 -dm_forest_initial_refinement 2 -dm_forest_minimum_refinement 0 -dm_plex_convert_type p4est -cells 2,2
1415:   test:
1416:     suffix: tensor_plex_3d
1417:     requires: hdf5
1418:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 1 -dim 3 -dm_refine_hierarchy 1 -cells 2,2,2
1419:   test:
1420:     suffix: tensor_p4est_3d
1421:     requires: hdf5 p4est
1422:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 1 -dm_forest_initial_refinement 1 -dm_forest_minimum_refinement 0 -dim 3 -dm_plex_convert_type p8est -cells 2,2,2
1423:   # Full solve tensor: AMR
1424:   test:
1425:     suffix: amr_0
1426:     requires: hdf5
1427:     nsize: 5
1428:     args: -run_type test -petscpartitioner_type simple -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 1 -dm_refine 1 -cells 2,2
1429:   test:
1430:     suffix: amr_1
1431:     requires: hdf5 p4est !complex
1432:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 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
1433:   test:
1434:     suffix: p4est_test_q2_conformal_serial
1435:     requires: hdf5 p4est
1436:     args: -run_type test -interpolate 1 -petscspace_order 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -cells 2,2
1437:   test:
1438:     suffix: p4est_test_q2_conformal_parallel
1439:     requires: hdf5 p4est
1440:     nsize: 7
1441:     args: -run_type test -interpolate 1 -petscspace_order 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type simple -cells 2,2
1442:   test:
1443:     suffix: p4est_test_q2_nonconformal_serial
1444:     requires: hdf5 p4est
1445:     filter: grep -v "CG or CGNE: variant"
1446:     args: -run_type test -interpolate 1 -petscspace_order 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
1447:   test:
1448:     suffix: p4est_test_q2_nonconformal_parallel
1449:     requires: hdf5 p4est
1450:     filter: grep -v "CG or CGNE: variant"
1451:     nsize: 7
1452:     args: -run_type test -interpolate 1 -petscspace_order 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
1453:   test:
1454:     suffix: p4est_exact_q2_conformal_serial
1455:     requires: hdf5 p4est
1456:     args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1457:   test:
1458:     suffix: p4est_exact_q2_conformal_parallel
1459:     requires: hdf5 p4est
1460:     nsize: 7
1461:     args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1462:   test:
1463:     suffix: p4est_exact_q2_nonconformal_serial
1464:     requires: hdf5 p4est
1465:     args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1466:   test:
1467:     suffix: p4est_exact_q2_nonconformal_parallel
1468:     requires: hdf5 p4est
1469:     nsize: 7
1470:     args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1471:   test:
1472:     suffix: p4est_full_q2_nonconformal_serial
1473:     requires: hdf5 p4est
1474:     filter: grep -v "CG or CGNE: variant"
1475:     args: -run_type full -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type cg -fas_coarse_pc_type jacobi -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1476:   test:
1477:     suffix: p4est_full_q2_nonconformal_parallel
1478:     requires: p4est hdf5
1479:     filter: grep -v "CG or CGNE: variant"
1480:     nsize: 7
1481:     args: -run_type full -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type cg -fas_coarse_pc_type jacobi -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1482:   test:
1483:     suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1484:     requires: hdf5 p4est
1485:     filter: grep -v "CG or CGNE: variant" | sed -e "s/BDDC: Graph max count: 9223372036854775807/BDDC: Graph max count: 2147483647/g"
1486:     nsize: 7
1487:     args: -run_type full -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -dm_mat_type is -pc_type bddc -ksp_type cg -fas_coarse_pc_type bddc -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1488:   test:
1489:     suffix: p4est_full_q2_nonconformal_parallel_bddc
1490:     requires: hdf5 p4est
1491:     filter: grep -v "CG or CGNE: variant" | sed -e "s/BDDC: Graph max count: 9223372036854775807/BDDC: Graph max count: 2147483647/g"
1492:     nsize: 7
1493:     args: -run_type full -interpolate 1 -petscspace_order 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 -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
1494:   test:
1495:     suffix: p4est_fas_q2_conformal_serial
1496:     requires: hdf5 p4est broken
1497:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 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 -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1498:   test:
1499:     suffix: p4est_fas_q2_nonconformal_serial
1500:     requires: hdf5 p4est broken
1501:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 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 -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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
1502:   test:
1503:     suffix: fas_newton_0_p4est
1504:     requires: hdf5 p4est
1505:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 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 -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -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

1507: TEST*/