Actual source code: ex12.c

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

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

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

 12: Information on refinement:

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

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

 23: typedef enum {NEUMANN, DIRICHLET, NONE} BCType;
 24: typedef enum {RUN_FULL, RUN_EXACT, RUN_TEST, RUN_PERF} RunType;
 25: typedef enum {COEFF_NONE, COEFF_ANALYTIC, COEFF_FIELD, COEFF_NONLINEAR, COEFF_CIRCLE, COEFF_CROSS} 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:   PetscBool      checkksp;          /* Whether to check the KSPSolve for runType == RUN_TEST */
 55: } AppCtx;

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

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

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

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

 75:   so that

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

 79:   For Neumann conditions, we have

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

 86:   Which we can express as

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

205:   so that

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

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

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

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

229:   so that

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

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

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

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

254:   so that

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

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

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

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

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

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

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

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

318:   so that

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

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

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

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

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

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

372:   so that

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

376:   For Neumann conditions, we have

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

385:   Which we can express as

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

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

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

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

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

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

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

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

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

494: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
495: {
496:   DM             plex;
497:   DMLabel        label;

501:   DMCreateLabel(dm, name);
502:   DMGetLabel(dm, name, &label);
503:   DMConvert(dm, DMPLEX, &plex);
504:   DMPlexMarkBoundaryFaces(plex, 1, label);
505:   DMDestroy(&plex);
506:   return(0);
507: }

509: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
510: {
511:   PetscInt       dim             = user->dim;
512:   const char    *filename        = user->filename;
513:   PetscBool      interpolate     = user->interpolate;
514:   PetscReal      refinementLimit = user->refinementLimit;
515:   size_t         len;

519:   PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
520:   PetscStrlen(filename, &len);
521:   if (!len) {
522:     PetscInt d;

524:     if (user->periodicity[0] || user->periodicity[1] || user->periodicity[2]) for (d = 0; d < dim; ++d) user->cells[d] = PetscMax(user->cells[d], 3);
525:     DMPlexCreateBoxMesh(comm, dim, user->simplex, user->cells, NULL, NULL, user->periodicity, interpolate, dm);
526:     PetscObjectSetName((PetscObject) *dm, "Mesh");
527:   } else {
528:     DMPlexCreateFromFile(comm, filename, interpolate, dm);
529:     DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
530:   }
531:   {
532:     PetscPartitioner part;
533:     DM               refinedMesh     = NULL;
534:     DM               distributedMesh = NULL;

536:     /* Refine mesh using a volume constraint */
537:     if (refinementLimit > 0.0) {
538:       DMPlexSetRefinementLimit(*dm, refinementLimit);
539:       DMRefine(*dm, comm, &refinedMesh);
540:       if (refinedMesh) {
541:         const char *name;

543:         PetscObjectGetName((PetscObject) *dm,         &name);
544:         PetscObjectSetName((PetscObject) refinedMesh,  name);
545:         DMDestroy(dm);
546:         *dm  = refinedMesh;
547:       }
548:     }
549:     /* Distribute mesh over processes */
550:     DMPlexGetPartitioner(*dm,&part);
551:     PetscPartitionerSetFromOptions(part);
552:     DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
553:     if (distributedMesh) {
554:       DMDestroy(dm);
555:       *dm  = distributedMesh;
556:     }
557:   }
558:   if (interpolate) {
559:     if (user->bcType == NEUMANN) {
560:       DMLabel   label;

562:       DMCreateLabel(*dm, "boundary");
563:       DMGetLabel(*dm, "boundary", &label);
564:       DMPlexMarkBoundaryFaces(*dm, 1, label);
565:     } else if (user->bcType == DIRICHLET) {
566:       PetscBool hasLabel;

568:       DMHasLabel(*dm,"marker",&hasLabel);
569:       if (!hasLabel) {CreateBCLabel(*dm, "marker");}
570:     }
571:   }
572:   {
573:     char      convType[256];
574:     PetscBool flg;

576:     PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
577:     PetscOptionsFList("-dm_plex_convert_type","Convert DMPlex to another format","ex12",DMList,DMPLEX,convType,256,&flg);
578:     PetscOptionsEnd();
579:     if (flg) {
580:       DM dmConv;

582:       DMConvert(*dm,convType,&dmConv);
583:       if (dmConv) {
584:         DMDestroy(dm);
585:         *dm  = dmConv;
586:       }
587:     }
588:   }
589:   DMLocalizeCoordinates(*dm); /* needed for periodic */
590:   DMSetFromOptions(*dm);
591:   DMViewFromOptions(*dm, NULL, "-dm_view");
592:   if (user->viewHierarchy) {
593:     DM       cdm = *dm;
594:     PetscInt i   = 0;
595:     char     buf[256];

597:     while (cdm) {
598:       DMSetUp(cdm);
599:       DMGetCoarseDM(cdm, &cdm);
600:       ++i;
601:     }
602:     cdm = *dm;
603:     while (cdm) {
604:       PetscViewer       viewer;
605:       PetscBool   isHDF5, isVTK;

607:       --i;
608:       PetscViewerCreate(comm,&viewer);
609:       PetscViewerSetType(viewer,PETSCVIEWERHDF5);
610:       PetscViewerSetOptionsPrefix(viewer,"hierarchy_");
611:       PetscViewerSetFromOptions(viewer);
612:       PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&isHDF5);
613:       PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERVTK,&isVTK);
614:       if (isHDF5) {
615:         PetscSNPrintf(buf, 256, "ex12-%d.h5", i);
616:       } else if (isVTK) {
617:         PetscSNPrintf(buf, 256, "ex12-%d.vtu", i);
618:         PetscViewerPushFormat(viewer,PETSC_VIEWER_VTK_VTU);
619:       } else {
620:         PetscSNPrintf(buf, 256, "ex12-%d", i);
621:       }
622:       PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
623:       PetscViewerFileSetName(viewer,buf);
624:       DMView(cdm, viewer);
625:       PetscViewerDestroy(&viewer);
626:       DMGetCoarseDM(cdm, &cdm);
627:     }
628:   }
629:   PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
630:   return(0);
631: }

