Actual source code: ex12.c

petsc-master 2018-01-17
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, check, 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:   /* Solver */
 52:   PC            pcmg;              /* This is needed for error monitoring */
 53: } AppCtx;

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

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

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

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

 73:   so that

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

 77:   For Neumann conditions, we have

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

 84:   Which we can express as

 86:     \nabla u \cdot  \hat n|_\Gamma = {2 x, 2 y} \cdot \hat n = 2 (x + y)
 87: */
 88: static PetscErrorCode quadratic_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 89: {
 90:   *u = x[0]*x[0] + x[1]*x[1];
 91:   return 0;
 92: }

 94: static void quadratic_u_field_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
 95:                                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
 96:                                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
 97:                                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
 98: {
 99:   uexact[0] = a[0];
100: }

102: static PetscErrorCode circle_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
103: {
104:   const PetscReal alpha   = 500.;
105:   const PetscReal radius2 = PetscSqr(0.15);
106:   const PetscReal r2      = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5);
107:   const PetscReal xi      = alpha*(radius2 - r2);

109:   *u = PetscTanhScalar(xi) + 1.0;
110:   return 0;
111: }

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

118:   *u = PetscSinScalar(alpha*xy) * (alpha*PetscAbsScalar(xy) < 2*PETSC_PI ? 1 : 0.01);
119:   return 0;
120: }

122: static void f0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
123:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
124:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
125:                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
126: {
127:   f0[0] = 4.0;
128: }

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

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

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

151:   f0[0] = PetscSinScalar(alpha*xy) * (alpha*PetscAbsScalar(xy) < 2*PETSC_PI ? 1 : 0.01);
152: }

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

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

172: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
173: static void f1_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
174:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
175:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
176:                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
177: {
178:   PetscInt d;
179:   for (d = 0; d < dim; ++d) f1[d] = u_x[d];
180: }

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

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

196:     u = sin(2 pi x)
197:     f = -4 pi^2 sin(2 pi x)

199:   so that

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

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

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

220:     u = sin(2 pi x) sin(2 pi y)
221:     f = -8 pi^2 sin(2 pi x)

223:   so that

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

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

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

244:     u  = x^2 + y^2
245:     f  = 6 (x + y)
246:     nu = (x + y)

248:   so that

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

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

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

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

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

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

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

308:     u  = x^2 + y^2
309:     f  = 16 (x^2 + y^2)
310:     nu = 1/2 |grad u|^2

312:   so that

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

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

