Actual source code: ex12.c

petsc-master 2019-07-18
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 -info_exclude null,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} 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:   /* Solver */
 53:   PC             pcmg;              /* This is needed for error monitoring */
 54: } AppCtx;

 56: static PetscErrorCode zero(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 57: {
 58:   u[0] = 0.0;
 59:   return 0;
 60: }

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

 68: /*
 69:   In 2D for Dirichlet conditions, we use exact solution:

 71:     u = x^2 + y^2
 72:     f = 4

 74:   so that

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

 78:   For Neumann conditions, we have

 80:     -\nabla u \cdot -\hat y |_{y=0} =  (2y)|_{y=0} =  0 (bottom)
 81:     -\nabla u \cdot  \hat y |_{y=1} = -(2y)|_{y=1} = -2 (top)
 82:     -\nabla u \cdot -\hat x |_{x=0} =  (2x)|_{x=0} =  0 (left)
 83:     -\nabla u \cdot  \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)

 85:   Which we can express as

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

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

 91:     \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
 92: */
 93: static PetscErrorCode quadratic_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 94: {
 95:   *u = x[0]*x[0] + x[1]*x[1];
 96:   return 0;
 97: }

 99: static void quadratic_u_field_2d(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[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
103: {
104:   uexact[0] = a[0];
105: }

107: static PetscErrorCode circle_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
108: {
109:   const PetscReal alpha   = 500.;
110:   const PetscReal radius2 = PetscSqr(0.15);
111:   const PetscReal r2      = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5);
112:   const PetscReal xi      = alpha*(radius2 - r2);

114:   *u = PetscTanhScalar(xi) + 1.0;
115:   return 0;
116: }

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

123:   *u = PetscSinReal(alpha*xy) * (alpha*PetscAbsReal(xy) < 2*PETSC_PI ? 1 : 0.01);
124:   return 0;
125: }

127: static void f0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
128:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
129:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
130:                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
131: {
132:   f0[0] = 4.0;
133: }

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

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

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

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

159: static void f0_bd_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
160:                     const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
161:                     const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
162:                     PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
163: {
164:   PetscInt d;
165:   for (d = 0, f0[0] = 0.0; d < dim; ++d) f0[0] += -n[d]*2.0*x[d];
166: }

168: static void f1_bd_zero(PetscInt dim, PetscInt Nf, PetscInt NfAux,
169:                        const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
170:                        const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
171:                        PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
172: {
173:   PetscInt comp;
174:   for (comp = 0; comp < dim; ++comp) f1[comp] = 0.0;
175: }

177: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
178: static void f1_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 f1[])
182: {
183:   PetscInt d;
184:   for (d = 0; d < dim; ++d) f1[d] = u_x[d];
185: }

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

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

201:     u = sin(2 pi x)
202:     f = -4 pi^2 sin(2 pi x)

204:   so that

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

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

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

225:     u = sin(2 pi x) sin(2 pi y)
226:     f = -8 pi^2 sin(2 pi x)

228:   so that

230:     -\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
231: */
232: static PetscErrorCode xytrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
233: {
234:   *u = PetscSinReal(2.0*PETSC_PI*x[0])*PetscSinReal(2.0*PETSC_PI*x[1]);
235:   return 0;
236: }

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

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

249:     u  = x^2 + y^2
250:     f  = 6 (x + y)
251:     nu = (x + y)

253:   so that

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

263: void f0_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
264:                    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
265:                    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
266:                    PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
267: {
268:   f0[0] = 6.0*(x[0] + x[1]);
269: }

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

281: void f1_field_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
282:                 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
283:                 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
284:                 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
285: {
286:   PetscInt d;
287:   for (d = 0; d < dim; ++d) f1[d] = a[0]*u_x[d];
288: }

290: /* < \nabla v, \nabla u + {\nabla u}^T >
291:    This just gives \nabla u, give the perdiagonal for the transpose */
292: void g3_analytic_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
293:                     const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
294:                     const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
295:                     PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
296: {
297:   PetscInt d;
298:   for (d = 0; d < dim; ++d) g3[d*dim+d] = x[0] + x[1];
299: }

301: void g3_field_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
302:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
303:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
304:                  PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
305: {
306:   PetscInt d;
307:   for (d = 0; d < dim; ++d) g3[d*dim+d] = a[0];
308: }

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

313:     u  = x^2 + y^2
314:     f  = 16 (x^2 + y^2)
315:     nu = 1/2 |grad u|^2

317:   so that

319:     -\div \nu \grad u + f = -16 (x^2 + y^2) + 16 (x^2 + y^2) = 0
320: */
321: void f0_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
322:                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
323:                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
324:                              PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
325: {
326:   f0[0] = 16.0*(x[0]*x[0] + x[1]*x[1]);
327: }

329: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
330: void f1_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
331:                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
332:                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
333:                              PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
334: {
335:   PetscScalar nu = 0.0;
336:   PetscInt    d;
337:   for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
338:   for (d = 0; d < dim; ++d) f1[d] = 0.5*nu*u_x[d];
339: }

341: /*
342:   grad (u + eps w) - grad u = eps grad w

344:   1/2 |grad (u + eps w)|^2 grad (u + eps w) - 1/2 |grad u|^2 grad u
345: = 1/2 (|grad u|^2 + 2 eps <grad u,grad w>) (grad u + eps grad w) - 1/2 |grad u|^2 grad u
346: = 1/2 (eps |grad u|^2 grad w + 2 eps <grad u,grad w> grad u)
347: = eps (1/2 |grad u|^2 grad w + grad u <grad u,grad w>)
348: */
349: void g3_analytic_nonlinear_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
350:                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
351:                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
352:                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
353: {
354:   PetscScalar nu = 0.0;
355:   PetscInt    d, e;
356:   for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
357:   for (d = 0; d < dim; ++d) {
358:     g3[d*dim+d] = 0.5*nu;
359:     for (e = 0; e < dim; ++e) {
360:       g3[d*dim+e] += u_x[d]*u_x[e];
361:     }
362:   }
363: }

365: /*
366:   In 3D for Dirichlet conditions we use exact solution:

368:     u = 2/3 (x^2 + y^2 + z^2)
369:     f = 4

371:   so that

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

375:   For Neumann conditions, we have

377:     -\nabla u \cdot -\hat z |_{z=0} =  (2z)|_{z=0} =  0 (bottom)
378:     -\nabla u \cdot  \hat z |_{z=1} = -(2z)|_{z=1} = -2 (top)
379:     -\nabla u \cdot -\hat y |_{y=0} =  (2y)|_{y=0} =  0 (front)
380:     -\nabla u \cdot  \hat y |_{y=1} = -(2y)|_{y=1} = -2 (back)
381:     -\nabla u \cdot -\hat x |_{x=0} =  (2x)|_{x=0} =  0 (left)
382:     -\nabla u \cdot  \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)

384:   Which we can express as

386:     \nabla u \cdot  \hat n|_\Gamma = {2 x, 2 y, 2z} \cdot \hat n = 2 (x + y + z)
387: */
388: static PetscErrorCode quadratic_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
389: {
390:   *u = 2.0*(x[0]*x[0] + x[1]*x[1] + x[2]*x[2])/3.0;
391:   return 0;
392: }

394: static void quadratic_u_field_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
395:                                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
396:                                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
397:                                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
398: {
399:   uexact[0] = a[0];
400: }

402: static void bd_integral_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
403:                            const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
404:                            const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
405:                            PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar *uint)
406: {
407:   uint[0] = u[0];
408: }