633: static PetscErrorCode SetupProblem(DM dm, AppCtx *user)
634: {
635:   PetscDS        prob;
636:   const PetscInt id = 1;

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

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

723:   DMCreateLocalVector(dmAux, &nu);
724:   DMProjectFunctionLocal(dmAux, 0.0, matFuncs, NULL, INSERT_ALL_VALUES, nu);
725:   PetscObjectCompose((PetscObject) dm, "A", (PetscObject) nu);
726:   VecDestroy(&nu);
727:   return(0);
728: }

730: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
731: {
732:   PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
733:   Vec            uexact;
734:   PetscInt       dim;

738:   DMGetDimension(dm, &dim);
739:   if (dim == 2) bcFuncs[0] = quadratic_u_2d;
740:   else          bcFuncs[0] = quadratic_u_3d;
741:   DMCreateLocalVector(dmAux, &uexact);
742:   DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
743:   PetscObjectCompose((PetscObject) dm, "A", (PetscObject) uexact);
744:   VecDestroy(&uexact);
745:   return(0);
746: }

748: static PetscErrorCode SetupAuxDM(DM dm, PetscFE feAux, AppCtx *user)
749: {
750:   DM             dmAux, coordDM;

754:   /* MUST call DMGetCoordinateDM() in order to get p4est setup if present */
755:   DMGetCoordinateDM(dm, &coordDM);
756:   if (!feAux) return(0);
757:   DMClone(dm, &dmAux);
758:   PetscObjectCompose((PetscObject) dm, "dmAux", (PetscObject) dmAux);
759:   DMSetCoordinateDM(dmAux, coordDM);
760:   DMSetField(dmAux, 0, NULL, (PetscObject) feAux);
761:   DMCreateDS(dmAux);
762:   if (user->fieldBC) {SetupBC(dm, dmAux, user);}
763:   else               {SetupMaterial(dm, dmAux, user);}
764:   DMDestroy(&dmAux);
765:   return(0);
766: }