336: /*
337:   grad (u + eps w) - grad u = eps grad w

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

360: /*
361:   In 3D for Dirichlet conditions we use exact solution:

363:     u = 2/3 (x^2 + y^2 + z^2)
364:     f = 4

366:   so that

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

370:   For Neumann conditions, we have

372:     -\nabla u \cdot -\hat z |_{z=0} =  (2z)|_{z=0} =  0 (bottom)
373:     -\nabla u \cdot  \hat z |_{z=1} = -(2z)|_{z=1} = -2 (top)
374:     -\nabla u \cdot -\hat y |_{y=0} =  (2y)|_{y=0} =  0 (front)
375:     -\nabla u \cdot  \hat y |_{y=1} = -(2y)|_{y=1} = -2 (back)
376:     -\nabla u \cdot -\hat x |_{x=0} =  (2x)|_{x=0} =  0 (left)
377:     -\nabla u \cdot  \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)

379:   Which we can express as

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

389: static void quadratic_u_field_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
390:                                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
391:                                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
392:                                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
393: {
394:   uexact[0] = a[0];
395: }

397: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
398: {
399:   const char    *bcTypes[3]  = {"neumann", "dirichlet", "none"};
400:   const char    *runTypes[4] = {"full", "exact", "test", "perf"};
401:   const char    *coeffTypes[6] = {"none", "analytic", "field", "nonlinear", "circle", "cross"};
402:   PetscInt       bd, bc, run, coeff, n;
403:   PetscBool      flg;

407:   options->debug               = 0;
408:   options->runType             = RUN_FULL;
409:   options->dim                 = 2;
410:   options->periodicity[0]      = DM_BOUNDARY_NONE;
411:   options->periodicity[1]      = DM_BOUNDARY_NONE;
412:   options->periodicity[2]      = DM_BOUNDARY_NONE;
413:   options->cells[0]            = 2;
414:   options->cells[1]            = 2;
415:   options->cells[2]            = 2;
416:   options->filename[0]         = '\0';
417:   options->interpolate         = PETSC_FALSE;
418:   options->refinementLimit     = 0.0;
419:   options->bcType              = DIRICHLET;
420:   options->variableCoefficient = COEFF_NONE;
421:   options->fieldBC             = PETSC_FALSE;
422:   options->jacobianMF          = PETSC_FALSE;
423:   options->showInitial         = PETSC_FALSE;
424:   options->showSolution        = PETSC_FALSE;
425:   options->restart             = PETSC_FALSE;
426:   options->check               = PETSC_FALSE;
427:   options->viewHierarchy       = PETSC_FALSE;
428:   options->simplex             = PETSC_TRUE;
429:   options->quiet               = PETSC_FALSE;
430:   options->nonzInit            = PETSC_FALSE;

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

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

439:   PetscOptionsInt("-dim", "The topological mesh dimension", "ex12.c", options->dim, &options->dim, NULL);
440:   bd = options->periodicity[0];
441:   PetscOptionsEList("-x_periodicity", "The x-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[0]], &bd, NULL);
442:   options->periodicity[0] = (DMBoundaryType) bd;
443:   bd = options->periodicity[1];
444:   PetscOptionsEList("-y_periodicity", "The y-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[1]], &bd, NULL);
445:   options->periodicity[1] = (DMBoundaryType) bd;
446:   bd = options->periodicity[2];
447:   PetscOptionsEList("-z_periodicity", "The z-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[2]], &bd, NULL);
448:   options->periodicity[2] = (DMBoundaryType) bd;
449:   n = 3;
450:   PetscOptionsIntArray("-cells", "The initial mesh division", "ex12.c", options->cells, &n, NULL);
451:   PetscOptionsString("-f", "Mesh filename to read", "ex12.c", options->filename, options->filename, sizeof(options->filename), &flg);
452:   PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex12.c", options->interpolate, &options->interpolate, NULL);
453:   PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex12.c", options->refinementLimit, &options->refinementLimit, NULL);
454:   bc   = options->bcType;
455:   PetscOptionsEList("-bc_type","Type of boundary condition","ex12.c",bcTypes,3,bcTypes[options->bcType],&bc,NULL);
456:   options->bcType = (BCType) bc;
457:   coeff = options->variableCoefficient;
458:   PetscOptionsEList("-variable_coefficient","Type of variable coefficent","ex12.c",coeffTypes,6,coeffTypes[options->variableCoefficient],&coeff,NULL);
459:   options->variableCoefficient = (CoeffType) coeff;

461:   PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL);
462:   PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL);
463:   PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL);
464:   PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL);
465:   PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL);
466:   PetscOptionsBool("-check", "Compare with default integration routines", "ex12.c", options->check, &options->check, NULL);
467:   PetscOptionsBool("-dm_view_hierarchy", "View the coarsened hierarchy", "ex12.c", options->viewHierarchy, &options->viewHierarchy, NULL);
468:   PetscOptionsBool("-simplex", "Simplicial (true) or tensor (false) mesh", "ex12.c", options->simplex, &options->simplex, NULL);
469:   PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL);
470:   PetscOptionsBool("-nonzero_initial_guess", "nonzero intial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL);
471:   PetscOptionsEnd();
472:   PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
473:   return(0);
474: }

476: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
477: {
478:   DMLabel        label;

482:   DMCreateLabel(dm, name);
483:   DMGetLabel(dm, name, &label);
484:   DMPlexMarkBoundaryFaces(dm, label);
485:   DMPlexLabelComplete(dm, label);
486:   return(0);
487: }

489: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
490: {
491:   PetscInt       dim             = user->dim;
492:   const char    *filename        = user->filename;
493:   PetscBool      interpolate     = user->interpolate;
494:   PetscReal      refinementLimit = user->refinementLimit;
495:   size_t         len;

499:   PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
500:   PetscStrlen(filename, &len);
501:   if (!len) {
502:     PetscInt d;

504:     if (user->periodicity[0] || user->periodicity[1] || user->periodicity[2]) for (d = 0; d < dim; ++d) user->cells[d] = PetscMax(user->cells[d], 3);
505:     DMPlexCreateBoxMesh(comm, dim, user->simplex, user->cells, NULL, NULL, user->periodicity, interpolate, dm);
506:     PetscObjectSetName((PetscObject) *dm, "Mesh");
507:   } else {
508:     DMPlexCreateFromFile(comm, filename, interpolate, dm);
509:     DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
510:   }
511:   {
512:     PetscPartitioner part;
513:     DM               refinedMesh     = NULL;
514:     DM               distributedMesh = NULL;

516:     /* Refine mesh using a volume constraint */
517:     if (refinementLimit > 0.0) {
518:       DMPlexSetRefinementLimit(*dm, refinementLimit);
519:       DMRefine(*dm, comm, &refinedMesh);
520:       if (refinedMesh) {
521:         const char *name;

523:         PetscObjectGetName((PetscObject) *dm,         &name);
524:         PetscObjectSetName((PetscObject) refinedMesh,  name);
525:         DMDestroy(dm);
526:         *dm  = refinedMesh;
527:       }
528:     }
529:     /* Distribute mesh over processes */
530:     DMPlexGetPartitioner(*dm,&part);
531:     PetscPartitionerSetFromOptions(part);
532:     DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
533:     if (distributedMesh) {
534:       DMDestroy(dm);
535:       *dm  = distributedMesh;
536:     }
537:   }
538:   if (user->bcType == NEUMANN) {
539:     DMLabel   label;

541:     DMCreateLabel(*dm, "boundary");
542:     DMGetLabel(*dm, "boundary", &label);
543:     DMPlexMarkBoundaryFaces(*dm, label);
544:   } else if (user->bcType == DIRICHLET) {
545:     PetscBool hasLabel;

547:     DMHasLabel(*dm,"marker",&hasLabel);
548:     if (!hasLabel) {CreateBCLabel(*dm, "marker");}
549:   }
550:   {
551:     char      convType[256];
552:     PetscBool flg;

554:     PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
555:     PetscOptionsFList("-dm_plex_convert_type","Convert DMPlex to another format","ex12",DMList,DMPLEX,convType,256,&flg);
556:     PetscOptionsEnd();
557:     if (flg) {
558:       DM dmConv;

560:       DMConvert(*dm,convType,&dmConv);
561:       if (dmConv) {
562:         DMDestroy(dm);
563:         *dm  = dmConv;
564:       }
565:     }
566:   }
567:   DMLocalizeCoordinates(*dm); /* needed for periodic */
568:   DMSetFromOptions(*dm);
569:   DMViewFromOptions(*dm, NULL, "-dm_view");
570:   if (user->viewHierarchy) {
571:     DM       cdm = *dm;
572:     PetscInt i   = 0;
573:     char     buf[256];

575:     while (cdm) {
576:       DMSetUp(cdm);
577:       DMGetCoarseDM(cdm, &cdm);
578:       ++i;
579:     }
580:     cdm = *dm;
581:     while (cdm) {
582:       PetscViewer       viewer;
583:       PetscBool   isHDF5, isVTK;

585:       --i;
586:       PetscViewerCreate(comm,&viewer);
587:       PetscViewerSetType(viewer,PETSCVIEWERHDF5);
588:       PetscViewerSetOptionsPrefix(viewer,"hierarchy_");
589:       PetscViewerSetFromOptions(viewer);
590:       PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&isHDF5);
591:       PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERVTK,&isVTK);
592:       if (isHDF5) {
593:         PetscSNPrintf(buf, 256, "ex12-%d.h5", i);
594:       }
595:       else if (isVTK) {
596:         PetscSNPrintf(buf, 256, "ex12-%d.vtu", i);
597:         PetscViewerPushFormat(viewer,PETSC_VIEWER_VTK_VTU);
598:       }
599:       else {
600:         PetscSNPrintf(buf, 256, "ex12-%d", i);
601:       }
602:       PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
603:       PetscViewerFileSetName(viewer,buf);
604:       DMView(cdm, viewer);
605:       PetscViewerDestroy(&viewer);
606:       DMGetCoarseDM(cdm, &cdm);
607:     }
608:   }
609:   PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
610:   return(0);
611: }