410: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
411: {
412:   const char    *bcTypes[3]  = {"neumann", "dirichlet", "none"};
413:   const char    *runTypes[4] = {"full", "exact", "test", "perf"};
414:   const char    *coeffTypes[6] = {"none", "analytic", "field", "nonlinear", "circle", "cross"};
415:   PetscInt       bd, bc, run, coeff, n;
416:   PetscBool      flg;

420:   options->debug               = 0;
421:   options->runType             = RUN_FULL;
422:   options->dim                 = 2;
423:   options->periodicity[0]      = DM_BOUNDARY_NONE;
424:   options->periodicity[1]      = DM_BOUNDARY_NONE;
425:   options->periodicity[2]      = DM_BOUNDARY_NONE;
426:   options->cells[0]            = 2;
427:   options->cells[1]            = 2;
428:   options->cells[2]            = 2;
429:   options->filename[0]         = '\0';
430:   options->interpolate         = PETSC_FALSE;
431:   options->refinementLimit     = 0.0;
432:   options->bcType              = DIRICHLET;
433:   options->variableCoefficient = COEFF_NONE;
434:   options->fieldBC             = PETSC_FALSE;
435:   options->jacobianMF          = PETSC_FALSE;
436:   options->showInitial         = PETSC_FALSE;
437:   options->showSolution        = PETSC_FALSE;
438:   options->restart             = PETSC_FALSE;
439:   options->viewHierarchy       = PETSC_FALSE;
440:   options->simplex             = PETSC_TRUE;
441:   options->quiet               = PETSC_FALSE;
442:   options->nonzInit            = PETSC_FALSE;
443:   options->bdIntegral          = PETSC_FALSE;

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

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

452:   PetscOptionsInt("-dim", "The topological mesh dimension", "ex12.c", options->dim, &options->dim, NULL);
453:   bd = options->periodicity[0];
454:   PetscOptionsEList("-x_periodicity", "The x-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[0]], &bd, NULL);
455:   options->periodicity[0] = (DMBoundaryType) bd;
456:   bd = options->periodicity[1];
457:   PetscOptionsEList("-y_periodicity", "The y-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[1]], &bd, NULL);
458:   options->periodicity[1] = (DMBoundaryType) bd;
459:   bd = options->periodicity[2];
460:   PetscOptionsEList("-z_periodicity", "The z-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[2]], &bd, NULL);
461:   options->periodicity[2] = (DMBoundaryType) bd;
462:   n = 3;
463:   PetscOptionsIntArray("-cells", "The initial mesh division", "ex12.c", options->cells, &n, NULL);
464:   PetscOptionsString("-f", "Mesh filename to read", "ex12.c", options->filename, options->filename, sizeof(options->filename), &flg);
465:   PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex12.c", options->interpolate, &options->interpolate, NULL);
466:   PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex12.c", options->refinementLimit, &options->refinementLimit, NULL);
467:   bc   = options->bcType;
468:   PetscOptionsEList("-bc_type","Type of boundary condition","ex12.c",bcTypes,3,bcTypes[options->bcType],&bc,NULL);
469:   options->bcType = (BCType) bc;
470:   coeff = options->variableCoefficient;
471:   PetscOptionsEList("-variable_coefficient","Type of variable coefficent","ex12.c",coeffTypes,6,coeffTypes[options->variableCoefficient],&coeff,NULL);
472:   options->variableCoefficient = (CoeffType) coeff;

474:   PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL);
475:   PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL);
476:   PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL);
477:   PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL);
478:   PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL);
479:   PetscOptionsBool("-dm_view_hierarchy", "View the coarsened hierarchy", "ex12.c", options->viewHierarchy, &options->viewHierarchy, NULL);
480:   PetscOptionsBool("-simplex", "Simplicial (true) or tensor (false) mesh", "ex12.c", options->simplex, &options->simplex, NULL);
481:   PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL);
482:   PetscOptionsBool("-nonzero_initial_guess", "nonzero intial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL);
483:   PetscOptionsBool("-bd_integral", "Compute the integral of the solution on the boundary", "ex12.c", options->bdIntegral, &options->bdIntegral, NULL);
484:   PetscOptionsEnd();
485:   PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
486:   return(0);
487: }

489: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
490: {
491:   DMLabel        label;

495:   DMCreateLabel(dm, name);
496:   DMGetLabel(dm, name, &label);
497:   DMPlexMarkBoundaryFaces(dm, 1, label);
498:   DMPlexLabelComplete(dm, label);
499:   return(0);
500: }