768: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
769: {
770:   DM             cdm = dm;
771:   const PetscInt dim = user->dim;
772:   PetscFE        fe, feAux = NULL;
773:   PetscBool      simplex   = user->simplex;
774:   MPI_Comm       comm;

778:   /* Create finite element for each field and auxiliary field */
779:   PetscObjectGetComm((PetscObject) dm, &comm);
780:   PetscFECreateDefault(comm, dim, 1, simplex, NULL, -1, &fe);
781:   PetscObjectSetName((PetscObject) fe, "potential");
782:   if (user->variableCoefficient == COEFF_FIELD) {
783:     PetscFECreateDefault(comm, dim, 1, simplex, "mat_", -1, &feAux);
784:     PetscFECopyQuadrature(fe, feAux);
785:   } else if (user->fieldBC) {
786:     PetscFECreateDefault(comm, dim, 1, simplex, "bc_", -1, &feAux);
787:     PetscFECopyQuadrature(fe, feAux);
788:   }
789:   /* Set discretization and boundary conditions for each mesh */
790:   DMSetField(dm, 0, NULL, (PetscObject) fe);
791:   DMCreateDS(dm);
792:   SetupProblem(dm, user);
793:   while (cdm) {
794:     SetupAuxDM(cdm, feAux, user);
795:     if (user->bcType == DIRICHLET && user->interpolate) {
796:       PetscBool hasLabel;

798:       DMHasLabel(cdm, "marker", &hasLabel);
799:       if (!hasLabel) {CreateBCLabel(cdm, "marker");}
800:     }
801:     DMCopyDisc(dm, cdm);
802:     DMGetCoarseDM(cdm, &cdm);
803:   }
804:   PetscFEDestroy(&fe);
805:   PetscFEDestroy(&feAux);
806:   return(0);
807: }

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

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

814:   Collective on KSP

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

822:   Level: intermediate

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

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

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

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

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

893:   Collective on SNES

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

901:   Level: intermediate

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

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

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

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

961:   PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
962:   SetupDiscretization(dm, &user);

964:   DMCreateGlobalVector(dm, &u);
965:   PetscObjectSetName((PetscObject) u, "potential");

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

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

981:     userJ.dm   = dm;
982:     userJ.J    = J;
983:     userJ.user = &user;

985:     DMCreateLocalVector(dm, &userJ.u);
986:     if (user.fieldBC) {DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);}
987:     else              {DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);}
988:     MatShellSetContext(A, &userJ);
989:   } else {
990:     A = J;
991:   }

993:   nullSpace = NULL;
994:   if (user.bcType != DIRICHLET) {
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:     VecViewFromOptions(u, NULL, "-guess_vec_view");
1068:     SNESSolve(snes, NULL, u);
1069:     SNESGetSolution(snes, &u);
1070:     SNESGetDM(snes, &dm);

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

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

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

1110:       SNESComputeJacobian(snes, u, A, A);
1111:       VecDuplicate(u, &b);
1112:       VecSet(r, 0.0);
1113:       SNESComputeFunction(snes, r, b);
1114:       MatMult(A, u, r);
1115:       VecAXPY(r, 1.0, 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", (double)res);
1121:       /* check solver */
1122:       if (user.checkksp) {
1123:         KSP ksp;

1125:         if (nullSpace) {
1126:           MatNullSpaceRemove(nullSpace, u);
1127:         }
1128:         SNESComputeJacobian(snes, u, A, J);
1129:         MatMult(A, u, b);
1130:         SNESGetKSP(snes, &ksp);
1131:         KSPSetOperators(ksp, A, J);
1132:         KSPSolve(ksp, b, r);
1133:         VecAXPY(r, -1.0, u);
1134:         VecNorm(r, NORM_2, &res);
1135:         PetscPrintf(PETSC_COMM_WORLD, "KSP Error: %g\n", (double)res);
1136:       }
1137:       VecDestroy(&b);
1138:     }
1139:   }
1140:   VecViewFromOptions(u, NULL, "-vec_view");

1142:   if (user.bdIntegral) {
1143:     DMLabel   label;
1144:     PetscInt  id = 1;
1145:     PetscScalar bdInt = 0.0;
1146:     PetscReal   exact = 3.3333333333;

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

1154:   MatNullSpaceDestroy(&nullSpace);
1155:   if (user.jacobianMF) {VecDestroy(&userJ.u);}
1156:   if (A != J) {MatDestroy(&A);}
1157:   MatDestroy(&J);
1158:   VecDestroy(&u);
1159:   SNESDestroy(&snes);
1160:   DMDestroy(&dm);
1161:   PetscFree2(user.exactFuncs, user.exactFields);
1162:   PetscFinalize();
1163:   return ierr;
1164: }