613: static PetscErrorCode SetupProblem(PetscDS prob, AppCtx *user)
614: {
615:   const PetscInt id = 1;

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

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

699:   DMCreateLocalVector(dmAux, &nu);
700:   DMProjectFunctionLocal(dmAux, 0.0, matFuncs, NULL, INSERT_ALL_VALUES, nu);
701:   PetscObjectCompose((PetscObject) dm, "A", (PetscObject) nu);
702:   VecDestroy(&nu);
703:   return(0);
704: }

706: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
707: {
708:   PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
709:   Vec            uexact;
710:   PetscInt       dim;

714:   DMGetDimension(dm, &dim);
715:   if (dim == 2) bcFuncs[0] = quadratic_u_2d;
716:   else          bcFuncs[0] = quadratic_u_3d;
717:   DMCreateLocalVector(dmAux, &uexact);
718:   DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
719:   PetscObjectCompose((PetscObject) dm, "A", (PetscObject) uexact);
720:   VecDestroy(&uexact);
721:   return(0);
722: }

724: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
725: {
726:   DM             cdm   = dm;
727:   const PetscInt dim   = user->dim;
728:   PetscFE        feAux = NULL;
729:   PetscFE        feCh  = NULL;
730:   PetscFE        fe;
731:   PetscDS        prob, probAux = NULL, probCh = NULL;
732:   PetscBool      simplex = user->simplex;

736:   /* Create finite element */
737:   PetscFECreateDefault(dm, dim, 1, simplex, NULL, -1, &fe);
738:   PetscObjectSetName((PetscObject) fe, "potential");
739:   if (user->variableCoefficient == COEFF_FIELD) {
740:     PetscQuadrature q;

742:     PetscFECreateDefault(dm, dim, 1, simplex, "mat_", -1, &feAux);
743:     PetscFEGetQuadrature(fe, &q);
744:     PetscFESetQuadrature(feAux, q);
745:     PetscDSCreate(PetscObjectComm((PetscObject)dm),&probAux);
746:     PetscDSSetDiscretization(probAux, 0, (PetscObject) feAux);
747:   } else if (user->fieldBC) {
748:     PetscQuadrature q;

750:     PetscFECreateDefault(dm, dim, 1, simplex, "bc_", -1, &feAux);
751:     PetscFEGetQuadrature(fe, &q);
752:     PetscFESetQuadrature(feAux, q);
753:     PetscDSCreate(PetscObjectComm((PetscObject)dm),&probAux);
754:     PetscDSSetDiscretization(probAux, 0, (PetscObject) feAux);
755:   }
756:   if (user->check) {
757:     PetscFECreateDefault(dm, dim, 1, simplex, "ch_", -1, &feCh);
758:     PetscDSCreate(PetscObjectComm((PetscObject)dm),&probCh);
759:     PetscDSSetDiscretization(probCh, 0, (PetscObject) feCh);
760:   }
761:   /* Set discretization and boundary conditions for each mesh */
762:   DMGetDS(dm, &prob);
763:   PetscDSSetDiscretization(prob, 0, (PetscObject) fe);
764:   SetupProblem(prob, user);
765:   while (cdm) {
766:     DM coordDM;

768:     DMSetDS(cdm,prob);
769:     DMGetCoordinateDM(cdm,&coordDM);
770:     if (feAux) {
771:       DM      dmAux;

773:       DMClone(cdm, &dmAux);
774:       DMSetCoordinateDM(dmAux, coordDM);
775:       DMSetDS(dmAux, probAux);
776:       PetscObjectCompose((PetscObject) dm, "dmAux", (PetscObject) dmAux);
777:       if (user->fieldBC) {SetupBC(cdm, dmAux, user);}
778:       else               {SetupMaterial(cdm, dmAux, user);}
779:       DMDestroy(&dmAux);
780:     }
781:     if (feCh) {
782:       DM      dmCh;

784:       DMClone(cdm, &dmCh);
785:       DMSetCoordinateDM(dmCh, coordDM);
786:       DMSetDS(dmCh, probCh);
787:       PetscObjectCompose((PetscObject) dm, "dmCh", (PetscObject) dmCh);
788:       DMDestroy(&dmCh);
789:     }
790:     if (user->bcType == DIRICHLET) {
791:       PetscBool hasLabel;

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

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

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

811:   Collective on KSP

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

819:   Level: intermediate

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

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

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

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

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

894:   Collective on SNES

896:   Input Parameters:
897: + snes  - the SNES
898: . its   - iteration number
899: . rnorm - 2-norm of residual
900: - ctx   - user context

902:   Level: intermediate

904: .keywords: SNES, nonlinear, default, monitor, norm
905: .seealso: SNESMonitorDefault(), SNESMonitorSet(), SNESMonitorSolution()
906: @*/
907: static PetscErrorCode SNESMonitorError(SNES snes, PetscInt its, PetscReal rnorm, void *ctx)
908: {
909:   AppCtx        *user = (AppCtx *) ctx;
910:   DM             dm;
911:   Vec            u, r;
912:   PetscInt       level = -1;
913:   PetscBool      hasLevel;
914: #if defined(PETSC_HAVE_HDF5)
915:   PetscViewer    viewer;
916: #endif
917:   char           buf[256];

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

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

956:   PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
957:   ProcessOptions(PETSC_COMM_WORLD, &user);
958:   SNESCreate(PETSC_COMM_WORLD, &snes);
959:   CreateMesh(PETSC_COMM_WORLD, &user, &dm);
960:   SNESSetDM(snes, dm);
961:   DMSetApplicationContext(dm, &user);

963:   PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
964:   SetupDiscretization(dm, &user);

966:   DMCreateGlobalVector(dm, &u);
967:   PetscObjectSetName((PetscObject) u, "potential");

969:   DMCreateMatrix(dm, &J);
970:   if (user.jacobianMF) {
971:     PetscInt M, m, N, n;

973:     MatGetSize(J, &M, &N);
974:     MatGetLocalSize(J, &m, &n);
975:     MatCreate(PETSC_COMM_WORLD, &A);
976:     MatSetSizes(A, m, n, M, N);
977:     MatSetType(A, MATSHELL);
978:     MatSetUp(A);
979: #if 0
980:     MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction);
981: #endif

983:     userJ.dm   = dm;
984:     userJ.J    = J;
985:     userJ.user = &user;

987:     DMCreateLocalVector(dm, &userJ.u);
988:     if (user.fieldBC) {DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);}
989:     else              {DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);}
990:     MatShellSetContext(A, &userJ);
991:   } else {
992:     A = J;
993:   }
994:   if (user.bcType == NEUMANN) {
995:     MatNullSpaceCreate(PetscObjectComm((PetscObject) dm), PETSC_TRUE, 0, NULL, &nullSpace);
996:     MatSetNullSpace(A, nullSpace);
997:   }

999:   DMPlexSetSNESLocalFEM(dm,&user,&user,&user);
1000:   SNESSetJacobian(snes, A, J, NULL, NULL);

1002:   SNESSetFromOptions(snes);

1004:   if (user.fieldBC) {DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u);}
1005:   else              {DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);}
1006:   if (user.restart) {
1007: #if defined(PETSC_HAVE_HDF5)
1008:     PetscViewer viewer;

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

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

1059:     if (user.nonzInit) initialGuess[0] = ecks;
1060:     if (user.runType == RUN_FULL) {
1061:       DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
1062:     }
1063:     if (user.debug) {
1064:       PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1065:       VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1066:     }
1067:     SNESSolve(snes, NULL, u);
1068:     SNESGetSolution(snes, &u);
1069:     SNESGetDM(snes, &dm);

1071:     if (user.showSolution) {
1072:       PetscPrintf(PETSC_COMM_WORLD, "Solution\n");
1073:       VecChop(u, 3.0e-9);
1074:       VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1075:     }
1076:     VecViewFromOptions(u, NULL, "-vec_view");
1077:   } else if (user.runType == RUN_PERF) {
1078:     Vec       r;
1079:     PetscReal res = 0.0;

1081:     SNESGetFunction(snes, &r, NULL, NULL);
1082:     SNESComputeFunction(snes, u, r);
1083:     PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1084:     VecChop(r, 1.0e-10);
1085:     VecNorm(r, NORM_2, &res);
1086:     PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1087:   } else {
1088:     Vec       r;
1089:     PetscReal res = 0.0;

1091:     /* Check discretization error */
1092:     SNESGetFunction(snes, &r, NULL, NULL);
1093:     PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1094:     if (!user.quiet) {VecView(u, PETSC_VIEWER_STDOUT_WORLD);}
1095:     DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
1096:     if (error < 1.0e-11) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < 1.0e-11\n");}
1097:     else                 {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", error);}
1098:     /* Check residual */
1099:     SNESComputeFunction(snes, u, r);
1100:     PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1101:     VecChop(r, 1.0e-10);
1102:     if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1103:     VecNorm(r, NORM_2, &res);
1104:     PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1105:     /* Check Jacobian */
1106:     {
1107:       Vec b;

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

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

1137: /*TEST
1138:   build:
1139:     requires: hdf5 triangle
1140:   # 2D serial P1 test 0-4
1141:   test:
1142:     suffix: 0
1143:     requires: hdf5 triangle
1144:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 0 -petscspace_order 1 -show_initial -dm_plex_print_fem 1

1146:   test:
1147:     suffix: 1
1148:     requires: hdf5 triangle
1149:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1

1151:   test:
1152:     suffix: 2
1153:     requires: hdf5 triangle
1154:     args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1

1156:   test:
1157:     suffix: 3
1158:     requires: hdf5 triangle
1159:     args: -run_type test -refinement_limit 0.0    -bc_type neumann   -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail

1161:   test:
1162:     suffix: 4
1163:     requires: hdf5 triangle
1164:     args: -run_type test -refinement_limit 0.0625 -bc_type neumann   -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1

1166:   # 2D serial P2 test 5-8
1167:   test:
1168:     suffix: 5
1169:     requires: hdf5 triangle
1170:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1

1172:   test:
1173:     suffix: 6
1174:     requires: hdf5 triangle
1175:     args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1

1177:   test:
1178:     suffix: 7
1179:     requires: hdf5 triangle
1180:     args: -run_type test -refinement_limit 0.0    -bc_type neumann   -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail

1182:   test:
1183:     suffix: 8
1184:     requires: hdf5 triangle
1185:     args: -run_type test -refinement_limit 0.0625 -bc_type neumann   -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail

1187:   # 3D serial P1 test 9-12
1188:   test:
1189:     suffix: 9
1190:     requires: hdf5 ctetgen
1191:     args: -run_type test -dim 3 -refinement_limit 0.0    -bc_type dirichlet -interpolate 0 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1

1193:   test:
1194:     suffix: 10
1195:     requires: hdf5 ctetgen
1196:     args: -run_type test -dim 3 -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1

1198:   test:
1199:     suffix: 11
1200:     requires: hdf5 ctetgen
1201:     args: -run_type test -dim 3 -refinement_limit 0.0125 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1

1203:   test:
1204:     suffix: 12
1205:     requires: hdf5 ctetgen
1206:     args: -run_type test -dim 3 -refinement_limit 0.0    -bc_type neumann   -interpolate 1 -petscspace_order 1 -snes_fd -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1

1208:   # Analytic variable coefficient 13-20
1209:   test:
1210:     suffix: 13
1211:     requires: hdf5 triangle
1212:     args: -run_type test -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1213:   test:
1214:     suffix: 14
1215:     requires: hdf5 triangle
1216:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1217:   test:
1218:     suffix: 15
1219:     requires: hdf5 triangle
1220:     args: -run_type test -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1221:   test:
1222:     suffix: 16
1223:     requires: hdf5 triangle
1224:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1225:   test:
1226:     suffix: 17
1227:     requires: hdf5 ctetgen
1228:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1230:   test:
1231:     suffix: 18
1232:     requires: hdf5 ctetgen
1233:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1235:   test:
1236:     suffix: 19
1237:     requires: hdf5 ctetgen
1238:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1240:   test:
1241:     suffix: 20
1242:     requires: hdf5 ctetgen
1243:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1245:   # P1 variable coefficient 21-28
1246:   test:
1247:     suffix: 21
1248:     requires: hdf5 triangle
1249:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1

1251:   test:
1252:     suffix: 22
1253:     requires: hdf5 triangle
1254:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1

1256:   test:
1257:     suffix: 23
1258:     requires: hdf5 triangle
1259:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1

1261:   test:
1262:     suffix: 24
1263:     requires: hdf5 triangle
1264:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1

1266:   test:
1267:     suffix: 25
1268:     requires: hdf5 ctetgen
1269:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1271:   test:
1272:     suffix: 26
1273:     requires: hdf5 ctetgen
1274:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1276:   test:
1277:     suffix: 27
1278:     requires: hdf5 ctetgen
1279:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1281:   test:
1282:     suffix: 28
1283:     requires: hdf5 ctetgen
1284:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1286:   # P0 variable coefficient 29-36
1287:   test:
1288:     suffix: 29
1289:     requires: hdf5 triangle
1290:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1

1292:   test:
1293:     suffix: 30
1294:     requires: hdf5 triangle
1295:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1

1297:   test:
1298:     suffix: 31
1299:     requires: hdf5 triangle
1300:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1

1302:   test:
1303:     requires: hdf5 triangle
1304:     suffix: 32
1305:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1

1307:   test:
1308:     requires: hdf5 ctetgen
1309:     suffix: 33
1310:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1312:   test:
1313:     suffix: 34
1314:     requires: hdf5 ctetgen
1315:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1317:   test:
1318:     suffix: 35
1319:     requires: hdf5 ctetgen
1320:     args: -run_type test -dim 3 -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1322:   test:
1323:     suffix: 36
1324:     requires: hdf5 ctetgen
1325:     args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field    -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1327:   # Full solve 39-44
1328:   test:
1329:     suffix: 39
1330:     requires: hdf5 triangle !single
1331:     args: -run_type full -refinement_limit 0.015625 -interpolate 1 -petscspace_order 2 -pc_type gamg -ksp_rtol 1.0e-10 -ksp_monitor_short -ksp_converged_reason -snes_monitor_short -snes_converged_reason ::ascii_info_detail

1333:   test:
1334:     suffix: 40
1335:     requires: hdf5 triangle !single
1336:     args: -run_type full -refinement_limit 0.015625 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 2 -pc_type svd -ksp_rtol 1.0e-10 -snes_monitor_short -snes_converged_reason ::ascii_info_detail

1338:   test:
1339:     suffix: 41
1340:     requires: hdf5 triangle !single
1341:     args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1343:   test:
1344:     suffix: 42
1345:     requires: hdf5 triangle !single
1346:     args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 3 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1348:   test:
1349:     suffix: 43
1350:     requires: hdf5 triangle !single
1351:     nsize: 2
1352:     args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1354:   test:
1355:     suffix: 44
1356:     requires: hdf5 triangle !single
1357:     nsize: 2
1358:     args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 3 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1360:   # Restarting
1361:   testset:
1362:     suffix: restart
1363:     requires: hdf5 triangle !complex
1364:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_order 1
1365:     test:
1366:       args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1367:     test:
1368:       args: -f sol.h5 -restart

1370:   # Periodicity
1371:   test:
1372:     suffix: periodic_0
1373:     requires: hdf5 triangle
1374:     args: -run_type full -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_order 1 -snes_converged_reason ::ascii_info_detail

1376:   # 2D serial P1 test with field bc
1377:   test:
1378:     suffix: field_bc_p1_0
1379:     requires: hdf5 triangle
1380:     args: -run_type test              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1

1382:   test:
1383:     suffix: field_bc_p1_1
1384:     requires: hdf5 triangle
1385:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1

1387:   test:
1388:     suffix: field_bc_p1_2
1389:     requires: hdf5 triangle
1390:     args: -run_type test              -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1

1392:   test:
1393:     suffix: field_bc_p1_3
1394:     requires: hdf5 triangle
1395:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1

1397:   # 3D serial P1 test with field bc
1398:   test:
1399:     suffix: field_bc_p1_4
1400:     requires: hdf5 ctetgen
1401:     args: -run_type test -dim 3              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1403:   test:
1404:     suffix: field_bc_p1_5
1405:     requires: hdf5 ctetgen
1406:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1408:   test:
1409:     suffix: field_bc_p1_6
1410:     requires: hdf5 ctetgen
1411:     args: -run_type test -dim 3              -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1413:   test:
1414:     suffix: field_bc_p1_7
1415:     requires: hdf5 ctetgen
1416:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1418:   # 2D serial P2 test with field bc
1419:   test:
1420:     suffix: field_bc_p2_0
1421:     requires: hdf5 triangle
1422:     args: -run_type test              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1

1424:   test:
1425:     suffix: field_bc_p2_1
1426:     requires: hdf5 triangle
1427:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1

1429:   test:
1430:     suffix: field_bc_p2_2
1431:     requires: hdf5 triangle
1432:     args: -run_type test              -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1

1434:   test:
1435:     suffix: field_bc_p2_3
1436:     requires: hdf5 triangle
1437:     args: -run_type test -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1

1439:   # 3D serial P2 test with field bc
1440:   test:
1441:     suffix: field_bc_p2_4
1442:     requires: hdf5 ctetgen
1443:     args: -run_type test -dim 3              -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1445:   test:
1446:     suffix: field_bc_p2_5
1447:     requires: hdf5 ctetgen
1448:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1450:   test:
1451:     suffix: field_bc_p2_6
1452:     requires: hdf5 ctetgen
1453:     args: -run_type test -dim 3              -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1455:   test:
1456:     suffix: field_bc_p2_7
1457:     requires: hdf5 ctetgen
1458:     args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type neumann   -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1

1460:   # Full solve simplex: Convergence
1461:   test:
1462:     suffix: tet_conv_p1_r0
1463:     requires: hdf5 ctetgen
1464:     args: -run_type full -dim 3 -dm_refine 0 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1

1466:   test:
1467:     suffix: tet_conv_p1_r2
1468:     requires: hdf5 ctetgen
1469:     args: -run_type full -dim 3 -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1

1471:   test:
1472:     suffix: tet_conv_p1_r3
1473:     requires: hdf5 ctetgen
1474:     args: -run_type full -dim 3 -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1

1476:   test:
1477:     suffix: tet_conv_p2_r0
1478:     requires: hdf5 ctetgen
1479:     args: -run_type full -dim 3 -dm_refine 0 -bc_type dirichlet -interpolate 1 -petscspace_order 2 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1

1481:   test:
1482:     suffix: tet_conv_p2_r2
1483:     requires: hdf5 ctetgen
1484:     args: -run_type full -dim 3 -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_order 2 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1

1486:   # Full solve simplex: BDDC
1487:   test:
1488:     suffix: tri_bddc
1489:     requires: hdf5 triangle !single
1490:     nsize: 5
1491:     args: -run_type full -petscpartitioner_type simple -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -dm_mat_type is -pc_type bddc -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0

1493:   # Full solve simplex: ASM
1494:   test:
1495:     suffix: tri_q2q1_asm_lu
1496:     requires: hdf5 triangle !single
1497:     args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type asm -pc_asm_type restrict -pc_asm_blocks 4 -sub_pc_type lu -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0

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

1504:   test:
1505:     suffix: tri_q2q1_asm_sor
1506:     requires: hdf5 triangle !single
1507:     args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type asm -pc_asm_type restrict -pc_asm_blocks 4 -sub_pc_type sor -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0

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

1514:   # Full solve simplex: FAS
1515:   test:
1516:     suffix: fas_newton_0
1517:     requires: hdf5 triangle !single
1518:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1520:   test:
1521:     suffix: fas_newton_1
1522:     requires: hdf5 triangle !single
1523:     args: -run_type full -dm_refine_hierarchy 3 -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 3 -ksp_rtol 1.0e-10 -fas_coarse_pc_type lu -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1525:   test:
1526:     suffix: fas_ngs_0
1527:     requires: hdf5 triangle !single
1528:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -dm_refine_hierarchy 1 -snes_view -fas_levels_1_snes_type ngs -fas_levels_1_snes_monitor_short

1530:   test:
1531:     suffix: fas_newton_coarse_0
1532:     requires: hdf5 pragmatic triangle
1533:     TODO: broken
1534:     args: -run_type full -dm_refine 2 -dm_plex_hash_location -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1536:   test:
1537:     suffix: mg_newton_coarse_0
1538:     requires: hdf5 triangle pragmatic
1539:     args: -run_type full -dm_refine 3 -interpolate 1 -petscspace_order 1 -snes_monitor_short -ksp_monitor_true_residual -snes_linesearch_type basic -snes_converged_reason ::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

1541:   test:
1542:     suffix: mg_newton_coarse_1
1543:     requires: hdf5 triangle pragmatic
1544:     args: -run_type full -dm_refine 5 -interpolate 1 -petscspace_order 1 -dm_coarsen_hierarchy 5 -dm_plex_hash_location -dm_plex_separate_marker -dm_plex_coarsen_bd_label marker -dm_plex_remesh_bd -ksp_type richardson -ksp_rtol 1.0e-12 -pc_type mg -pc_mg_levels 3 -mg_levels_ksp_max_it 2 -snes_converged_reason ::ascii_info_detail -snes_monitor -ksp_monitor_true_residual -mg_levels_ksp_monitor_true_residual -dm_view -ksp_view

1546:   test:
1547:     suffix: mg_newton_coarse_2
1548:     requires: hdf5 triangle pragmatic
1549:     args: -run_type full -dm_refine 5 -interpolate 1 -petscspace_order 1 -dm_coarsen_hierarchy 5 -dm_plex_hash_location -dm_plex_separate_marker -dm_plex_remesh_bd -ksp_type richardson -ksp_rtol 1.0e-12 -pc_type mg -pc_mg_levels 3 -mg_levels_ksp_max_it 2 -snes_converged_reason ::ascii_info_detail -snes_monitor -ksp_monitor_true_residual -mg_levels_ksp_monitor_true_residual -dm_view -ksp_view

1551:   # Full solve tensor
1552:   test:
1553:     suffix: tensor_plex_2d
1554:     requires: hdf5
1555:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 1 -dm_refine_hierarchy 2 -cells 2,2

1557:   test:
1558:     suffix: tensor_p4est_2d
1559:     requires: hdf5 p4est
1560:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 1 -dm_forest_initial_refinement 2 -dm_forest_minimum_refinement 0 -dm_plex_convert_type p4est -cells 2,2

1562:   test:
1563:     suffix: tensor_plex_3d
1564:     requires: hdf5
1565:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 1 -dim 3 -dm_refine_hierarchy 1 -cells 2,2,2

1567:   test:
1568:     suffix: tensor_p4est_3d
1569:     requires: hdf5 p4est
1570:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 1 -dm_forest_initial_refinement 1 -dm_forest_minimum_refinement 0 -dim 3 -dm_plex_convert_type p8est -cells 2,2,2

1572:   test:
1573:     suffix: p4est_test_q2_conformal_serial
1574:     requires: hdf5 p4est
1575:     args: -run_type test -interpolate 1 -petscspace_order 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -cells 2,2

1577:   test:
1578:     suffix: p4est_test_q2_conformal_parallel
1579:     requires: hdf5 p4est
1580:     nsize: 7
1581:     args: -run_type test -interpolate 1 -petscspace_order 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type simple -cells 2,2

1583:   test:
1584:     suffix: p4est_test_q2_nonconformal_serial
1585:     requires: hdf5 p4est
1586:     args: -run_type test -interpolate 1 -petscspace_order 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2

1588:   test:
1589:     suffix: p4est_test_q2_nonconformal_parallel
1590:     requires: hdf5 p4est
1591:     nsize: 7
1592:     args: -run_type test -interpolate 1 -petscspace_order 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -cells 2,2

1594:   test:
1595:     suffix: p4est_exact_q2_conformal_serial
1596:     requires: hdf5 p4est !single
1597:     args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1599:   test:
1600:     suffix: p4est_exact_q2_conformal_parallel
1601:     requires: hdf5 p4est !single
1602:     nsize: 7
1603:     args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1605:   test:
1606:     suffix: p4est_exact_q2_nonconformal_serial
1607:     requires: hdf5 p4est
1608:     args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1610:   test:
1611:     suffix: p4est_exact_q2_nonconformal_parallel
1612:     requires: hdf5 p4est
1613:     nsize: 7
1614:     args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1616:   test:
1617:     suffix: p4est_full_q2_nonconformal_serial
1618:     requires: hdf5 p4est !single
1619:     filter: grep -v "variant HERMITIAN"
1620:     args: -run_type full -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type cg -fas_coarse_pc_type jacobi -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1622:   test:
1623:     suffix: p4est_full_q2_nonconformal_parallel
1624:     requires: p4est hdf5 !single
1625:     filter: grep -v "variant HERMITIAN"
1626:     nsize: 7
1627:     args: -run_type full -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type cg -fas_coarse_pc_type jacobi -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1629:   test:
1630:     suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1631:     requires: hdf5 p4est
1632:     filter: grep -v "variant HERMITIAN"
1633:     nsize: 7
1634:     args: -run_type full -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -dm_mat_type is -pc_type bddc -ksp_type cg -fas_coarse_pc_type bddc -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1636:   test:
1637:     suffix: p4est_full_q2_nonconformal_parallel_bddc
1638:     requires: hdf5 p4est
1639:     filter: grep -v "variant HERMITIAN"
1640:     nsize: 7
1641:     args: -run_type full -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type newtonls -dm_mat_type is -pc_type bddc -ksp_type cg -snes_monitor_short -snes_linesearch_type basic -snes_converged_reason ::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

1643:   test:
1644:     suffix: p4est_fas_q2_conformal_serial
1645:     requires: hdf5 p4est
1646:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type gmres -fas_coarse_pc_type svd -fas_coarse_ksp_type gmres -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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
1647:     TODO: broken

1649:   test:
1650:     suffix: p4est_fas_q2_nonconformal_serial
1651:     requires: hdf5 p4est broken
1652:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type gmres -fas_coarse_pc_type jacobi -fas_coarse_ksp_type gmres -fas_coarse_ksp_monitor_true_residual -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1654:   test:
1655:     suffix: fas_newton_0_p4est
1656:     requires: hdf5 p4est !single
1657:     args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason ::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

1659:   # Full solve simplicial AMR
1660:   test:
1661:     suffix: tri_p1_adapt_0
1662:     requires: pragmatic
1663:     args: -run_type exact -dim 2 -dm_refine 5 -bc_type dirichlet -interpolate 1 -petscspace_order 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

1665:   test:
1666:     suffix: tri_p1_adapt_1
1667:     requires: pragmatic
1668:     args: -run_type exact -dim 2 -dm_refine 5 -bc_type dirichlet -interpolate 1 -petscspace_order 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

1670:   test:
1671:     suffix: tri_p1_adapt_analytic_0
1672:     requires: pragmatic
1673:     args: -run_type exact -dim 2 -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -adaptor_monitor -dm_view -dm_adapt_iter_view

1675:   # Full solve tensor AMR
1676:   test:
1677:     suffix: quad_q1_adapt_0
1678:     requires: p4est
1679:     args: -run_type exact -dim 2 -simplex 0 -dm_plex_convert_type p4est -bc_type dirichlet -interpolate 1 -petscspace_order 1 -variable_coefficient circle -snes_converged_reason ::ascii_info_detail -pc_type lu -dm_forest_initial_refinement 4   -snes_adapt_initial -dm_view
1680:     filter: grep -v DM_

1682:   test:
1683:     suffix: amr_0
1684:     requires: hdf5
1685:     nsize: 5
1686:     args: -run_type test -petscpartitioner_type simple -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 1 -dm_refine 1 -cells 2,2

1688:   test:
1689:     suffix: amr_1
1690:     requires: hdf5 p4est !complex
1691:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_order 1 -dm_plex_convert_type p4est -dm_p4est_refine_pattern center -dm_forest_maximum_refinement 5 -dm_view vtk:amr.vtu:vtk_vtu -vec_view vtk:amr.vtu:vtk_vtu:append -cells 2,2

1693:   test:
1694:     suffix: p4est_solve_bddc
1695:     requires: p4est
1696:     args: -run_type full -variable_coefficient nonlinear -nonzero_initial_guess 1 -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type newtonls -dm_mat_type is -pc_type bddc -ksp_type cg -snes_monitor_short -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
1697:     nsize: 4

1699:   test:
1700:     suffix: p4est_solve_fas
1701:     requires: p4est
1702:     args: -run_type full -variable_coefficient nonlinear -nonzero_initial_guess 1 -interpolate 1 -petscspace_order 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
1703:     nsize: 4
1704:     TODO: identical machine two runs produce slightly different solver trackers

1706:   test:
1707:     suffix: p4est_convergence_test_1
1708:     requires: p4est
1709:     args:  -quiet -run_type test -interpolate 1 -petscspace_order 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
1710:     nsize: 4

1712:   test:
1713:     suffix: p4est_convergence_test_2
1714:     requires: p4est
1715:     args: -quiet -run_type test -interpolate 1 -petscspace_order 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

1717:   test:
1718:     suffix: p4est_convergence_test_3
1719:     requires: p4est
1720:     args: -quiet -run_type test -interpolate 1 -petscspace_order 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

1722:   test:
1723:     suffix: p4est_convergence_test_4
1724:     requires: p4est
1725:     args: -quiet -run_type test -interpolate 1 -petscspace_order 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
1726:     timeoutfactor: 5

1728: TEST*/