502: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
503: {
504:   PetscInt       dim             = user->dim;
505:   const char    *filename        = user->filename;
506:   PetscBool      interpolate     = user->interpolate;
507:   PetscReal      refinementLimit = user->refinementLimit;
508:   size_t         len;

512:   PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
513:   PetscStrlen(filename, &len);
514:   if (!len) {
515:     PetscInt d;

517:     if (user->periodicity[0] || user->periodicity[1] || user->periodicity[2]) for (d = 0; d < dim; ++d) user->cells[d] = PetscMax(user->cells[d], 3);
518:     DMPlexCreateBoxMesh(comm, dim, user->simplex, user->cells, NULL, NULL, user->periodicity, interpolate, dm);
519:     PetscObjectSetName((PetscObject) *dm, "Mesh");
520:   } else {
521:     DMPlexCreateFromFile(comm, filename, interpolate, dm);
522:     DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
523:   }
524:   {
525:     PetscPartitioner part;
526:     DM               refinedMesh     = NULL;
527:     DM               distributedMesh = NULL;

529:     /* Refine mesh using a volume constraint */
530:     if (refinementLimit > 0.0) {
531:       DMPlexSetRefinementLimit(*dm, refinementLimit);
532:       DMRefine(*dm, comm, &refinedMesh);
533:       if (refinedMesh) {
534:         const char *name;

536:         PetscObjectGetName((PetscObject) *dm,         &name);
537:         PetscObjectSetName((PetscObject) refinedMesh,  name);
538:         DMDestroy(dm);
539:         *dm  = refinedMesh;
540:       }
541:     }
542:     /* Distribute mesh over processes */
543:     DMPlexGetPartitioner(*dm,&part);
544:     PetscPartitionerSetFromOptions(part);
545:     DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
546:     if (distributedMesh) {
547:       DMDestroy(dm);
548:       *dm  = distributedMesh;
549:     }
550:   }
551:   if (user->bcType == NEUMANN) {
552:     DMLabel   label;

554:     DMCreateLabel(*dm, "boundary");
555:     DMGetLabel(*dm, "boundary", &label);
556:     DMPlexMarkBoundaryFaces(*dm, 1, label);
557:   } else if (user->bcType == DIRICHLET) {
558:     PetscBool hasLabel;

560:     DMHasLabel(*dm,"marker",&hasLabel);
561:     if (!hasLabel) {CreateBCLabel(*dm, "marker");}
562:   }
563:   {
564:     char      convType[256];
565:     PetscBool flg;

567:     PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
568:     PetscOptionsFList("-dm_plex_convert_type","Convert DMPlex to another format","ex12",DMList,DMPLEX,convType,256,&flg);
569:     PetscOptionsEnd();
570:     if (flg) {
571:       DM dmConv;

573:       DMConvert(*dm,convType,&dmConv);
574:       if (dmConv) {
575:         DMDestroy(dm);
576:         *dm  = dmConv;
577:       }
578:     }
579:   }
580:   DMLocalizeCoordinates(*dm); /* needed for periodic */
581:   DMSetFromOptions(*dm);
582:   DMViewFromOptions(*dm, NULL, "-dm_view");
583:   if (user->viewHierarchy) {
584:     DM       cdm = *dm;
585:     PetscInt i   = 0;
586:     char     buf[256];

588:     while (cdm) {
589:       DMSetUp(cdm);
590:       DMGetCoarseDM(cdm, &cdm);
591:       ++i;
592:     }
593:     cdm = *dm;
594:     while (cdm) {
595:       PetscViewer       viewer;
596:       PetscBool   isHDF5, isVTK;

598:       --i;
599:       PetscViewerCreate(comm,&viewer);
600:       PetscViewerSetType(viewer,PETSCVIEWERHDF5);
601:       PetscViewerSetOptionsPrefix(viewer,"hierarchy_");
602:       PetscViewerSetFromOptions(viewer);
603:       PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&isHDF5);
604:       PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERVTK,&isVTK);
605:       if (isHDF5) {
606:         PetscSNPrintf(buf, 256, "ex12-%d.h5", i);
607:       } else if (isVTK) {
608:         PetscSNPrintf(buf, 256, "ex12-%d.vtu", i);
609:         PetscViewerPushFormat(viewer,PETSC_VIEWER_VTK_VTU);
610:       } else {
611:         PetscSNPrintf(buf, 256, "ex12-%d", i);
612:       }
613:       PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
614:       PetscViewerFileSetName(viewer,buf);
615:       DMView(cdm, viewer);
616:       PetscViewerDestroy(&viewer);
617:       DMGetCoarseDM(cdm, &cdm);
618:     }
619:   }
620:   PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
621:   return(0);
622: }

624: static PetscErrorCode SetupProblem(DM dm, AppCtx *user)
625: {
626:   PetscDS        prob;
627:   const PetscInt id = 1;

631:   DMGetDS(dm, &prob);
632:   switch (user->variableCoefficient) {
633:   case COEFF_NONE:
634:     if (user->periodicity[0]) {
635:       if (user->periodicity[1]) {
636:         PetscDSSetResidual(prob, 0, f0_xytrig_u, f1_u);
637:         PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
638:       } else {
639:         PetscDSSetResidual(prob, 0, f0_xtrig_u,  f1_u);
640:         PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
641:       }
642:     } else {
643:       PetscDSSetResidual(prob, 0, f0_u, f1_u);
644:       PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
645:     }
646:     break;
647:   case COEFF_ANALYTIC:
648:     PetscDSSetResidual(prob, 0, f0_analytic_u, f1_analytic_u);
649:     PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_analytic_uu);
650:     break;
651:   case COEFF_FIELD:
652:     PetscDSSetResidual(prob, 0, f0_analytic_u, f1_field_u);
653:     PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_field_uu);
654:     break;
655:   case COEFF_NONLINEAR:
656:     PetscDSSetResidual(prob, 0, f0_analytic_nonlinear_u, f1_analytic_nonlinear_u);
657:     PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_analytic_nonlinear_uu);
658:     break;
659:   case COEFF_CIRCLE:
660:     PetscDSSetResidual(prob, 0, f0_circle_u, f1_u);
661:     PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
662:     break;
663:   case COEFF_CROSS:
664:     PetscDSSetResidual(prob, 0, f0_cross_u, f1_u);
665:     PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
666:     break;
667:   default: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid variable coefficient type %d", user->variableCoefficient);
668:   }
669:   switch (user->dim) {
670:   case 2:
671:     switch (user->variableCoefficient) {
672:     case COEFF_CIRCLE:
673:       user->exactFuncs[0]  = circle_u_2d;break;
674:     case COEFF_CROSS:
675:       user->exactFuncs[0]  = cross_u_2d;break;
676:     default:
677:       if (user->periodicity[0]) {
678:         if (user->periodicity[1]) {
679:           user->exactFuncs[0] = xytrig_u_2d;
680:         } else {
681:           user->exactFuncs[0] = xtrig_u_2d;
682:         }
683:       } else {
684:         user->exactFuncs[0]  = quadratic_u_2d;
685:         user->exactFields[0] = quadratic_u_field_2d;
686:       }
687:     }
688:     if (user->bcType == NEUMANN) {PetscDSSetBdResidual(prob, 0, f0_bd_u, f1_bd_zero);}
689:     break;
690:   case 3:
691:     user->exactFuncs[0]  = quadratic_u_3d;
692:     user->exactFields[0] = quadratic_u_field_3d;
693:     if (user->bcType == NEUMANN) {PetscDSSetBdResidual(prob, 0, f0_bd_u, f1_bd_zero);}
694:     break;
695:   default:
696:     SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %d", user->dim);
697:   }
698:   PetscDSAddBoundary(prob, user->bcType == DIRICHLET ? (user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL) : DM_BC_NATURAL,
699:                             "wall", user->bcType == DIRICHLET ? "marker" : "boundary", 0, 0, NULL,
700:                             user->fieldBC ? (void (*)(void)) user->exactFields[0] : (void (*)(void)) user->exactFuncs[0], 1, &id, user);
701:   PetscDSSetExactSolution(prob, 0, user->exactFuncs[0], user);
702:   return(0);
703: }

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