1166: /*TEST
1167:   # 2D serial P1 test 0-4
1168:   test:
1169:     suffix: 2d_p1_0
1170:     requires: triangle
1171:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1173:   test:
1174:     suffix: 2d_p1_1
1175:     requires: triangle
1176:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1178:   test:
1179:     suffix: 2d_p1_2
1180:     requires: triangle
1181:     args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1183:   test:
1184:     suffix: 2d_p1_neumann_0
1185:     requires: triangle
1186:     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

1188:   test:
1189:     suffix: 2d_p1_neumann_1
1190:     requires: triangle
1191:     args: -run_type test -refinement_limit 0.0625 -bc_type neumann   -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1193:   # 2D serial P2 test 5-8
1194:   test:
1195:     suffix: 2d_p2_0
1196:     requires: triangle
1197:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1199:   test:
1200:     suffix: 2d_p2_1
1201:     requires: triangle
1202:     args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1204:   test:
1205:     suffix: 2d_p2_neumann_0
1206:     requires: triangle
1207:     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

1209:   test:
1210:     suffix: 2d_p2_neumann_1
1211:     requires: triangle
1212:     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

1214:   test:
1215:     suffix: bd_int_0
1216:     requires: triangle
1217:     args: -run_type test -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet

1219:   test:
1220:     suffix: bd_int_1
1221:     requires: triangle
1222:     args: -run_type test -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet

1224:   # 3D serial P1 test 9-12
1225:   test:
1226:     suffix: 3d_p1_0
1227:     requires: ctetgen
1228:     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

1230:   test:
1231:     suffix: 3d_p1_1
1232:     requires: ctetgen
1233:     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

1235:   test:
1236:     suffix: 3d_p1_2
1237:     requires: ctetgen
1238:     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

1240:   test:
1241:     suffix: 3d_p1_neumann_0
1242:     requires: ctetgen
1243:     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

1245:   # Analytic variable coefficient 13-20
1246:   test:
1247:     suffix: 13
1248:     requires: triangle
1249:     args: -run_type test -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1250:   test:
1251:     suffix: 14
1252:     requires: triangle
1253:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1254:   test:
1255:     suffix: 15
1256:     requires: triangle
1257:     args: -run_type test -refinement_limit 0.0    -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1258:   test:
1259:     suffix: 16
1260:     requires: triangle
1261:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1262:   test:
1263:     suffix: 17
1264:     requires: ctetgen
1265:     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

1267:   test:
1268:     suffix: 18
1269:     requires: ctetgen
1270:     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

1272:   test:
1273:     suffix: 19
1274:     requires: ctetgen
1275:     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

1277:   test:
1278:     suffix: 20
1279:     requires: ctetgen
1280:     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

1282:   # P1 variable coefficient 21-28
1283:   test:
1284:     suffix: 21
1285:     requires: triangle
1286:     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

1288:   test:
1289:     suffix: 22
1290:     requires: triangle
1291:     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

1293:   test:
1294:     suffix: 23
1295:     requires: triangle
1296:     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

1298:   test:
1299:     suffix: 24
1300:     requires: triangle
1301:     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

1303:   test:
1304:     suffix: 25
1305:     requires: ctetgen
1306:     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

1308:   test:
1309:     suffix: 26
1310:     requires: ctetgen
1311:     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

1313:   test:
1314:     suffix: 27
1315:     requires: ctetgen
1316:     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

1318:   test:
1319:     suffix: 28
1320:     requires: ctetgen
1321:     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

1323:   # P0 variable coefficient 29-36
1324:   test:
1325:     suffix: 29
1326:     requires: triangle
1327:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1329:   test:
1330:     suffix: 30
1331:     requires: triangle
1332:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1334:   test:
1335:     suffix: 31
1336:     requires: triangle
1337:     args: -run_type test -refinement_limit 0.0    -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1339:   test:
1340:     requires: triangle
1341:     suffix: 32
1342:     args: -run_type test -refinement_limit 0.0625 -variable_coefficient field    -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1344:   test:
1345:     requires: ctetgen
1346:     suffix: 33
1347:     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

1349:   test:
1350:     suffix: 34
1351:     requires: ctetgen
1352:     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

1354:   test:
1355:     suffix: 35
1356:     requires: ctetgen
1357:     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

1359:   test:
1360:     suffix: 36
1361:     requires: ctetgen
1362:     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

1364:   # Full solve 39-44
1365:   test:
1366:     suffix: 39
1367:     requires: triangle !single
1368:     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
1369:   test:
1370:     suffix: 40
1371:     requires: triangle !single
1372:     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
1373:   test:
1374:     suffix: 41
1375:     requires: triangle !single
1376:     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
1377:   test:
1378:     suffix: 42
1379:     requires: triangle !single
1380:     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
1381:   test:
1382:     suffix: 43
1383:     requires: triangle !single
1384:     nsize: 2
1385:     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

1387:   test:
1388:     suffix: 44
1389:     requires: triangle !single
1390:     nsize: 2
1391:     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

1393:   # These tests use a loose tolerance just to exercise the PtAP operations for MATIS and multiple PCBDDC setup calls inside PCMG
1394:   testset:
1395:     requires: triangle !single
1396:     nsize: 3
1397:     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
1398:     test:
1399:       suffix: gmg_bddc
1400:       filter: sed -e "s/CONVERGED_FNORM_RELATIVE iterations 3/CONVERGED_FNORM_RELATIVE iterations 4/g"
1401:       args: -mg_levels_pc_type jacobi
1402:     test:
1403:       filter: sed -e "s/iterations [0-4]/iterations 4/g"
1404:       suffix: gmg_bddc_lev
1405:       args: -mg_levels_pc_type bddc

1407:   # Restarting
1408:   testset:
1409:     suffix: restart
1410:     requires: hdf5 triangle !complex
1411:     args: -run_type test -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 1
1412:     test:
1413:       args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1414:     test:
1415:       args: -f sol.h5 -restart

1417:   # Periodicity
1418:   test:
1419:     suffix: periodic_0
1420:     requires: triangle
1421:     args: -run_type full -refinement_limit 0.0    -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -snes_converged_reason ::ascii_info_detail

1423:   test:
1424:     requires: !complex
1425:     suffix: periodic_1
1426:     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

1428:   # 2D serial P1 test with field bc
1429:   test:
1430:     suffix: field_bc_2d_p1_0
1431:     requires: triangle
1432:     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

1434:   test:
1435:     suffix: field_bc_2d_p1_1
1436:     requires: triangle
1437:     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

1439:   test:
1440:     suffix: field_bc_2d_p1_neumann_0
1441:     requires: triangle
1442:     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

1444:   test:
1445:     suffix: field_bc_2d_p1_neumann_1
1446:     requires: triangle
1447:     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

1449:   # 3D serial P1 test with field bc
1450:   test:
1451:     suffix: field_bc_3d_p1_0
1452:     requires: ctetgen
1453:     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

1455:   test:
1456:     suffix: field_bc_3d_p1_1
1457:     requires: ctetgen
1458:     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

1460:   test:
1461:     suffix: field_bc_3d_p1_neumann_0
1462:     requires: ctetgen
1463:     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

1465:   test:
1466:     suffix: field_bc_3d_p1_neumann_1
1467:     requires: ctetgen
1468:     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

1470:   # 2D serial P2 test with field bc
1471:   test:
1472:     suffix: field_bc_2d_p2_0
1473:     requires: triangle
1474:     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

1476:   test:
1477:     suffix: field_bc_2d_p2_1
1478:     requires: triangle
1479:     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

1481:   test:
1482:     suffix: field_bc_2d_p2_neumann_0
1483:     requires: triangle
1484:     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

1486:   test:
1487:     suffix: field_bc_2d_p2_neumann_1
1488:     requires: triangle
1489:     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

1491:   # 3D serial P2 test with field bc
1492:   test:
1493:     suffix: field_bc_3d_p2_0
1494:     requires: ctetgen
1495:     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

1497:   test:
1498:     suffix: field_bc_3d_p2_1
1499:     requires: ctetgen
1500:     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

1502:   test:
1503:     suffix: field_bc_3d_p2_neumann_0
1504:     requires: ctetgen
1505:     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

1507:   test:
1508:     suffix: field_bc_3d_p2_neumann_1
1509:     requires: ctetgen
1510:     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

1512:   # Full solve simplex: Convergence
1513:   test:
1514:     suffix: tet_conv_p1_r0
1515:     requires: ctetgen
1516:     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
1517:   test:
1518:     suffix: tet_conv_p1_r2
1519:     requires: ctetgen
1520:     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
1521:   test:
1522:     suffix: tet_conv_p1_r3
1523:     requires: ctetgen
1524:     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
1525:   test:
1526:     suffix: tet_conv_p2_r0
1527:     requires: ctetgen
1528:     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
1529:   test:
1530:     suffix: tet_conv_p2_r2
1531:     requires: ctetgen
1532:     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

1534:   # Full solve simplex: PCBDDC
1535:   test:
1536:     suffix: tri_bddc
1537:     requires: triangle !single
1538:     nsize: 5
1539:     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

1541:   # Full solve simplex: PCBDDC
1542:   test:
1543:     suffix: tri_parmetis_bddc
1544:     requires: triangle !single parmetis
1545:     nsize: 4
1546:     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

1548:   testset:
1549:     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
1550:     nsize: 5
1551:     output_file: output/ex12_quad_bddc.out
1552:     filter: sed -e "s/aijcusparse/aij/g" -e "s/aijviennacl/aij/g" -e "s/factorization: cusparse/factorization: petsc/g"
1553:     test:
1554:       requires: !single
1555:       suffix: quad_bddc
1556:     test:
1557:       requires: !single cuda
1558:       suffix: quad_bddc_cuda
1559:       args: -matis_localmat_type aijcusparse -pc_bddc_dirichlet_pc_factor_mat_solver_type cusparse -pc_bddc_neumann_pc_factor_mat_solver_type cusparse
1560:     test:
1561:       requires: !single viennacl
1562:       suffix: quad_bddc_viennacl
1563:       args: -matis_localmat_type aijviennacl

1565:   # Full solve simplex: ASM
1566:   test:
1567:     suffix: tri_q2q1_asm_lu
1568:     requires: triangle !single
1569:     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

1571:   test:
1572:     suffix: tri_q2q1_msm_lu
1573:     requires: triangle !single
1574:     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

1576:   test:
1577:     suffix: tri_q2q1_asm_sor
1578:     requires: triangle !single
1579:     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

1581:   test:
1582:     suffix: tri_q2q1_msm_sor
1583:     requires: triangle !single
1584:     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

1586:   # Full solve simplex: FAS
1587:   test:
1588:     suffix: fas_newton_0
1589:     requires: triangle !single
1590:     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

1592:   test:
1593:     suffix: fas_newton_1
1594:     requires: triangle !single
1595:     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
1596:     filter: sed -e "s/total number of linear solver iterations=14/total number of linear solver iterations=15/g"

1598:   test:
1599:     suffix: fas_ngs_0
1600:     requires: triangle !single
1601:     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

1603:   test:
1604:     suffix: fas_newton_coarse_0
1605:     requires: pragmatic triangle
1606:     TODO: broken
1607:     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

1609:   test:
1610:     suffix: mg_newton_coarse_0
1611:     requires: triangle pragmatic
1612:     TODO: broken
1613:     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

1615:   test:
1616:     suffix: mg_newton_coarse_1
1617:     requires: triangle pragmatic
1618:     TODO: broken
1619:     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

1621:   test:
1622:     suffix: mg_newton_coarse_2
1623:     requires: triangle pragmatic
1624:     TODO: broken
1625:     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

1627:   # Full solve tensor
1628:   test:
1629:     suffix: tensor_plex_2d
1630:     args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_refine_hierarchy 2 -cells 2,2

1632:   test:
1633:     suffix: tensor_p4est_2d
1634:     requires: p4est
1635:     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

1637:   test:
1638:     suffix: tensor_plex_3d
1639:     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

1641:   test:
1642:     suffix: tensor_p4est_3d
1643:     requires: p4est
1644:     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

1646:   test:
1647:     suffix: p4est_test_q2_conformal_serial
1648:     requires: p4est
1649:     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

1651:   test:
1652:     suffix: p4est_test_q2_conformal_parallel
1653:     requires: p4est
1654:     nsize: 7
1655:     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

1657:   test:
1658:     suffix: p4est_test_q2_conformal_parallel_parmetis
1659:     requires: parmetis p4est
1660:     nsize: 4
1661:     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

1663:   test:
1664:     suffix: p4est_test_q2_nonconformal_serial
1665:     requires: p4est
1666:     filter: grep -v "CG or CGNE: variant"
1667:     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

1669:   test:
1670:     suffix: p4est_test_q2_nonconformal_parallel
1671:     requires: p4est
1672:     filter: grep -v "CG or CGNE: variant"
1673:     nsize: 7
1674:     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

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

1682:   test:
1683:     suffix: p4est_exact_q2_conformal_serial
1684:     requires: p4est !single !complex !__float128
1685:     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

1687:   test:
1688:     suffix: p4est_exact_q2_conformal_parallel
1689:     requires: p4est !single !complex !__float128
1690:     nsize: 4
1691:     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

1693:   test:
1694:     suffix: p4est_exact_q2_conformal_parallel_parmetis
1695:     requires: parmetis p4est !single
1696:     nsize: 4
1697:     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

1699:   test:
1700:     suffix: p4est_exact_q2_nonconformal_serial
1701:     requires: p4est
1702:     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

1704:   test:
1705:     suffix: p4est_exact_q2_nonconformal_parallel
1706:     requires: p4est
1707:     nsize: 7
1708:     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

1710:   test:
1711:     suffix: p4est_exact_q2_nonconformal_parallel_parmetis
1712:     requires: parmetis p4est
1713:     nsize: 4
1714:     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

1716:   test:
1717:     suffix: p4est_full_q2_nonconformal_serial
1718:     requires: p4est !single
1719:     filter: grep -v "variant HERMITIAN"
1720:     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

1722:   test:
1723:     suffix: p4est_full_q2_nonconformal_parallel
1724:     requires: p4est !single
1725:     filter: grep -v "variant HERMITIAN"
1726:     nsize: 7
1727:     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

1729:   test:
1730:     suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1731:     requires: p4est !single
1732:     filter: grep -v "variant HERMITIAN"
1733:     nsize: 7
1734:     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

1736:   test:
1737:     suffix: p4est_full_q2_nonconformal_parallel_bddc
1738:     requires: p4est !single
1739:     filter: grep -v "variant HERMITIAN"
1740:     nsize: 7
1741:     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

1743:   test:
1744:     TODO: broken
1745:     suffix: p4est_fas_q2_conformal_serial
1746:     requires: p4est !complex !__float128
1747:     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

1749:   test:
1750:     TODO: broken
1751:     suffix: p4est_fas_q2_nonconformal_serial
1752:     requires: p4est
1753:     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

1755:   test:
1756:     suffix: fas_newton_0_p4est
1757:     requires: p4est !single !__float128
1758:     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

1760:   # Full solve simplicial AMR
1761:   test:
1762:     suffix: tri_p1_adapt_0
1763:     requires: pragmatic
1764:     TODO: broken
1765:     args: -run_type exact -dim 2 -dm_refine 5 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -variable_coefficient circle -snes_converged_reason ::ascii_info_detail -pc_type lu -adaptor_refinement_factor 1.0 -dm_view -dm_adapt_view -snes_adapt_initial 1

1767:   test:
1768:     suffix: tri_p1_adapt_1
1769:     requires: pragmatic
1770:     TODO: broken
1771:     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

1773:   test:
1774:     suffix: tri_p1_adapt_analytic_0
1775:     requires: pragmatic
1776:     TODO: broken
1777:     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

1779:   # Full solve tensor AMR
1780:   test:
1781:     suffix: quad_q1_adapt_0
1782:     requires: p4est
1783:     args: -run_type exact -dim 2 -simplex 0 -dm_plex_convert_type p4est -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -variable_coefficient circle -snes_converged_reason ::ascii_info_detail -pc_type lu -dm_forest_initial_refinement 4   -snes_adapt_initial 1 -dm_view
1784:     filter: grep -v DM_

1786:   test:
1787:     suffix: amr_0
1788:     nsize: 5
1789:     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

1791:   test:
1792:     suffix: amr_1
1793:     requires: p4est !complex
1794:     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

1796:   test:
1797:     suffix: p4est_solve_bddc
1798:     requires: p4est !complex
1799:     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
1800:     nsize: 4

1802:   test:
1803:     suffix: p4est_solve_fas
1804:     requires: p4est
1805:     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
1806:     nsize: 4
1807:     TODO: identical machine two runs produce slightly different solver trackers

1809:   test:
1810:     suffix: p4est_convergence_test_1
1811:     requires: p4est
1812:     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
1813:     nsize: 4

1815:   test:
1816:     suffix: p4est_convergence_test_2
1817:     requires: p4est
1818:     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

1820:   test:
1821:     suffix: p4est_convergence_test_3
1822:     requires: p4est
1823:     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

1825:   test:
1826:     suffix: p4est_convergence_test_4
1827:     requires: p4est
1828:     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
1829:     timeoutfactor: 5

1831:   # Serial tests with GLVis visualization
1832:   test:
1833:     suffix: glvis_2d_tet_p1
1834:     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
1835:   test:
1836:     suffix: glvis_2d_tet_p2
1837:     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
1838:   test:
1839:     suffix: glvis_2d_hex_p1
1840:     args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -simplex 0 -dm_refine 1
1841:   test:
1842:     suffix: glvis_2d_hex_p2
1843:     args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -simplex 0 -dm_refine 1
1844:   test:
1845:     suffix: glvis_2d_hex_p2_p4est
1846:     requires: p4est
1847:     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
1848:   test:
1849:     suffix: glvis_2d_tet_p0
1850:     args: -run_type exact  -interpolate 1 -guess_vec_view glvis: -nonzero_initial_guess 1 -f ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -petscspace_degree 0
1851:   test:
1852:     suffix: glvis_2d_hex_p0
1853:     args: -run_type exact  -interpolate 1 -guess_vec_view glvis: -nonzero_initial_guess 1 -cells 5,7  -simplex 0 -petscspace_degree 0

1855:   # PCHPDDM tests
1856:   testset:
1857:     nsize: 4
1858:     requires: hpddm slepc !single
1859:     args: -run_type test -run_test_check_ksp -quiet -petscspace_degree 1 -interpolate 1 -petscpartitioner_type simple -bc_type none -simplex 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 2 -pc_hpddm_coarse_p 1 -pc_hpddm_coarse_pc_type svd -ksp_rtol 1.e-10 -pc_hpddm_levels_1_st_pc_factor_shift_type INBLOCKS -ksp_converged_reason
1860:     test:
1861:       suffix: quad_singular_hpddm
1862:       args: -cells 6,7
1863:     test:
1864:       requires: p4est
1865:       suffix: p4est_singular_2d_hpddm
1866:       args: -dm_plex_convert_type p4est -dm_forest_minimum_refinement 1 -dm_forest_initial_refinement 3 -dm_forest_maximum_refinement 3
1867:     test:
1868:       requires: p4est
1869:       suffix: p4est_nc_singular_2d_hpddm
1870:       args: -dm_plex_convert_type p4est -dm_forest_minimum_refinement 1 -dm_forest_initial_refinement 1 -dm_forest_maximum_refinement 3 -dm_p4est_refine_pattern hash
1871:   testset:
1872:     nsize: 4
1873:     requires: hpddm slepc triangle !single
1874:     args: -run_type full -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 4 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-1
1875:     test:
1876:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1877:       suffix: tri_hpddm_reuse_baij
1878:     test:
1879:       requires: !complex
1880:       suffix: tri_hpddm_reuse
1881:   testset:
1882:     nsize: 4
1883:     requires: hpddm slepc !single
1884:     args: -run_type full -petscpartitioner_type simple -cells 7,5 -dm_refine 2 -simplex 0 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 4 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-1
1885:     test:
1886:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1887:       suffix: quad_hpddm_reuse_baij
1888:     test:
1889:       requires: !complex
1890:       suffix: quad_hpddm_reuse
1891:   testset:
1892:     nsize: 4
1893:     requires: hpddm slepc !single
1894:     args: -run_type full -petscpartitioner_type simple -cells 7,5 -dm_refine 2 -simplex 0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_threshold 0.1 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-1
1895:     test:
1896:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1897:       suffix: quad_hpddm_reuse_threshold_baij
1898:     test:
1899:       requires: !complex
1900:       suffix: quad_hpddm_reuse_threshold
1901:   testset:
1902:     nsize: 4
1903:     requires: hpddm slepc parmetis !single
1904:     filter: sed -e "s/linear solver iterations=17/linear solver iterations=16/g"
1905:     args: -run_type full -petscpartitioner_type parmetis -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -snes_converged_reason ::ascii_info_detail -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type icc -pc_hpddm_levels_1_eps_nev 20 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-10 -f ${PETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -pc_hpddm_levels_1_sub_pc_factor_levels 3 -variable_coefficient circle -dm_plex_gmsh_periodic 0
1906:     test:
1907:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1908:       suffix: tri_parmetis_hpddm_baij
1909:     test:
1910:       requires: !complex
1911:       suffix: tri_parmetis_hpddm
1912: TEST*/