712:   DMCreateLocalVector(dmAux, &nu);
713:   DMProjectFunctionLocal(dmAux, 0.0, matFuncs, NULL, INSERT_ALL_VALUES, nu);
714:   PetscObjectCompose((PetscObject) dm, "A", (PetscObject) nu);
715:   VecDestroy(&nu);
716:   return(0);
717: }

719: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
720: {
721:   PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
722:   Vec            uexact;
723:   PetscInt       dim;

727:   DMGetDimension(dm, &dim);
728:   if (dim == 2) bcFuncs[0] = quadratic_u_2d;
729:   else          bcFuncs[0] = quadratic_u_3d;
730:   DMCreateLocalVector(dmAux, &uexact);
731:   DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
732:   PetscObjectCompose((PetscObject) dm, "A", (PetscObject) uexact);
733:   VecDestroy(&uexact);
734:   return(0);
735: }

737: static PetscErrorCode SetupAuxDM(DM dm, PetscFE feAux, AppCtx *user)
738: {
739:   DM             dmAux, coordDM;

743:   /* MUST call DMGetCoordinateDM() in order to get p4est setup if present */
744:   DMGetCoordinateDM(dm, &coordDM);
745:   if (!feAux) return(0);
746:   DMClone(dm, &dmAux);
747:   PetscObjectCompose((PetscObject) dm, "dmAux", (PetscObject) dmAux);
748:   DMSetCoordinateDM(dmAux, coordDM);
749:   DMSetField(dmAux, 0, NULL, (PetscObject) feAux);
750:   DMCreateDS(dmAux);
751:   if (user->fieldBC) {SetupBC(dm, dmAux, user);}
752:   else               {SetupMaterial(dm, dmAux, user);}
753:   DMDestroy(&dmAux);
754:   return(0);
755: }

757: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
758: {
759:   DM             cdm = dm;
760:   const PetscInt dim = user->dim;
761:   PetscFE        fe, feAux = NULL;
762:   PetscBool      simplex   = user->simplex;
763:   MPI_Comm       comm;

767:   /* Create finite element for each field and auxiliary field */
768:   PetscObjectGetComm((PetscObject) dm, &comm);
769:   PetscFECreateDefault(comm, dim, 1, simplex, NULL, -1, &fe);
770:   PetscObjectSetName((PetscObject) fe, "potential");
771:   if (user->variableCoefficient == COEFF_FIELD) {
772:     PetscQuadrature q;

774:     PetscFECreateDefault(comm, dim, 1, simplex, "mat_", -1, &feAux);
775:     PetscFEGetQuadrature(fe, &q);
776:     PetscFESetQuadrature(feAux, q);
777:   } else if (user->fieldBC) {
778:     PetscQuadrature q;

780:     PetscFECreateDefault(comm, dim, 1, simplex, "bc_", -1, &feAux);
781:     PetscFEGetQuadrature(fe, &q);
782:     PetscFESetQuadrature(feAux, q);
783:   }
784:   /* Set discretization and boundary conditions for each mesh */
785:   DMSetField(dm, 0, NULL, (PetscObject) fe);
786:   DMCreateDS(dm);
787:   SetupProblem(dm, user);
788:   while (cdm) {
789:     DMCopyDisc(dm, cdm);
790:     SetupAuxDM(cdm, feAux, user);
791:     if (user->bcType == DIRICHLET) {
792:       PetscBool hasLabel;

794:       DMHasLabel(cdm, "marker", &hasLabel);
795:       if (!hasLabel) {CreateBCLabel(cdm, "marker");}
796:     }
797:     DMGetCoarseDM(cdm, &cdm);
798:   }
799:   PetscFEDestroy(&fe);
800:   PetscFEDestroy(&feAux);
801:   return(0);
802: }

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

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

809:   Collective on KSP

811:   Input Parameters:
812: + ksp   - the KSP
813: . its   - iteration number
814: . rnorm - 2-norm, preconditioned residual value (may be estimated).
815: - ctx   - monitor context

817:   Level: intermediate

819: .seealso: KSPMonitorSet(), KSPMonitorTrueResidualNorm(), KSPMonitorDefault()
820: @*/
821: static PetscErrorCode KSPMonitorError(KSP ksp, PetscInt its, PetscReal rnorm, void *ctx)
822: {
823:   AppCtx        *user = (AppCtx *) ctx;
824:   DM             dm;
825:   Vec            du = NULL, r;
826:   PetscInt       level = 0;
827:   PetscBool      hasLevel;
828: #if defined(PETSC_HAVE_HDF5)
829:   PetscViewer    viewer;
830:   char           buf[256];
831: #endif

835:   KSPGetDM(ksp, &dm);
836:   /* Calculate solution */
837:   {
838:     PC        pc = user->pcmg; /* The MG PC */
839:     DM        fdm = NULL,  cdm = NULL;
840:     KSP       fksp, cksp;
841:     Vec       fu,   cu = NULL;
842:     PetscInt  levels, l;

844:     KSPBuildSolution(ksp, NULL, &du);
845:     PetscObjectComposedDataGetInt((PetscObject) ksp, PetscMGLevelId, level, hasLevel);
846:     PCMGGetLevels(pc, &levels);
847:     PCMGGetSmoother(pc, levels-1, &fksp);
848:     KSPBuildSolution(fksp, NULL, &fu);
849:     for (l = levels-1; l > level; --l) {
850:       Mat R;
851:       Vec s;

853:       PCMGGetSmoother(pc, l-1, &cksp);
854:       KSPGetDM(cksp, &cdm);
855:       DMGetGlobalVector(cdm, &cu);
856:       PCMGGetRestriction(pc, l, &R);
857:       PCMGGetRScale(pc, l, &s);
858:       MatRestrict(R, fu, cu);
859:       VecPointwiseMult(cu, cu, s);
860:       if (l < levels-1) {DMRestoreGlobalVector(fdm, &fu);}
861:       fdm  = cdm;
862:       fu   = cu;
863:     }
864:     if (levels-1 > level) {
865:       VecAXPY(du, 1.0, cu);
866:       DMRestoreGlobalVector(cdm, &cu);
867:     }
868:   }
869:   /* Calculate error */
870:   DMGetGlobalVector(dm, &r);
871:   DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
872:   VecAXPY(r,-1.0,du);
873:   PetscObjectSetName((PetscObject) r, "solution error");
874:   /* View error */
875: #if defined(PETSC_HAVE_HDF5)
876:   PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
877:   PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
878:   VecView(r, viewer);
879:   PetscViewerDestroy(&viewer);
880: #endif
881:   DMRestoreGlobalVector(dm, &r);
882:   return(0);
883: }

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

888:   Collective on SNES

890:   Input Parameters:
891: + snes  - the SNES
892: . its   - iteration number
893: . rnorm - 2-norm of residual
894: - ctx   - user context

896:   Level: intermediate

898: .seealso: SNESMonitorDefault(), SNESMonitorSet(), SNESMonitorSolution()
899: @*/
900: static PetscErrorCode SNESMonitorError(SNES snes, PetscInt its, PetscReal rnorm, void *ctx)
901: {
902:   AppCtx        *user = (AppCtx *) ctx;
903:   DM             dm;
904:   Vec            u, r;
905:   PetscInt       level = -1;
906:   PetscBool      hasLevel;
907: #if defined(PETSC_HAVE_HDF5)
908:   PetscViewer    viewer;
909: #endif
910:   char           buf[256];

914:   SNESGetDM(snes, &dm);
915:   /* Calculate error */
916:   SNESGetSolution(snes, &u);
917:   DMGetGlobalVector(dm, &r);
918:   PetscObjectSetName((PetscObject) r, "solution error");
919:   DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
920:   VecAXPY(r, -1.0, u);
921:   /* View error */
922:   PetscObjectComposedDataGetInt((PetscObject) snes, PetscMGLevelId, level, hasLevel);
923:   PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
924: #if defined(PETSC_HAVE_HDF5)
925:   PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
926:   VecView(r, viewer);
927:   PetscViewerDestroy(&viewer);
928:   /* Cleanup */
929:   DMRestoreGlobalVector(dm, &r);
930:   return(0);
931: #else
932:   SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"You need to configure with --download-hdf5");
933: #endif
934: }

936: int main(int argc, char **argv)
937: {
938:   DM             dm;          /* Problem specification */
939:   SNES           snes;        /* nonlinear solver */
940:   Vec            u;           /* solution vector */
941:   Mat            A,J;         /* Jacobian matrix */
942:   MatNullSpace   nullSpace;   /* May be necessary for Neumann conditions */
943:   AppCtx         user;        /* user-defined work context */
944:   JacActionCtx   userJ;       /* context for Jacobian MF action */
945:   PetscReal      error = 0.0; /* L_2 error in the solution */
946:   PetscBool      isFAS;

949:   PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
950:   ProcessOptions(PETSC_COMM_WORLD, &user);
951:   SNESCreate(PETSC_COMM_WORLD, &snes);
952:   CreateMesh(PETSC_COMM_WORLD, &user, &dm);
953:   SNESSetDM(snes, dm);
954:   DMSetApplicationContext(dm, &user);

956:   PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
957:   SetupDiscretization(dm, &user);

959:   DMCreateGlobalVector(dm, &u);
960:   PetscObjectSetName((PetscObject) u, "potential");

962:   DMCreateMatrix(dm, &J);
963:   if (user.jacobianMF) {
964:     PetscInt M, m, N, n;

966:     MatGetSize(J, &M, &N);
967:     MatGetLocalSize(J, &m, &n);
968:     MatCreate(PETSC_COMM_WORLD, &A);
969:     MatSetSizes(A, m, n, M, N);
970:     MatSetType(A, MATSHELL);
971:     MatSetUp(A);
972: #if 0
973:     MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction);
974: #endif

976:     userJ.dm   = dm;
977:     userJ.J    = J;
978:     userJ.user = &user;

980:     DMCreateLocalVector(dm, &userJ.u);
981:     if (user.fieldBC) {DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);}
982:     else              {DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);}
983:     MatShellSetContext(A, &userJ);
984:   } else {
985:     A = J;
986:   }
987:   if (user.bcType == NEUMANN) {
988:     MatNullSpaceCreate(PetscObjectComm((PetscObject) dm), PETSC_TRUE, 0, NULL, &nullSpace);
989:     MatSetNullSpace(A, nullSpace);
990:   }

992:   DMPlexSetSNESLocalFEM(dm,&user,&user,&user);
993:   SNESSetJacobian(snes, A, J, NULL, NULL);

995:   SNESSetFromOptions(snes);

997:   if (user.fieldBC) {DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u);}
998:   else              {DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);}
999:   if (user.restart) {
1000: #if defined(PETSC_HAVE_HDF5)
1001:     PetscViewer viewer;

1003:     PetscViewerCreate(PETSC_COMM_WORLD, &viewer);
1004:     PetscViewerSetType(viewer, PETSCVIEWERHDF5);
1005:     PetscViewerFileSetMode(viewer, FILE_MODE_READ);
1006:     PetscViewerFileSetName(viewer, user.filename);
1007:     PetscViewerHDF5PushGroup(viewer, "/fields");
1008:     VecLoad(u, viewer);
1009:     PetscViewerHDF5PopGroup(viewer);
1010:     PetscViewerDestroy(&viewer);
1011: #endif
1012:   }
1013:   if (user.showInitial) {
1014:     Vec lv;
1015:     DMGetLocalVector(dm, &lv);
1016:     DMGlobalToLocalBegin(dm, u, INSERT_VALUES, lv);
1017:     DMGlobalToLocalEnd(dm, u, INSERT_VALUES, lv);
1018:     DMPrintLocalVec(dm, "Local function", 1.0e-10, lv);
1019:     DMRestoreLocalVector(dm, &lv);
1020:   }
1021:   if (user.viewHierarchy) {
1022:     SNES      lsnes;
1023:     KSP       ksp;
1024:     PC        pc;
1025:     PetscInt  numLevels, l;
1026:     PetscBool isMG;

1028:     PetscObjectTypeCompare((PetscObject) snes, SNESFAS, &isFAS);
1029:     if (isFAS) {
1030:       SNESFASGetLevels(snes, &numLevels);
1031:       for (l = 0; l < numLevels; ++l) {
1032:         SNESFASGetCycleSNES(snes, l, &lsnes);
1033:         SNESMonitorSet(lsnes, SNESMonitorError, &user, NULL);
1034:       }
1035:     } else {
1036:       SNESGetKSP(snes, &ksp);
1037:       KSPGetPC(ksp, &pc);
1038:       PetscObjectTypeCompare((PetscObject) pc, PCMG, &isMG);
1039:       if (isMG) {
1040:         user.pcmg = pc;
1041:         PCMGGetLevels(pc, &numLevels);
1042:         for (l = 0; l < numLevels; ++l) {
1043:           PCMGGetSmootherDown(pc, l, &ksp);
1044:           KSPMonitorSet(ksp, KSPMonitorError, &user, NULL);
1045:         }
1046:       }
1047:     }
1048:   }
1049:   if (user.runType == RUN_FULL || user.runType == RUN_EXACT) {
1050:     PetscErrorCode (*initialGuess[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {zero};

1052:     if (user.nonzInit) initialGuess[0] = ecks;
1053:     if (user.runType == RUN_FULL) {
1054:       DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
1055:     }
1056:     if (user.debug) {
1057:       PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1058:       VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1059:     }
1060:     SNESSolve(snes, NULL, u);
1061:     SNESGetSolution(snes, &u);
1062:     SNESGetDM(snes, &dm);

1064:     if (user.showSolution) {
1065:       PetscPrintf(PETSC_COMM_WORLD, "Solution\n");
1066:       VecChop(u, 3.0e-9);
1067:       VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1068:     }
1069:     VecViewFromOptions(u, NULL, "-vec_view");
1070:   } else if (user.runType == RUN_PERF) {
1071:     Vec       r;
1072:     PetscReal res = 0.0;

1074:     SNESGetFunction(snes, &r, NULL, NULL);
1075:     SNESComputeFunction(snes, u, r);
1076:     PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1077:     VecChop(r, 1.0e-10);
1078:     VecNorm(r, NORM_2, &res);
1079:     PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1080:   } else {
1081:     Vec       r;
1082:     PetscReal res = 0.0, tol = 1.0e-11;

1084:     /* Check discretization error */
1085:     SNESGetFunction(snes, &r, NULL, NULL);
1086:     PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1087:     if (!user.quiet) {VecView(u, PETSC_VIEWER_STDOUT_WORLD);}
1088:     DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
1089:     if (error < tol) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < %2.1e\n", tol);}
1090:     else             {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", error);}
1091:     /* Check residual */
1092:     SNESComputeFunction(snes, u, r);
1093:     PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1094:     VecChop(r, 1.0e-10);
1095:     if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1096:     VecNorm(r, NORM_2, &res);
1097:     PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1098:     /* Check Jacobian */
1099:     {
1100:       Vec b;

1102:       SNESComputeJacobian(snes, u, A, A);
1103:       VecDuplicate(u, &b);
1104:       VecSet(r, 0.0);
1105:       SNESComputeFunction(snes, r, b);
1106:       MatMult(A, u, r);
1107:       VecAXPY(r, 1.0, b);
1108:       VecDestroy(&b);
1109:       PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n");
1110:       VecChop(r, 1.0e-10);
1111:       if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1112:       VecNorm(r, NORM_2, &res);
1113:       PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", res);
1114:     }
1115:   }
1116:   VecViewFromOptions(u, NULL, "-vec_view");

1118:   if (user.bdIntegral) {
1119:     DMLabel   label;
1120:     PetscInt  id = 1;
1121:     PetscScalar bdInt = 0.0;
1122:     PetscReal   exact = 3.3333333333;

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

1130:   if (user.bcType == NEUMANN) {MatNullSpaceDestroy(&nullSpace);}
1131:   if (user.jacobianMF) {VecDestroy(&userJ.u);}
1132:   if (A != J) {MatDestroy(&A);}
1133:   MatDestroy(&J);
1134:   VecDestroy(&u);
1135:   SNESDestroy(&snes);
1136:   DMDestroy(&dm);
1137:   PetscFree2(user.exactFuncs, user.exactFields);
1138:   PetscFinalize();
1139:   return ierr;
1140: }

1142: /*TEST
1143:   # 2D serial P1 test 0-4
1144:   test:
1145:     suffix: 2d_p1_0
1146:     requires: triangle
1147:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1149:   test:
1150:     suffix: 2d_p1_1
1151:     requires: triangle
1152:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1154:   test:
1155:     suffix: 2d_p1_2
1156:     requires: triangle
1157:     args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1159:   test:
1160:     suffix: 2d_p1_neumann_0
1161:     requires: triangle
1162:     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

1164:   test:
1165:     suffix: 2d_p1_neumann_1
1166:     requires: triangle
1167:     args: -run_type test -refinement_limit 0.0625 -bc_type neumann   -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1169:   # 2D serial P2 test 5-8
1170:   test:
1171:     suffix: 2d_p2_0
1172:     requires: triangle
1173:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1175:   test:
1176:     suffix: 2d_p2_1
1177:     requires: triangle
1178:     args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1180:   test:
1181:     suffix: 2d_p2_neumann_0
1182:     requires: triangle
1183:     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

1185:   test:
1186:     suffix: 2d_p2_neumann_1
1187:     requires: triangle
1188:     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

1190:   test:
1191:     suffix: bd_int_0
1192:     requires: triangle
1193:     args: -run_type test -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet

1195:   test:
1196:     suffix: bd_int_1
1197:     requires: triangle
1198:     args: -run_type test -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet

1200:   # 3D serial P1 test 9-12
1201:   test:
1202:     suffix: 3d_p1_0
1203:     requires: ctetgen
1204:     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

1206:   test:
1207:     suffix: 3d_p1_1
1208:     requires: ctetgen
1209:     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

1211:   test:
1212:     suffix: 3d_p1_2
1213:     requires: ctetgen
1214:     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

1216:   test:
1217:     suffix: 3d_p1_neumann_0
1218:     requires: ctetgen
1219:     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

1221:   # Analytic variable coefficient 13-20
1222:   test:
1223:     suffix: 13
1224:     requires: triangle
1225:     args: -run_type test -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1226:   test:
1227:     suffix: 14
1228:     requires: triangle
1229:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1230:   test:
1231:     suffix: 15
1232:     requires: triangle
1233:     args: -run_type test -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1234:   test:
1235:     suffix: 16
1236:     requires: triangle
1237:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1238:   test:
1239:     suffix: 17
1240:     requires: ctetgen
1241:     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

1243:   test:
1244:     suffix: 18
1245:     requires: ctetgen
1246:     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

1248:   test:
1249:     suffix: 19
1250:     requires: ctetgen
1251:     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

1253:   test:
1254:     suffix: 20
1255:     requires: ctetgen
1256:     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

1258:   # P1 variable coefficient 21-28
1259:   test:
1260:     suffix: 21
1261:     requires: triangle
1262:     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

1264:   test:
1265:     suffix: 22
1266:     requires: triangle
1267:     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

1269:   test:
1270:     suffix: 23
1271:     requires: triangle
1272:     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

1274:   test:
1275:     suffix: 24
1276:     requires: triangle
1277:     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

1279:   test:
1280:     suffix: 25
1281:     requires: ctetgen
1282:     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

1284:   test:
1285:     suffix: 26
1286:     requires: ctetgen
1287:     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

1289:   test:
1290:     suffix: 27
1291:     requires: ctetgen
1292:     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

1294:   test:
1295:     suffix: 28
1296:     requires: ctetgen
1297:     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

1299:   # P0 variable coefficient 29-36
1300:   test:
1301:     suffix: 29
1302:     requires: triangle
1303:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1305:   test:
1306:     suffix: 30
1307:     requires: triangle
1308:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1310:   test:
1311:     suffix: 31
1312:     requires: triangle
1313:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1315:   test:
1316:     requires: triangle
1317:     suffix: 32
1318:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1320:   test:
1321:     requires: ctetgen
1322:     suffix: 33
1323:     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

1325:   test:
1326:     suffix: 34
1327:     requires: ctetgen
1328:     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

1330:   test:
1331:     suffix: 35
1332:     requires: ctetgen
1333:     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

1335:   test:
1336:     suffix: 36
1337:     requires: ctetgen
1338:     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

1340:   # Full solve 39-44
1341:   test:
1342:     suffix: 39
1343:     requires: triangle !single
1344:     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
1345:   test:
1346:     suffix: 40
1347:     requires: triangle !single
1348:     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
1349:   test:
1350:     suffix: 41
1351:     requires: triangle !single
1352:     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
1353:   test:
1354:     suffix: 42
1355:     requires: triangle !single
1356:     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
1357:   test:
1358:     suffix: 43
1359:     requires: triangle !single
1360:     nsize: 2
1361:     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

1363:   test:
1364:     suffix: 44
1365:     requires: triangle !single
1366:     nsize: 2
1367:     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

1369:   # These tests use a loose tolerance just to exercise the PtAP operations for MATIS and multiple BDDC setup calls inside PCMG
1370:   testset:
1371:     requires: triangle !single
1372:     nsize: 3
1373:     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
1374:     test:
1375:       suffix: gmg_bddc
1376:       filter: sed -e "s/CONVERGED_FNORM_RELATIVE iterations 3/CONVERGED_FNORM_RELATIVE iterations 4/g"
1377:       args: -mg_levels_pc_type jacobi
1378:     test:
1379:       filter: sed -e "s/iterations [0-4]/iterations 4/g"
1380:       suffix: gmg_bddc_lev
1381:       args: -mg_levels_pc_type bddc

1383:   # Restarting
1384:   testset:
1385:     suffix: restart
1386:     requires: hdf5 triangle !complex
1387:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 1
1388:     test:
1389:       args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1390:     test:
1391:       args: -f sol.h5 -restart

1393:   # Periodicity
1394:   test:
1395:     suffix: periodic_0
1396:     requires: triangle
1397:     args: -run_type full -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -snes_converged_reason ::ascii_info_detail

1399:   test:
1400:     requires: !complex
1401:     suffix: periodic_1
1402:     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

1404:   # 2D serial P1 test with field bc
1405:   test:
1406:     suffix: field_bc_2d_p1_0
1407:     requires: triangle
1408:     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

1410:   test:
1411:     suffix: field_bc_2d_p1_1
1412:     requires: triangle
1413:     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

1415:   test:
1416:     suffix: field_bc_2d_p1_neumann_0
1417:     requires: triangle
1418:     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

1420:   test:
1421:     suffix: field_bc_2d_p1_neumann_1
1422:     requires: triangle
1423:     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

1425:   # 3D serial P1 test with field bc
1426:   test:
1427:     suffix: field_bc_3d_p1_0
1428:     requires: ctetgen
1429:     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

1431:   test:
1432:     suffix: field_bc_3d_p1_1
1433:     requires: ctetgen
1434:     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

1436:   test:
1437:     suffix: field_bc_3d_p1_neumann_0
1438:     requires: ctetgen
1439:     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

1441:   test:
1442:     suffix: field_bc_3d_p1_neumann_1
1443:     requires: ctetgen
1444:     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

1446:   # 2D serial P2 test with field bc
1447:   test:
1448:     suffix: field_bc_2d_p2_0
1449:     requires: triangle
1450:     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

1452:   test:
1453:     suffix: field_bc_2d_p2_1
1454:     requires: triangle
1455:     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

1457:   test:
1458:     suffix: field_bc_2d_p2_neumann_0
1459:     requires: triangle
1460:     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

1462:   test:
1463:     suffix: field_bc_2d_p2_neumann_1
1464:     requires: triangle
1465:     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

1467:   # 3D serial P2 test with field bc
1468:   test:
1469:     suffix: field_bc_3d_p2_0
1470:     requires: ctetgen
1471:     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

1473:   test:
1474:     suffix: field_bc_3d_p2_1
1475:     requires: ctetgen
1476:     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

1478:   test:
1479:     suffix: field_bc_3d_p2_neumann_0
1480:     requires: ctetgen
1481:     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

1483:   test:
1484:     suffix: field_bc_3d_p2_neumann_1
1485:     requires: ctetgen
1486:     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

1488:   # Full solve simplex: Convergence
1489:   test:
1490:     suffix: tet_conv_p1_r0
1491:     requires: ctetgen
1492:     args: -run_type full -dim 3 -dm_refine 0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1
1493:   test:
1494:     suffix: tet_conv_p1_r2
1495:     requires: ctetgen
1496:     args: -run_type full -dim 3 -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1
1497:   test:
1498:     suffix: tet_conv_p1_r3
1499:     requires: ctetgen
1500:     args: -run_type full -dim 3 -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1
1501:   test:
1502:     suffix: tet_conv_p2_r0
1503:     requires: ctetgen
1504:     args: -run_type full -dim 3 -dm_refine 0 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1
1505:   test:
1506:     suffix: tet_conv_p2_r2
1507:     requires: ctetgen
1508:     args: -run_type full -dim 3 -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1

1510:   # Full solve simplex: BDDC
1511:   test:
1512:     suffix: tri_bddc
1513:     requires: triangle !single
1514:     nsize: 5
1515:     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

1517:   # Full solve simplex: BDDC
1518:   test:
1519:     suffix: tri_bddc_parmetis
1520:     requires: triangle !single parmetis
1521:     nsize: 4
1522:     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

1524:   testset:
1525:     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
1526:     nsize: 5
1527:     output_file: output/ex12_quad_bddc.out
1528:     filter: sed -e "s/aijcusparse/aij/g" -e "s/aijviennacl/aij/g" -e "s/factorization: cusparse/factorization: petsc/g"
1529:     test:
1530:       requires: !single
1531:       suffix: quad_bddc
1532:     test:
1533:       requires: !single cuda
1534:       suffix: quad_bddc_cuda
1535:       args: -matis_localmat_type aijcusparse -pc_bddc_dirichlet_pc_factor_mat_solver_type cusparse -pc_bddc_neumann_pc_factor_mat_solver_type cusparse
1536:     test:
1537:       requires: !single viennacl
1538:       suffix: quad_bddc_viennacl
1539:       args: -matis_localmat_type aijviennacl

1541:   # Full solve simplex: ASM
1542:   test:
1543:     suffix: tri_q2q1_asm_lu
1544:     requires: triangle !single
1545:     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

1547:   test:
1548:     suffix: tri_q2q1_msm_lu
1549:     requires: triangle !single
1550:     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

1552:   test:
1553:     suffix: tri_q2q1_asm_sor
1554:     requires: triangle !single
1555:     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

1557:   test:
1558:     suffix: tri_q2q1_msm_sor
1559:     requires: triangle !single
1560:     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

1562:   # Full solve simplex: FAS
1563:   test:
1564:     suffix: fas_newton_0
1565:     requires: triangle !single
1566:     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

1568:   test:
1569:     suffix: fas_newton_1
1570:     requires: triangle !single
1571:     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

1573:   test:
1574:     suffix: fas_ngs_0
1575:     requires: triangle !single
1576:     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

1578:   test:
1579:     suffix: fas_newton_coarse_0
1580:     requires: pragmatic triangle
1581:     TODO: broken
1582:     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

1584:   test:
1585:     suffix: mg_newton_coarse_0
1586:     requires: triangle pragmatic
1587:     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

1589:   test:
1590:     suffix: mg_newton_coarse_1
1591:     requires: triangle pragmatic
1592:     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

1594:   test:
1595:     suffix: mg_newton_coarse_2
1596:     requires: triangle pragmatic
1597:     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

1599:   # Full solve tensor
1600:   test:
1601:     suffix: tensor_plex_2d
1602:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_refine_hierarchy 2 -cells 2,2

1604:   test:
1605:     suffix: tensor_p4est_2d
1606:     requires: p4est
1607:     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

1609:   test:
1610:     suffix: tensor_plex_3d
1611:     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

1613:   test:
1614:     suffix: tensor_p4est_3d
1615:     requires: p4est
1616:     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

1618:   test:
1619:     suffix: p4est_test_q2_conformal_serial
1620:     requires: p4est
1621:     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

1623:   test:
1624:     suffix: p4est_test_q2_conformal_parallel
1625:     requires: p4est
1626:     nsize: 7
1627:     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

1629:   test:
1630:     suffix: p4est_test_q2_conformal_parallel_parmetis
1631:     requires: parmetis p4est
1632:     nsize: 4
1633:     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

1635:   test:
1636:     suffix: p4est_test_q2_nonconformal_serial
1637:     requires: p4est
1638:     filter: grep -v "CG or CGNE: variant"
1639:     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

1641:   test:
1642:     suffix: p4est_test_q2_nonconformal_parallel
1643:     requires: p4est
1644:     filter: grep -v "CG or CGNE: variant"
1645:     nsize: 7
1646:     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

1648:   test:
1649:     suffix: p4est_test_q2_nonconformal_parallel_parmetis
1650:     requires: parmetis p4est
1651:     nsize: 4
1652:     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

1654:   test:
1655:     suffix: p4est_exact_q2_conformal_serial
1656:     requires: p4est !single !complex !__float128
1657:     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

1659:   test:
1660:     suffix: p4est_exact_q2_conformal_parallel
1661:     requires: p4est !single !complex !__float128
1662:     nsize: 4
1663:     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

1665:   test:
1666:     suffix: p4est_exact_q2_conformal_parallel_parmetis
1667:     requires: parmetis p4est !single
1668:     nsize: 4
1669:     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

1671:   test:
1672:     suffix: p4est_exact_q2_nonconformal_serial
1673:     requires: p4est
1674:     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

1676:   test:
1677:     suffix: p4est_exact_q2_nonconformal_parallel
1678:     requires: p4est
1679:     nsize: 7
1680:     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

1682:   test:
1683:     suffix: p4est_exact_q2_nonconformal_parallel_parmetis
1684:     requires: parmetis p4est
1685:     nsize: 4
1686:     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

1688:   test:
1689:     suffix: p4est_full_q2_nonconformal_serial
1690:     requires: p4est !single
1691:     filter: grep -v "variant HERMITIAN"
1692:     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

1694:   test:
1695:     suffix: p4est_full_q2_nonconformal_parallel
1696:     requires: p4est !single
1697:     filter: grep -v "variant HERMITIAN"
1698:     nsize: 7
1699:     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

1701:   test:
1702:     suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1703:     requires: p4est !single
1704:     filter: grep -v "variant HERMITIAN"
1705:     nsize: 7
1706:     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

1708:   test:
1709:     suffix: p4est_full_q2_nonconformal_parallel_bddc
1710:     requires: p4est !single
1711:     filter: grep -v "variant HERMITIAN"
1712:     nsize: 7
1713:     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

1715:   test:
1716:     TODO: broken
1717:     suffix: p4est_fas_q2_conformal_serial
1718:     requires: p4est !complex !__float128
1719:     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

1721:   test:
1722:     TODO: broken
1723:     suffix: p4est_fas_q2_nonconformal_serial
1724:     requires: p4est
1725:     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

1727:   test:
1728:     suffix: fas_newton_0_p4est
1729:     requires: p4est !single !__float128
1730:     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

1732:   # Full solve simplicial AMR
1733:   test:
1734:     suffix: tri_p1_adapt_0
1735:     requires: pragmatic
1736:     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

1738:   test:
1739:     suffix: tri_p1_adapt_1
1740:     requires: pragmatic
1741:     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

1743:   test:
1744:     suffix: tri_p1_adapt_analytic_0
1745:     requires: pragmatic
1746:     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

1748:   # Full solve tensor AMR
1749:   test:
1750:     suffix: quad_q1_adapt_0
1751:     requires: p4est
1752:     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 -dm_view
1753:     filter: grep -v DM_

1755:   test:
1756:     suffix: amr_0
1757:     nsize: 5
1758:     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

1760:   test:
1761:     suffix: amr_1
1762:     requires: p4est !complex
1763:     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

1765:   test:
1766:     suffix: p4est_solve_bddc
1767:     requires: p4est !complex
1768:     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
1769:     nsize: 4

1771:   test:
1772:     suffix: p4est_solve_fas
1773:     requires: p4est
1774:     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
1775:     nsize: 4
1776:     TODO: identical machine two runs produce slightly different solver trackers

1778:   test:
1779:     suffix: p4est_convergence_test_1
1780:     requires: p4est
1781:     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
1782:     nsize: 4

1784:   test:
1785:     suffix: p4est_convergence_test_2
1786:     requires: p4est
1787:     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

1789:   test:
1790:     suffix: p4est_convergence_test_3
1791:     requires: p4est
1792:     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

1794:   test:
1795:     suffix: p4est_convergence_test_4
1796:     requires: p4est
1797:     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
1798:     timeoutfactor: 5

1800:   # Serial tests with GLVis visualization
1801:   test:
1802:     suffix: glvis_2d_tet_p1
1803:     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
1804:   test:
1805:     suffix: glvis_2d_tet_p2
1806:     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
1807:   test:
1808:     suffix: glvis_2d_hex_p1
1809:     args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -simplex 0 -dm_refine 1
1810:   test:
1811:     suffix: glvis_2d_hex_p2
1812:     args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -simplex 0 -dm_refine 1
1813:   test:
1814:     suffix: glvis_2d_hex_p2_p4est
1815:     requires: p4est
1816:     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

1818: TEST*/