Actual source code: ex12.c

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

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

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

 12: Information on refinement:

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

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

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

 27: typedef struct {
 28:   RunType        runType;           /* Whether to run tests, or solve the full problem */
 29:   PetscBool      jacobianMF;        /* Whether to calculate the Jacobian action on the fly */
 30:   PetscBool      showInitial, showSolution, restart, quiet, nonzInit;
 31:   /* Problem definition */
 32:   BCType         bcType;
 33:   CoeffType      variableCoefficient;
 34:   PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx);
 35:   PetscBool      fieldBC;
 36:   void           (**exactFields)(PetscInt, PetscInt, PetscInt,
 37:                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
 38:                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
 39:                                  PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]);
 40:   PetscBool      bdIntegral;        /* Compute the integral of the solution on the boundary */
 41:   /* Reproducing tests from SISC 40(3), pp. A1473-A1493, 2018 */
 42:   PetscInt       div;               /* Number of divisions */
 43:   PetscInt       k;                 /* Parameter for checkerboard coefficient */
 44:   PetscInt      *kgrid;             /* Random parameter grid */
 45:   PetscBool      rand;              /* Make random assignments */
 46:   /* Solver */
 47:   PC             pcmg;              /* This is needed for error monitoring */
 48:   PetscBool      checkksp;          /* Whether to check the KSPSolve for runType == RUN_TEST */
 49: } AppCtx;

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

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

 63: /*
 64:   In 2D for Dirichlet conditions, we use exact solution:

 66:     u = x^2 + y^2
 67:     f = 4

 69:   so that

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

 73:   For Neumann conditions, we have

 75:     -\nabla u \cdot -\hat y |_{y=0} =  (2y)|_{y=0} =  0 (bottom)
 76:     -\nabla u \cdot  \hat y |_{y=1} = -(2y)|_{y=1} = -2 (top)
 77:     -\nabla u \cdot -\hat x |_{x=0} =  (2x)|_{x=0} =  0 (left)
 78:     -\nabla u \cdot  \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)

 80:   Which we can express as

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

198:   so that

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

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

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

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

222:   so that

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

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

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

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

247:   so that

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

257: static PetscErrorCode checkerboardCoeff(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
258: {
259:   AppCtx  *user = (AppCtx *) ctx;
260:   PetscInt div  = user->div;
261:   PetscInt k    = user->k;
262:   PetscInt mask = 0, ind = 0, d;

265:   for (d = 0; d < dim; ++d) mask = (mask + (PetscInt) (x[d]*div)) % 2;
266:   if (user->kgrid) {
267:     for (d = 0; d < dim; ++d) {
268:       if (d > 0) ind *= dim;
269:       ind += (PetscInt) (x[d]*div);
270:     }
271:     k = user->kgrid[ind];
272:   }
273:   u[0] = mask ? 1.0 : PetscPowRealInt(10.0, -k);
274:   return(0);
275: }

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

285: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
286: void f1_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
287:                    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
288:                    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
289:                    PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
290: {
291:   PetscInt d;
292:   for (d = 0; d < dim; ++d) f1[d] = (x[0] + x[1])*u_x[d];
293: }

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

304: /* < \nabla v, \nabla u + {\nabla u}^T >
305:    This just gives \nabla u, give the perdiagonal for the transpose */
306: void g3_analytic_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
307:                     const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
308:                     const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
309:                     PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
310: {
311:   PetscInt d;
312:   for (d = 0; d < dim; ++d) g3[d*dim+d] = x[0] + x[1];
313: }

315: void g3_field_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
316:                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
317:                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
318:                  PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
319: {
320:   PetscInt d;
321:   for (d = 0; d < dim; ++d) g3[d*dim+d] = a[0];
322: }

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

327:     u  = x^2 + y^2
328:     f  = 16 (x^2 + y^2)
329:     nu = 1/2 |grad u|^2

331:   so that

333:     -\div \nu \grad u + f = -16 (x^2 + y^2) + 16 (x^2 + y^2) = 0
334: */
335: void f0_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
336:                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
337:                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
338:                              PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
339: {
340:   f0[0] = 16.0*(x[0]*x[0] + x[1]*x[1]);
341: }

343: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
344: void f1_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
345:                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
346:                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
347:                              PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
348: {
349:   PetscScalar nu = 0.0;
350:   PetscInt    d;
351:   for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
352:   for (d = 0; d < dim; ++d) f1[d] = 0.5*nu*u_x[d];
353: }

355: /*
356:   grad (u + eps w) - grad u = eps grad w

358:   1/2 |grad (u + eps w)|^2 grad (u + eps w) - 1/2 |grad u|^2 grad u
359: = 1/2 (|grad u|^2 + 2 eps <grad u,grad w>) (grad u + eps grad w) - 1/2 |grad u|^2 grad u
360: = 1/2 (eps |grad u|^2 grad w + 2 eps <grad u,grad w> grad u)
361: = eps (1/2 |grad u|^2 grad w + grad u <grad u,grad w>)
362: */
363: void g3_analytic_nonlinear_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
364:                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
365:                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
366:                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
367: {
368:   PetscScalar nu = 0.0;
369:   PetscInt    d, e;
370:   for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
371:   for (d = 0; d < dim; ++d) {
372:     g3[d*dim+d] = 0.5*nu;
373:     for (e = 0; e < dim; ++e) {
374:       g3[d*dim+e] += u_x[d]*u_x[e];
375:     }
376:   }
377: }

379: /*
380:   In 3D for Dirichlet conditions we use exact solution:

382:     u = 2/3 (x^2 + y^2 + z^2)
383:     f = 4

385:   so that

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

389:   For Neumann conditions, we have

391:     -\nabla u \cdot -\hat z |_{z=0} =  (2z)|_{z=0} =  0 (bottom)
392:     -\nabla u \cdot  \hat z |_{z=1} = -(2z)|_{z=1} = -2 (top)
393:     -\nabla u \cdot -\hat y |_{y=0} =  (2y)|_{y=0} =  0 (front)
394:     -\nabla u \cdot  \hat y |_{y=1} = -(2y)|_{y=1} = -2 (back)
395:     -\nabla u \cdot -\hat x |_{x=0} =  (2x)|_{x=0} =  0 (left)
396:     -\nabla u \cdot  \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)

398:   Which we can express as

400:     \nabla u \cdot  \hat n|_\Gamma = {2 x, 2 y, 2z} \cdot \hat n = 2 (x + y + z)
401: */
402: static PetscErrorCode quadratic_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
403: {
404:   *u = 2.0*(x[0]*x[0] + x[1]*x[1] + x[2]*x[2])/3.0;
405:   return 0;
406: }

408: static void quadratic_u_field_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
409:                                  const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
410:                                  const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
411:                                  PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
412: {
413:   uexact[0] = a[0];
414: }

416: static void bd_integral_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
417:                            const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
418:                            const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
419:                            PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar *uint)
420: {
421:   uint[0] = u[0];
422: }

424: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
425: {
426:   const char    *bcTypes[3]  = {"neumann", "dirichlet", "none"};
427:   const char    *runTypes[4] = {"full", "exact", "test", "perf"};
428:   const char    *coeffTypes[8] = {"none", "analytic", "field", "nonlinear", "circle", "cross", "checkerboard_0", "checkerboard_1"};
429:   PetscInt       bc, run, coeff;

433:   options->runType             = RUN_FULL;
434:   options->bcType              = DIRICHLET;
435:   options->variableCoefficient = COEFF_NONE;
436:   options->fieldBC             = PETSC_FALSE;
437:   options->jacobianMF          = PETSC_FALSE;
438:   options->showInitial         = PETSC_FALSE;
439:   options->showSolution        = PETSC_FALSE;
440:   options->restart             = PETSC_FALSE;
441:   options->quiet               = PETSC_FALSE;
442:   options->nonzInit            = PETSC_FALSE;
443:   options->bdIntegral          = PETSC_FALSE;
444:   options->checkksp            = PETSC_FALSE;
445:   options->div                 = 4;
446:   options->k                   = 1;
447:   options->kgrid               = NULL;
448:   options->rand                = PETSC_FALSE;

450:   PetscOptionsBegin(comm, "", "Poisson Problem Options", "DMPLEX");
451:   run  = options->runType;
452:   PetscOptionsEList("-run_type", "The run type", "ex12.c", runTypes, 4, runTypes[options->runType], &run, NULL);
453:   options->runType = (RunType) run;
454:   bc   = options->bcType;
455:   PetscOptionsEList("-bc_type","Type of boundary condition","ex12.c",bcTypes,3,bcTypes[options->bcType],&bc,NULL);
456:   options->bcType = (BCType) bc;
457:   coeff = options->variableCoefficient;
458:   PetscOptionsEList("-variable_coefficient","Type of variable coefficent","ex12.c",coeffTypes,8,coeffTypes[options->variableCoefficient],&coeff,NULL);
459:   options->variableCoefficient = (CoeffType) coeff;

461:   PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL);
462:   PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL);
463:   PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL);
464:   PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL);
465:   PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL);
466:   PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL);
467:   PetscOptionsBool("-nonzero_initial_guess", "nonzero initial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL);
468:   PetscOptionsBool("-bd_integral", "Compute the integral of the solution on the boundary", "ex12.c", options->bdIntegral, &options->bdIntegral, NULL);
469:   if (options->runType == RUN_TEST) {
470:     PetscOptionsBool("-run_test_check_ksp", "Check solution of KSP", "ex12.c", options->checkksp, &options->checkksp, NULL);
471:   }
472:   PetscOptionsInt("-div", "The number of division for the checkerboard coefficient", "ex12.c", options->div, &options->div, NULL);
473:   PetscOptionsInt("-k", "The exponent for the checkerboard coefficient", "ex12.c", options->k, &options->k, NULL);
474:   PetscOptionsBool("-k_random", "Assign random k values to checkerboard", "ex12.c", options->rand, &options->rand, NULL);
475:   PetscOptionsEnd();
476:   return(0);
477: }

479: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
480: {
481:   DM             plex;
482:   DMLabel        label;

486:   DMCreateLabel(dm, name);
487:   DMGetLabel(dm, name, &label);
488:   DMConvert(dm, DMPLEX, &plex);
489:   DMPlexMarkBoundaryFaces(plex, 1, label);
490:   DMDestroy(&plex);
491:   return(0);
492: }

494: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
495: {

499:   DMCreate(comm, dm);
500:   DMSetType(*dm, DMPLEX);
501:   DMSetFromOptions(*dm);
502:   {
503:     char      convType[256];
504:     PetscBool flg;

506:     PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
507:     PetscOptionsFList("-dm_plex_convert_type","Convert DMPlex to another format","ex12",DMList,DMPLEX,convType,256,&flg);
508:     PetscOptionsEnd();
509:     if (flg) {
510:       DM dmConv;

512:       DMConvert(*dm,convType,&dmConv);
513:       if (dmConv) {
514:         DMDestroy(dm);
515:         *dm  = dmConv;
516:       }
517:       DMSetFromOptions(*dm);
518:       DMSetUp(*dm);
519:     }
520:   }
521:   DMViewFromOptions(*dm, NULL, "-dm_view");
522:   if (user->rand) {
523:     PetscRandom r;
524:     PetscReal   val;
525:     PetscInt    dim, N, i;

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

544: static PetscErrorCode SetupProblem(DM dm, AppCtx *user)
545: {
546:   PetscDS         ds;
547:   DMLabel         label;
548:   PetscWeakForm   wf;
549:   const DMBoundaryType *periodicity;
550:   const PetscInt  id = 1;
551:   PetscInt        bd, dim;
552:   PetscErrorCode  ierr;

555:   DMGetDS(dm, &ds);
556:   DMGetDimension(dm, &dim);
557:   DMGetPeriodicity(dm, NULL, NULL, NULL, &periodicity);
558:   switch (user->variableCoefficient) {
559:   case COEFF_NONE:
560:     if (periodicity && periodicity[0]) {
561:       if (periodicity && periodicity[1]) {
562:         PetscDSSetResidual(ds, 0, f0_xytrig_u, f1_u);
563:         PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
564:       } else {
565:         PetscDSSetResidual(ds, 0, f0_xtrig_u,  f1_u);
566:         PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
567:       }
568:     } else {
569:       PetscDSSetResidual(ds, 0, f0_u, f1_u);
570:       PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
571:     }
572:     break;
573:   case COEFF_ANALYTIC:
574:     PetscDSSetResidual(ds, 0, f0_analytic_u, f1_analytic_u);
575:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_analytic_uu);
576:     break;
577:   case COEFF_FIELD:
578:     PetscDSSetResidual(ds, 0, f0_analytic_u, f1_field_u);
579:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_field_uu);
580:     break;
581:   case COEFF_NONLINEAR:
582:     PetscDSSetResidual(ds, 0, f0_analytic_nonlinear_u, f1_analytic_nonlinear_u);
583:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_analytic_nonlinear_uu);
584:     break;
585:   case COEFF_CIRCLE:
586:     PetscDSSetResidual(ds, 0, f0_circle_u, f1_u);
587:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
588:     break;
589:   case COEFF_CROSS:
590:     PetscDSSetResidual(ds, 0, f0_cross_u, f1_u);
591:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
592:     break;
593:   case COEFF_CHECKERBOARD_0:
594:     PetscDSSetResidual(ds, 0, f0_checkerboard_0_u, f1_field_u);
595:     PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_field_uu);
596:     break;
597:   default: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid variable coefficient type %d", user->variableCoefficient);
598:   }
599:   switch (dim) {
600:   case 2:
601:     switch (user->variableCoefficient) {
602:     case COEFF_CIRCLE:
603:       user->exactFuncs[0]  = circle_u_2d;break;
604:     case COEFF_CROSS:
605:       user->exactFuncs[0]  = cross_u_2d;break;
606:     case COEFF_CHECKERBOARD_0:
607:       user->exactFuncs[0]  = zero;break;
608:     default:
609:       if (periodicity && periodicity[0]) {
610:         if (periodicity && periodicity[1]) {
611:           user->exactFuncs[0] = xytrig_u_2d;
612:         } else {
613:           user->exactFuncs[0] = xtrig_u_2d;
614:         }
615:       } else {
616:         user->exactFuncs[0]  = quadratic_u_2d;
617:         user->exactFields[0] = quadratic_u_field_2d;
618:       }
619:     }
620:     if (user->bcType == NEUMANN) {
621:       DMGetLabel(dm, "boundary", &label);
622:       DMAddBoundary(dm, DM_BC_NATURAL, "wall", label, 1, &id, 0, 0, NULL, NULL, NULL, user, &bd);
623:       PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL);
624:       PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, 0, f0_bd_u, 0, NULL);
625:     }
626:     break;
627:   case 3:
628:     user->exactFuncs[0]  = quadratic_u_3d;
629:     user->exactFields[0] = quadratic_u_field_3d;
630:     if (user->bcType == NEUMANN) {
631:       DMGetLabel(dm, "boundary", &label);
632:       DMAddBoundary(dm, DM_BC_NATURAL, "wall", label, 1, &id, 0, 0, NULL, NULL, NULL, user, &bd);
633:       PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL);
634:       PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, 0, f0_bd_u, 0, NULL);
635:     }
636:     break;
637:   default:
638:     SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %d", dim);
639:   }
640:   /* Setup constants */
641:   switch (user->variableCoefficient) {
642:   case COEFF_CHECKERBOARD_0:
643:   {
644:     PetscScalar constants[2];

646:     constants[0] = user->div;
647:     constants[1] = user->k;
648:     PetscDSSetConstants(ds, 2, constants);
649:   }
650:   break;
651:   default: break;
652:   }
653:   PetscDSSetExactSolution(ds, 0, user->exactFuncs[0], user);
654:   /* Setup Boundary Conditions */
655:   if (user->bcType == DIRICHLET) {
656:     DMGetLabel(dm, "marker", &label);
657:     if (!label) {
658:       /* Right now, p4est cannot create labels immediately */
659:       PetscDSAddBoundaryByName(ds, user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL, "wall", "marker", 1, &id, 0, 0, NULL, user->fieldBC ? (void (*)(void)) user->exactFields[0] : (void (*)(void)) user->exactFuncs[0], NULL, user, NULL);
660:     } else {
661:       DMAddBoundary(dm, user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL, "wall", label, 1, &id, 0, 0, NULL, user->fieldBC ? (void (*)(void)) user->exactFields[0] : (void (*)(void)) user->exactFuncs[0], NULL, user, NULL);
662:     }
663:   }
664:   return(0);
665: }

667: static PetscErrorCode SetupMaterial(DM dm, DM dmAux, AppCtx *user)
668: {
669:   PetscErrorCode (*matFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {nu_2d};
670:   void            *ctx[1];
671:   Vec              nu;
672:   PetscErrorCode   ierr;

675:   ctx[0] = user;
676:   if (user->variableCoefficient == COEFF_CHECKERBOARD_0) {matFuncs[0] = checkerboardCoeff;}
677:   DMCreateLocalVector(dmAux, &nu);
678:   PetscObjectSetName((PetscObject) nu, "Coefficient");
679:   DMProjectFunctionLocal(dmAux, 0.0, matFuncs, ctx, INSERT_ALL_VALUES, nu);
680:   DMSetAuxiliaryVec(dm, NULL, 0, nu);
681:   VecDestroy(&nu);
682:   return(0);
683: }

685: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
686: {
687:   PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
688:   Vec            uexact;
689:   PetscInt       dim;

693:   DMGetDimension(dm, &dim);
694:   if (dim == 2) bcFuncs[0] = quadratic_u_2d;
695:   else          bcFuncs[0] = quadratic_u_3d;
696:   DMCreateLocalVector(dmAux, &uexact);
697:   DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
698:   DMSetAuxiliaryVec(dm, NULL, 0, uexact);
699:   VecDestroy(&uexact);
700:   return(0);
701: }

703: static PetscErrorCode SetupAuxDM(DM dm, PetscFE feAux, AppCtx *user)
704: {
705:   DM             dmAux, coordDM;

709:   /* MUST call DMGetCoordinateDM() in order to get p4est setup if present */
710:   DMGetCoordinateDM(dm, &coordDM);
711:   if (!feAux) return(0);
712:   DMClone(dm, &dmAux);
713:   DMSetCoordinateDM(dmAux, coordDM);
714:   DMSetField(dmAux, 0, NULL, (PetscObject) feAux);
715:   DMCreateDS(dmAux);
716:   if (user->fieldBC) {SetupBC(dm, dmAux, user);}
717:   else               {SetupMaterial(dm, dmAux, user);}
718:   DMDestroy(&dmAux);
719:   return(0);
720: }

722: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
723: {
724:   DM             plex, cdm = dm;
725:   PetscFE        fe, feAux = NULL;
726:   PetscBool      simplex;
727:   PetscInt       dim;
728:   MPI_Comm       comm;

732:   DMGetDimension(dm, &dim);
733:   DMConvert(dm, DMPLEX, &plex);
734:   DMPlexIsSimplex(plex, &simplex);
735:   DMDestroy(&plex);
736:   PetscObjectGetComm((PetscObject) dm, &comm);
737:   PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, NULL, -1, &fe);
738:   PetscObjectSetName((PetscObject) fe, "potential");
739:   if (user->variableCoefficient == COEFF_FIELD || user->variableCoefficient == COEFF_CHECKERBOARD_0) {
740:     PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, "mat_", -1, &feAux);
741:     PetscObjectSetName((PetscObject) feAux, "coefficient");
742:     PetscFECopyQuadrature(fe, feAux);
743:   } else if (user->fieldBC) {
744:     PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, "bc_", -1, &feAux);
745:     PetscFECopyQuadrature(fe, feAux);
746:   }
747:   /* Set discretization and boundary conditions for each mesh */
748:   DMSetField(dm, 0, NULL, (PetscObject) fe);
749:   DMCreateDS(dm);
750:   SetupProblem(dm, user);
751:   while (cdm) {
752:     SetupAuxDM(cdm, feAux, user);
753:     if (user->bcType == DIRICHLET) {
754:       PetscBool hasLabel;

756:       DMHasLabel(cdm, "marker", &hasLabel);
757:       if (!hasLabel) {CreateBCLabel(cdm, "marker");}
758:     }
759:     DMCopyDisc(dm, cdm);
760:     DMGetCoarseDM(cdm, &cdm);
761:   }
762:   PetscFEDestroy(&fe);
763:   PetscFEDestroy(&feAux);
764:   return(0);
765: }

767: int main(int argc, char **argv)
768: {
769:   DM             dm;          /* Problem specification */
770:   SNES           snes;        /* nonlinear solver */
771:   Vec            u;           /* solution vector */
772:   Mat            A,J;         /* Jacobian matrix */
773:   MatNullSpace   nullSpace;   /* May be necessary for Neumann conditions */
774:   AppCtx         user;        /* user-defined work context */
775:   JacActionCtx   userJ;       /* context for Jacobian MF action */
776:   PetscReal      error = 0.0; /* L_2 error in the solution */

779:   PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
780:   ProcessOptions(PETSC_COMM_WORLD, &user);
781:   SNESCreate(PETSC_COMM_WORLD, &snes);
782:   CreateMesh(PETSC_COMM_WORLD, &user, &dm);
783:   SNESSetDM(snes, dm);
784:   DMSetApplicationContext(dm, &user);

786:   PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
787:   SetupDiscretization(dm, &user);

789:   DMCreateGlobalVector(dm, &u);
790:   PetscObjectSetName((PetscObject) u, "potential");

792:   DMCreateMatrix(dm, &J);
793:   if (user.jacobianMF) {
794:     PetscInt M, m, N, n;

796:     MatGetSize(J, &M, &N);
797:     MatGetLocalSize(J, &m, &n);
798:     MatCreate(PETSC_COMM_WORLD, &A);
799:     MatSetSizes(A, m, n, M, N);
800:     MatSetType(A, MATSHELL);
801:     MatSetUp(A);
802: #if 0
803:     MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction);
804: #endif

806:     userJ.dm   = dm;
807:     userJ.J    = J;
808:     userJ.user = &user;

810:     DMCreateLocalVector(dm, &userJ.u);
811:     if (user.fieldBC) {DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);}
812:     else              {DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);}
813:     MatShellSetContext(A, &userJ);
814:   } else {
815:     A = J;
816:   }

818:   nullSpace = NULL;
819:   if (user.bcType != DIRICHLET) {
820:     MatNullSpaceCreate(PetscObjectComm((PetscObject) dm), PETSC_TRUE, 0, NULL, &nullSpace);
821:     MatSetNullSpace(A, nullSpace);
822:   }

824:   DMPlexSetSNESLocalFEM(dm,&user,&user,&user);
825:   SNESSetJacobian(snes, A, J, NULL, NULL);

827:   SNESSetFromOptions(snes);

829:   if (user.fieldBC) {DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u);}
830:   else              {DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);}
831:   if (user.restart) {
832: #if defined(PETSC_HAVE_HDF5)
833:     PetscViewer viewer;
834:     char        filename[PETSC_MAX_PATH_LEN];

836:     PetscOptionsGetString(NULL, NULL, "-dm_plex_filename", filename, sizeof(filename), NULL);
837:     PetscViewerCreate(PETSC_COMM_WORLD, &viewer);
838:     PetscViewerSetType(viewer, PETSCVIEWERHDF5);
839:     PetscViewerFileSetMode(viewer, FILE_MODE_READ);
840:     PetscViewerFileSetName(viewer, filename);
841:     PetscViewerHDF5PushGroup(viewer, "/fields");
842:     VecLoad(u, viewer);
843:     PetscViewerHDF5PopGroup(viewer);
844:     PetscViewerDestroy(&viewer);
845: #endif
846:   }
847:   if (user.showInitial) {
848:     Vec lv;
849:     DMGetLocalVector(dm, &lv);
850:     DMGlobalToLocalBegin(dm, u, INSERT_VALUES, lv);
851:     DMGlobalToLocalEnd(dm, u, INSERT_VALUES, lv);
852:     DMPrintLocalVec(dm, "Local function", 1.0e-10, lv);
853:     DMRestoreLocalVector(dm, &lv);
854:   }
855:   if (user.runType == RUN_FULL || user.runType == RUN_EXACT) {
856:     PetscErrorCode (*initialGuess[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {zero};

858:     if (user.nonzInit) initialGuess[0] = ecks;
859:     if (user.runType == RUN_FULL) {
860:       DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
861:     }
862:     VecViewFromOptions(u, NULL, "-guess_vec_view");
863:     SNESSolve(snes, NULL, u);
864:     SNESGetSolution(snes, &u);
865:     SNESGetDM(snes, &dm);

867:     if (user.showSolution) {
868:       PetscPrintf(PETSC_COMM_WORLD, "Solution\n");
869:       VecChop(u, 3.0e-9);
870:       VecView(u, PETSC_VIEWER_STDOUT_WORLD);
871:     }
872:   } else if (user.runType == RUN_PERF) {
873:     Vec       r;
874:     PetscReal res = 0.0;

876:     SNESGetFunction(snes, &r, NULL, NULL);
877:     SNESComputeFunction(snes, u, r);
878:     PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
879:     VecChop(r, 1.0e-10);
880:     VecNorm(r, NORM_2, &res);
881:     PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res);
882:   } else {
883:     Vec       r;
884:     PetscReal res = 0.0, tol = 1.0e-11;

886:     /* Check discretization error */
887:     SNESGetFunction(snes, &r, NULL, NULL);
888:     PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
889:     if (!user.quiet) {VecView(u, PETSC_VIEWER_STDOUT_WORLD);}
890:     DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
891:     if (error < tol) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < %2.1e\n", (double)tol);}
892:     else             {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", (double)error);}
893:     /* Check residual */
894:     SNESComputeFunction(snes, u, r);
895:     PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
896:     VecChop(r, 1.0e-10);
897:     if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
898:     VecNorm(r, NORM_2, &res);
899:     PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res);
900:     /* Check Jacobian */
901:     {
902:       Vec b;

904:       SNESComputeJacobian(snes, u, A, A);
905:       VecDuplicate(u, &b);
906:       VecSet(r, 0.0);
907:       SNESComputeFunction(snes, r, b);
908:       MatMult(A, u, r);
909:       VecAXPY(r, 1.0, b);
910:       PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n");
911:       VecChop(r, 1.0e-10);
912:       if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
913:       VecNorm(r, NORM_2, &res);
914:       PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", (double)res);
915:       /* check solver */
916:       if (user.checkksp) {
917:         KSP ksp;

919:         if (nullSpace) {
920:           MatNullSpaceRemove(nullSpace, u);
921:         }
922:         SNESComputeJacobian(snes, u, A, J);
923:         MatMult(A, u, b);
924:         SNESGetKSP(snes, &ksp);
925:         KSPSetOperators(ksp, A, J);
926:         KSPSolve(ksp, b, r);
927:         VecAXPY(r, -1.0, u);
928:         VecNorm(r, NORM_2, &res);
929:         PetscPrintf(PETSC_COMM_WORLD, "KSP Error: %g\n", (double)res);
930:       }
931:       VecDestroy(&b);
932:     }
933:   }
934:   VecViewFromOptions(u, NULL, "-vec_view");
935:   {
936:     Vec nu;

938:     DMGetAuxiliaryVec(dm, NULL, 0, &nu);
939:     if (nu) {VecViewFromOptions(nu, NULL, "-coeff_view");}
940:   }

942:   if (user.bdIntegral) {
943:     DMLabel   label;
944:     PetscInt  id = 1;
945:     PetscScalar bdInt = 0.0;
946:     PetscReal   exact = 3.3333333333;

948:     DMGetLabel(dm, "marker", &label);
949:     DMPlexComputeBdIntegral(dm, u, label, 1, &id, bd_integral_2d, &bdInt, NULL);
950:     PetscPrintf(PETSC_COMM_WORLD, "Solution boundary integral: %.4g\n", (double) PetscAbsScalar(bdInt));
951:     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);
952:   }

954:   MatNullSpaceDestroy(&nullSpace);
955:   if (user.jacobianMF) {VecDestroy(&userJ.u);}
956:   if (A != J) {MatDestroy(&A);}
957:   MatDestroy(&J);
958:   VecDestroy(&u);
959:   SNESDestroy(&snes);
960:   DMDestroy(&dm);
961:   PetscFree2(user.exactFuncs, user.exactFields);
962:   PetscFree(user.kgrid);
963:   PetscFinalize();
964:   return ierr;
965: }

967: /*TEST
968:   # 2D serial P1 test 0-4
969:   test:
970:     suffix: 2d_p1_0
971:     requires: triangle
972:     args: -run_type test -bc_type dirichlet -dm_plex_interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

974:   test:
975:     suffix: 2d_p1_1
976:     requires: triangle
977:     args: -run_type test -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

979:   test:
980:     suffix: 2d_p1_2
981:     requires: triangle
982:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

984:   test:
985:     suffix: 2d_p1_neumann_0
986:     requires: triangle
987:     args: -dm_coord_space 0 -run_type test -bc_type neumann -dm_plex_boundary_label boundary -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail

989:   test:
990:     suffix: 2d_p1_neumann_1
991:     requires: triangle
992:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -bc_type neumann -dm_plex_boundary_label boundary -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

994:   # 2D serial P2 test 5-8
995:   test:
996:     suffix: 2d_p2_0
997:     requires: triangle
998:     args: -run_type test -bc_type dirichlet -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1000:   test:
1001:     suffix: 2d_p2_1
1002:     requires: triangle
1003:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -bc_type dirichlet -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1005:   test:
1006:     suffix: 2d_p2_neumann_0
1007:     requires: triangle
1008:     args: -dm_coord_space 0 -run_type test -bc_type neumann -dm_plex_boundary_label boundary -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail

1010:   test:
1011:     suffix: 2d_p2_neumann_1
1012:     requires: triangle
1013:     args: -dm_coord_space 0 -run_type test -dm_refine_volume_limit_pre 0.0625 -bc_type neumann -dm_plex_boundary_label boundary -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail

1015:   test:
1016:     suffix: bd_int_0
1017:     requires: triangle
1018:     args: -run_type test -bc_type dirichlet -petscspace_degree 2 -bd_integral -dm_view -quiet

1020:   test:
1021:     suffix: bd_int_1
1022:     requires: triangle
1023:     args: -run_type test -dm_refine 2 -bc_type dirichlet -petscspace_degree 2 -bd_integral -dm_view -quiet

1025:   # 3D serial P1 test 9-12
1026:   test:
1027:     suffix: 3d_p1_0
1028:     requires: ctetgen
1029:     args: -run_type test -dm_plex_dim 3 -bc_type dirichlet -dm_plex_interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view

1031:   test:
1032:     suffix: 3d_p1_1
1033:     requires: ctetgen
1034:     args: -run_type test -dm_plex_dim 3 -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view

1036:   test:
1037:     suffix: 3d_p1_2
1038:     requires: ctetgen
1039:     args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view

1041:   test:
1042:     suffix: 3d_p1_neumann_0
1043:     requires: ctetgen
1044:     args: -run_type test -dm_plex_dim 3 -bc_type neumann -dm_plex_boundary_label boundary -petscspace_degree 1 -snes_fd -show_initial -dm_plex_print_fem 1 -dm_view

1046:   # Analytic variable coefficient 13-20
1047:   test:
1048:     suffix: 13
1049:     requires: triangle
1050:     args: -run_type test -variable_coefficient analytic -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1051:   test:
1052:     suffix: 14
1053:     requires: triangle
1054:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient analytic -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1055:   test:
1056:     suffix: 15
1057:     requires: triangle
1058:     args: -run_type test -variable_coefficient analytic -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1059:   test:
1060:     suffix: 16
1061:     requires: triangle
1062:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient analytic -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1063:   test:
1064:     suffix: 17
1065:     requires: ctetgen
1066:     args: -run_type test -dm_plex_dim 3 -variable_coefficient analytic -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1068:   test:
1069:     suffix: 18
1070:     requires: ctetgen
1071:     args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient analytic -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1073:   test:
1074:     suffix: 19
1075:     requires: ctetgen
1076:     args: -run_type test -dm_plex_dim 3 -variable_coefficient analytic -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1078:   test:
1079:     suffix: 20
1080:     requires: ctetgen
1081:     args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient analytic -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1083:   # P1 variable coefficient 21-28
1084:   test:
1085:     suffix: 21
1086:     requires: triangle
1087:     args: -run_type test -variable_coefficient field -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1089:   test:
1090:     suffix: 22
1091:     requires: triangle
1092:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient field -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1094:   test:
1095:     suffix: 23
1096:     requires: triangle
1097:     args: -run_type test -variable_coefficient field -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1099:   test:
1100:     suffix: 24
1101:     requires: triangle
1102:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient field -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1104:   test:
1105:     suffix: 25
1106:     requires: ctetgen
1107:     args: -run_type test -dm_plex_dim 3 -variable_coefficient field -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1109:   test:
1110:     suffix: 26
1111:     requires: ctetgen
1112:     args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient field -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1114:   test:
1115:     suffix: 27
1116:     requires: ctetgen
1117:     args: -run_type test -dm_plex_dim 3 -variable_coefficient field -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1119:   test:
1120:     suffix: 28
1121:     requires: ctetgen
1122:     args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient field -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1124:   # P0 variable coefficient 29-36
1125:   test:
1126:     suffix: 29
1127:     requires: triangle
1128:     args: -run_type test -variable_coefficient field -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1130:   test:
1131:     suffix: 30
1132:     requires: triangle
1133:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient field -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1135:   test:
1136:     suffix: 31
1137:     requires: triangle
1138:     args: -run_type test -variable_coefficient field -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1140:   test:
1141:     requires: triangle
1142:     suffix: 32
1143:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient field -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1145:   test:
1146:     requires: ctetgen
1147:     suffix: 33
1148:     args: -run_type test -dm_plex_dim 3 -variable_coefficient field -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1150:   test:
1151:     suffix: 34
1152:     requires: ctetgen
1153:     args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient field -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1155:   test:
1156:     suffix: 35
1157:     requires: ctetgen
1158:     args: -run_type test -dm_plex_dim 3 -variable_coefficient field -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1160:   test:
1161:     suffix: 36
1162:     requires: ctetgen
1163:     args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient field -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1165:   # Full solve 39-44
1166:   test:
1167:     suffix: 39
1168:     requires: triangle !single
1169:     args: -run_type full -dm_refine_volume_limit_pre 0.015625 -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
1170:   test:
1171:     suffix: 40
1172:     requires: triangle !single
1173:     args: -run_type full -dm_refine_volume_limit_pre 0.015625 -variable_coefficient nonlinear -petscspace_degree 2 -pc_type svd -ksp_rtol 1.0e-10 -snes_monitor_short -snes_converged_reason ::ascii_info_detail
1174:   test:
1175:     suffix: 41
1176:     requires: triangle !single
1177:     args: -run_type full -dm_refine_volume_limit_pre 0.03125 -variable_coefficient nonlinear -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
1178:   test:
1179:     suffix: 42
1180:     requires: triangle !single
1181:     args: -run_type full -dm_refine_volume_limit_pre 0.0625 -variable_coefficient nonlinear -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 -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
1182:   test:
1183:     suffix: 43
1184:     requires: triangle !single
1185:     nsize: 2
1186:     args: -run_type full -dm_distribute -dm_refine_volume_limit_pre 0.03125 -variable_coefficient nonlinear -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

1188:   test:
1189:     suffix: 44
1190:     requires: triangle !single
1191:     nsize: 2
1192:     args: -run_type full -dm_distribute -dm_refine_volume_limit_pre 0.0625 -variable_coefficient nonlinear -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

1194:   # These tests use a loose tolerance just to exercise the PtAP operations for MATIS and multiple PCBDDC setup calls inside PCMG
1195:   testset:
1196:     requires: triangle !single
1197:     nsize: 3
1198:     args: -run_type full -dm_distribute -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
1199:     test:
1200:       suffix: gmg_bddc
1201:       filter: sed -e "s/CONVERGED_FNORM_RELATIVE iterations 3/CONVERGED_FNORM_RELATIVE iterations 4/g"
1202:       args: -mg_levels_pc_type jacobi
1203:     test:
1204:       filter: sed -e "s/iterations [0-4]/iterations 4/g"
1205:       suffix: gmg_bddc_lev
1206:       args: -mg_levels_pc_type bddc

1208:   # Restarting
1209:   testset:
1210:     suffix: restart
1211:     requires: hdf5 triangle !complex
1212:     args: -run_type test -bc_type dirichlet -petscspace_degree 1
1213:     test:
1214:       args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1215:     test:
1216:       args: -dm_plex_filename sol.h5 -restart

1218:   # Periodicity
1219:   test:
1220:     suffix: periodic_0
1221:     requires: triangle
1222:     args: -run_type full -bc_type dirichlet -petscspace_degree 1 -snes_converged_reason ::ascii_info_detail

1224:   test:
1225:     requires: !complex
1226:     suffix: periodic_1
1227:     args: -quiet -run_type test -dm_plex_simplex 0 -dm_plex_box_faces 3,3 -dm_plex_box_bd periodic,periodic -vec_view vtk:test.vtu:vtk_vtu -petscspace_degree 1 -dm_refine 1

1229:   # 2D serial P1 test with field bc
1230:   test:
1231:     suffix: field_bc_2d_p1_0
1232:     requires: triangle
1233:     args: -run_type test -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1235:   test:
1236:     suffix: field_bc_2d_p1_1
1237:     requires: triangle
1238:     args: -run_type test -dm_refine 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1240:   test:
1241:     suffix: field_bc_2d_p1_neumann_0
1242:     requires: triangle
1243:     args: -run_type test -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1245:   test:
1246:     suffix: field_bc_2d_p1_neumann_1
1247:     requires: triangle
1248:     args: -run_type test -dm_refine 1 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1250:   # 3D serial P1 test with field bc
1251:   test:
1252:     suffix: field_bc_3d_p1_0
1253:     requires: ctetgen
1254:     args: -run_type test -dm_plex_dim 3 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1256:   test:
1257:     suffix: field_bc_3d_p1_1
1258:     requires: ctetgen
1259:     args: -run_type test -dm_plex_dim 3 -dm_refine 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1261:   test:
1262:     suffix: field_bc_3d_p1_neumann_0
1263:     requires: ctetgen
1264:     args: -run_type test -dm_plex_dim 3 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1266:   test:
1267:     suffix: field_bc_3d_p1_neumann_1
1268:     requires: ctetgen
1269:     args: -run_type test -dm_plex_dim 3 -dm_refine 1 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1271:   # 2D serial P2 test with field bc
1272:   test:
1273:     suffix: field_bc_2d_p2_0
1274:     requires: triangle
1275:     args: -run_type test -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1277:   test:
1278:     suffix: field_bc_2d_p2_1
1279:     requires: triangle
1280:     args: -run_type test -dm_refine 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1282:   test:
1283:     suffix: field_bc_2d_p2_neumann_0
1284:     requires: triangle
1285:     args: -run_type test -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1287:   test:
1288:     suffix: field_bc_2d_p2_neumann_1
1289:     requires: triangle
1290:     args: -run_type test -dm_refine 1 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1292:   # 3D serial P2 test with field bc
1293:   test:
1294:     suffix: field_bc_3d_p2_0
1295:     requires: ctetgen
1296:     args: -run_type test -dm_plex_dim 3 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1298:   test:
1299:     suffix: field_bc_3d_p2_1
1300:     requires: ctetgen
1301:     args: -run_type test -dm_plex_dim 3 -dm_refine 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1303:   test:
1304:     suffix: field_bc_3d_p2_neumann_0
1305:     requires: ctetgen
1306:     args: -run_type test -dm_plex_dim 3 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1308:   test:
1309:     suffix: field_bc_3d_p2_neumann_1
1310:     requires: ctetgen
1311:     args: -run_type test -dm_plex_dim 3 -dm_refine 1 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1313:   # Full solve simplex: Convergence
1314:   test:
1315:     suffix: 3d_p1_conv
1316:     requires: ctetgen
1317:     args: -run_type full -dm_plex_dim 3 -dm_refine 1 -bc_type dirichlet -petscspace_degree 1 \
1318:       -snes_convergence_estimate -convest_num_refine 1 -pc_type lu

1320:   # Full solve simplex: PCBDDC
1321:   test:
1322:     suffix: tri_bddc
1323:     requires: triangle !single
1324:     nsize: 5
1325:     args: -run_type full -dm_distribute -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -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

1327:   # Full solve simplex: PCBDDC
1328:   test:
1329:     suffix: tri_parmetis_bddc
1330:     requires: triangle !single parmetis
1331:     nsize: 4
1332:     args: -run_type full -dm_distribute -petscpartitioner_type parmetis -dm_refine 2 -bc_type dirichlet -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

1334:   testset:
1335:     args: -run_type full -dm_distribute -dm_plex_simplex 0 -dm_plex_box_faces 3,3 -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -petscspace_degree 2 -dm_mat_type is -pc_type bddc -ksp_type gmres -snes_monitor_short -ksp_monitor_short -snes_view -petscspace_poly_tensor -pc_bddc_corner_selection -ksp_rtol 1.e-9 -pc_bddc_use_edges 0
1336:     nsize: 5
1337:     output_file: output/ex12_quad_bddc.out
1338:     filter: sed -e "s/aijcusparse/aij/g" -e "s/aijviennacl/aij/g" -e "s/factorization: cusparse/factorization: petsc/g"
1339:     test:
1340:       requires: !single
1341:       suffix: quad_bddc
1342:     test:
1343:       requires: !single cuda
1344:       suffix: quad_bddc_cuda
1345:       args: -matis_localmat_type aijcusparse -pc_bddc_dirichlet_pc_factor_mat_solver_type cusparse -pc_bddc_neumann_pc_factor_mat_solver_type cusparse
1346:     test:
1347:       requires: !single viennacl
1348:       suffix: quad_bddc_viennacl
1349:       args: -matis_localmat_type aijviennacl

1351:   # Full solve simplex: ASM
1352:   test:
1353:     suffix: tri_q2q1_asm_lu
1354:     requires: triangle !single
1355:     args: -run_type full -dm_refine 3 -bc_type dirichlet -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

1357:   test:
1358:     suffix: tri_q2q1_msm_lu
1359:     requires: triangle !single
1360:     args: -run_type full -dm_refine 3 -bc_type dirichlet -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

1362:   test:
1363:     suffix: tri_q2q1_asm_sor
1364:     requires: triangle !single
1365:     args: -run_type full -dm_refine 3 -bc_type dirichlet -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

1367:   test:
1368:     suffix: tri_q2q1_msm_sor
1369:     requires: triangle !single
1370:     args: -run_type full -dm_refine 3 -bc_type dirichlet -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

1372:   # Full solve simplex: FAS
1373:   test:
1374:     suffix: fas_newton_0
1375:     requires: triangle !single
1376:     args: -run_type full -variable_coefficient nonlinear -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

1378:   test:
1379:     suffix: fas_newton_1
1380:     requires: triangle !single
1381:     args: -run_type full -dm_refine_hierarchy 3 -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
1382:     filter: sed -e "s/total number of linear solver iterations=14/total number of linear solver iterations=15/g"

1384:   test:
1385:     suffix: fas_ngs_0
1386:     requires: triangle !single
1387:     args: -run_type full -variable_coefficient nonlinear -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

1389:   test:
1390:     suffix: fas_newton_coarse_0
1391:     requires: pragmatic triangle
1392:     TODO: broken
1393:     args: -run_type full -dm_refine 2 -dm_plex_hash_location -variable_coefficient nonlinear -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

1395:   test:
1396:     suffix: mg_newton_coarse_0
1397:     requires: triangle pragmatic
1398:     TODO: broken
1399:     args: -run_type full -dm_refine 3 -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

1401:   test:
1402:     suffix: mg_newton_coarse_1
1403:     requires: triangle pragmatic
1404:     TODO: broken
1405:     args: -run_type full -dm_refine 5 -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

1407:   test:
1408:     suffix: mg_newton_coarse_2
1409:     requires: triangle pragmatic
1410:     TODO: broken
1411:     args: -run_type full -dm_refine 5 -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

1413:   # Full solve tensor
1414:   test:
1415:     suffix: tensor_plex_2d
1416:     args: -run_type test -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_refine_hierarchy 2

1418:   test:
1419:     suffix: tensor_p4est_2d
1420:     requires: p4est
1421:     args: -run_type test -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_forest_initial_refinement 2 -dm_forest_minimum_refinement 0 -dm_plex_convert_type p4est

1423:   test:
1424:     suffix: tensor_plex_3d
1425:     args: -run_type test -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_plex_dim 3 -dm_refine_hierarchy 1 -dm_plex_box_faces 2,2,2

1427:   test:
1428:     suffix: tensor_p4est_3d
1429:     requires: p4est
1430:     args: -run_type test -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_forest_initial_refinement 1 -dm_forest_minimum_refinement 0 -dm_plex_dim 3 -dm_plex_convert_type p8est -dm_plex_box_faces 2,2,2

1432:   test:
1433:     suffix: p4est_test_q2_conformal_serial
1434:     requires: p4est
1435:     args: -run_type test -petscspace_degree 2 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2

1437:   test:
1438:     suffix: p4est_test_q2_conformal_parallel
1439:     requires: p4est
1440:     nsize: 7
1441:     args: -run_type test -dm_distribute -petscspace_degree 2 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type simple

1443:   test:
1444:     suffix: p4est_test_q2_conformal_parallel_parmetis
1445:     requires: parmetis p4est
1446:     nsize: 4
1447:     args: -run_type test -dm_distribute -petscspace_degree 2 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type parmetis

1449:   test:
1450:     suffix: p4est_test_q2_nonconformal_serial
1451:     requires: p4est
1452:     filter: grep -v "CG or CGNE: variant"
1453:     args: -run_type test -petscspace_degree 2 -dm_plex_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

1455:   test:
1456:     suffix: p4est_test_q2_nonconformal_parallel
1457:     requires: p4est
1458:     filter: grep -v "CG or CGNE: variant"
1459:     nsize: 7
1460:     args: -run_type test -dm_distribute -petscspace_degree 2 -dm_plex_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

1462:   test:
1463:     suffix: p4est_test_q2_nonconformal_parallel_parmetis
1464:     requires: parmetis p4est
1465:     nsize: 4
1466:     args: -run_type test -dm_distribute -petscspace_degree 2 -dm_plex_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

1468:   test:
1469:     suffix: p4est_exact_q2_conformal_serial
1470:     requires: p4est !single !complex !__float128
1471:     args: -run_type exact -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 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2

1473:   test:
1474:     suffix: p4est_exact_q2_conformal_parallel
1475:     requires: p4est !single !complex !__float128
1476:     nsize: 4
1477:     args: -run_type exact -dm_distribute -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 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2

1479:   test:
1480:     suffix: p4est_exact_q2_conformal_parallel_parmetis
1481:     requires: parmetis p4est !single
1482:     nsize: 4
1483:     args: -run_type exact -dm_distribute -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 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type parmetis

1485:   test:
1486:     suffix: p4est_exact_q2_nonconformal_serial
1487:     requires: p4est
1488:     args: -run_type exact -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 -dm_plex_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

1490:   test:
1491:     suffix: p4est_exact_q2_nonconformal_parallel
1492:     requires: p4est
1493:     nsize: 7
1494:     args: -run_type exact -dm_distribute -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 -dm_plex_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

1496:   test:
1497:     suffix: p4est_exact_q2_nonconformal_parallel_parmetis
1498:     requires: parmetis p4est
1499:     nsize: 4
1500:     args: -run_type exact -dm_distribute -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 -dm_plex_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

1502:   test:
1503:     suffix: p4est_full_q2_nonconformal_serial
1504:     requires: p4est !single
1505:     filter: grep -v "variant HERMITIAN"
1506:     args: -run_type full -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 -dm_plex_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

1508:   test:
1509:     suffix: p4est_full_q2_nonconformal_parallel
1510:     requires: p4est !single
1511:     filter: grep -v "variant HERMITIAN"
1512:     nsize: 7
1513:     args: -run_type full -dm_distribute -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 -dm_plex_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

1515:   test:
1516:     suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1517:     requires: p4est !single
1518:     filter: grep -v "variant HERMITIAN"
1519:     nsize: 7
1520:     args: -run_type full -dm_distribute -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 -dm_plex_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

1522:   test:
1523:     suffix: p4est_full_q2_nonconformal_parallel_bddc
1524:     requires: p4est !single
1525:     filter: grep -v "variant HERMITIAN"
1526:     nsize: 7
1527:     args: -run_type full -dm_distribute -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 -dm_plex_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

1529:   test:
1530:     TODO: broken
1531:     suffix: p4est_fas_q2_conformal_serial
1532:     requires: p4est !complex !__float128
1533:     args: -run_type full -variable_coefficient nonlinear -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 -dm_plex_simplex 0 -dm_refine_hierarchy 3

1535:   test:
1536:     TODO: broken
1537:     suffix: p4est_fas_q2_nonconformal_serial
1538:     requires: p4est
1539:     args: -run_type full -variable_coefficient nonlinear -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 -dm_plex_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

1541:   test:
1542:     suffix: fas_newton_0_p4est
1543:     requires: p4est !single !__float128
1544:     args: -run_type full -variable_coefficient nonlinear -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 -dm_plex_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

1546:   # Full solve simplicial AMR
1547:   test:
1548:     suffix: tri_p1_adapt_0
1549:     requires: pragmatic
1550:     args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient circle -snes_converged_reason ::ascii_info_detail -pc_type lu -dm_view -dm_adapt_view -snes_adapt_initial 1 -adaptor_target_num 4000 -adaptor_metric_h_min 0.0001 -adaptor_metric_h_max 0.05

1552:   test:
1553:     suffix: tri_p1_adapt_1
1554:     requires: pragmatic
1555:     args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient circle -snes_converged_reason ::ascii_info_detail -pc_type lu -dm_view -dm_adapt_iter_view -dm_adapt_view -snes_adapt_sequence 2 -adaptor_target_num 4000 -adaptor_metric_h_min 0.0001 -adaptor_metric_h_max 0.05

1557:   test:
1558:     suffix: tri_p1_adapt_analytic_0
1559:     requires: pragmatic
1560:     args: -run_type exact -dm_refine 3 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -adaptor_monitor -dm_view -dm_adapt_iter_view -adaptor_metric_h_min 0.0001 -adaptor_metric_h_max 0.05

1562:   # Full solve tensor AMR
1563:   test:
1564:     suffix: quad_q1_adapt_0
1565:     requires: p4est
1566:     args: -run_type exact -dm_plex_simplex 0 -dm_plex_convert_type p4est -bc_type dirichlet -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
1567:     filter: grep -v DM_

1569:   test:
1570:     suffix: amr_0
1571:     nsize: 5
1572:     args: -run_type test -dm_distribute -petscpartitioner_type simple -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_refine 1

1574:   test:
1575:     suffix: amr_1
1576:     requires: p4est !complex
1577:     args: -run_type test -dm_plex_simplex 0 -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

1579:   test:
1580:     suffix: p4est_solve_bddc
1581:     requires: p4est !complex
1582:     args: -run_type full -dm_distribute -variable_coefficient nonlinear -nonzero_initial_guess 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 -dm_plex_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
1583:     nsize: 4

1585:   test:
1586:     suffix: p4est_solve_fas
1587:     requires: p4est
1588:     args: -run_type full -dm_distribute -variable_coefficient nonlinear -nonzero_initial_guess 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 -dm_plex_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
1589:     nsize: 4
1590:     TODO: identical machine two runs produce slightly different solver trackers

1592:   test:
1593:     suffix: p4est_convergence_test_1
1594:     requires: p4est
1595:     args:  -quiet -run_type test -dm_distribute -petscspace_degree 1 -dm_plex_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
1596:     nsize: 4

1598:   test:
1599:     suffix: p4est_convergence_test_2
1600:     requires: p4est
1601:     args: -quiet -run_type test -petscspace_degree 1 -dm_plex_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

1603:   test:
1604:     suffix: p4est_convergence_test_3
1605:     requires: p4est
1606:     args: -quiet -run_type test -petscspace_degree 1 -dm_plex_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

1608:   test:
1609:     suffix: p4est_convergence_test_4
1610:     requires: p4est
1611:     args: -quiet -run_type test -petscspace_degree 1 -dm_plex_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
1612:     timeoutfactor: 5

1614:   # Serial tests with GLVis visualization
1615:   test:
1616:     suffix: glvis_2d_tet_p1
1617:     args: -quiet -run_type test -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_plex_boundary_label marker -dm_plex_gmsh_periodic 0 -dm_coord_space 0
1618:   test:
1619:     suffix: glvis_2d_tet_p2
1620:     args: -quiet -run_type test -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_plex_boundary_label marker -dm_plex_gmsh_periodic 0 -dm_coord_space 0
1621:   test:
1622:     suffix: glvis_2d_hex_p1
1623:     args: -quiet -run_type test -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -dm_plex_simplex 0 -dm_refine 1 -dm_coord_space 0
1624:   test:
1625:     suffix: glvis_2d_hex_p2
1626:     args: -quiet -run_type test -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -dm_plex_simplex 0 -dm_refine 1 -dm_coord_space 0
1627:   test:
1628:     suffix: glvis_2d_hex_p2_p4est
1629:     requires: p4est
1630:     args: -quiet -run_type test -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -dm_plex_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 -viewer_glvis_dm_plex_enable_ncmesh
1631:   test:
1632:     suffix: glvis_2d_tet_p0
1633:     args: -run_type exact -guess_vec_view glvis: -nonzero_initial_guess 1 -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_plex_boundary_label marker -petscspace_degree 0 -dm_coord_space 0
1634:   test:
1635:     suffix: glvis_2d_hex_p0
1636:     args: -run_type exact -guess_vec_view glvis: -nonzero_initial_guess 1 -dm_plex_box_faces 5,7 -dm_plex_simplex 0 -petscspace_degree 0 -dm_coord_space 0

1638:   # PCHPDDM tests
1639:   testset:
1640:     nsize: 4
1641:     requires: hpddm slepc !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1642:     args: -run_type test -run_test_check_ksp -quiet -petscspace_degree 1 -dm_distribute -petscpartitioner_type simple -bc_type none -dm_plex_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
1643:     test:
1644:       suffix: quad_singular_hpddm
1645:       args: -dm_plex_box_faces 6,7
1646:     test:
1647:       requires: p4est
1648:       suffix: p4est_singular_2d_hpddm
1649:       args: -dm_plex_convert_type p4est -dm_forest_minimum_refinement 1 -dm_forest_initial_refinement 3 -dm_forest_maximum_refinement 3
1650:     test:
1651:       requires: p4est
1652:       suffix: p4est_nc_singular_2d_hpddm
1653:       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
1654:   testset:
1655:     nsize: 4
1656:     requires: hpddm slepc triangle !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1657:     args: -run_type full -dm_distribute -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -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
1658:     test:
1659:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1660:       suffix: tri_hpddm_reuse_baij
1661:     test:
1662:       requires: !complex
1663:       suffix: tri_hpddm_reuse
1664:   testset:
1665:     nsize: 4
1666:     requires: hpddm slepc !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1667:     args: -run_type full -dm_distribute -petscpartitioner_type simple -dm_plex_box_faces 7,5 -dm_refine 2 -dm_plex_simplex 0 -bc_type dirichlet -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
1668:     test:
1669:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1670:       suffix: quad_hpddm_reuse_baij
1671:     test:
1672:       requires: !complex
1673:       suffix: quad_hpddm_reuse
1674:   testset:
1675:     nsize: 4
1676:     requires: hpddm slepc !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1677:     args: -run_type full -dm_distribute -petscpartitioner_type simple -dm_plex_box_faces 7,5 -dm_refine 2 -dm_plex_simplex 0 -bc_type dirichlet -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
1678:     test:
1679:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1680:       suffix: quad_hpddm_reuse_threshold_baij
1681:     test:
1682:       requires: !complex
1683:       suffix: quad_hpddm_reuse_threshold
1684:   testset:
1685:     nsize: 4
1686:     requires: hpddm slepc parmetis !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1687:     filter: sed -e "s/linear solver iterations=17/linear solver iterations=16/g"
1688:     args: -run_type full -dm_distribute -petscpartitioner_type parmetis -dm_refine 3 -bc_type dirichlet -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 -dm_plex_filename ${PETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_plex_boundary_label marker -pc_hpddm_levels_1_sub_pc_factor_levels 3 -variable_coefficient circle -dm_plex_gmsh_periodic 0
1689:     test:
1690:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1691:       filter: grep -v "      total: nonzeros=" | grep -v "      rows=" | sed -e "s/total number of linear solver iterations=1[5-7]/total number of linear solver iterations=16/g"
1692:       suffix: tri_parmetis_hpddm_baij
1693:     test:
1694:       filter: grep -v "      total: nonzeros=" | grep -v "      rows=" | sed -e "s/total number of linear solver iterations=1[5-7]/total number of linear solver iterations=16/g"
1695:       requires: !complex
1696:       suffix: tri_parmetis_hpddm

1698:   # 2D serial P1 tests for adaptive MG
1699:   test:
1700:     suffix: 2d_p1_adaptmg_0
1701:     requires: triangle bamg
1702:     args: -dm_refine_hierarchy 3 -dm_plex_box_faces 4,4 -bc_type dirichlet -petscspace_degree 1 \
1703:           -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
1704:           -snes_max_it 1 -ksp_converged_reason \
1705:           -ksp_rtol 1e-8 -pc_type mg
1706:   # -ksp_monitor_true_residual -ksp_converged_reason -mg_levels_ksp_monitor_true_residual -pc_mg_mesp_monitor -dm_adapt_interp_view_fine draw -dm_adapt_interp_view_coarse draw -draw_pause 1
1707:   test:
1708:     suffix: 2d_p1_adaptmg_1
1709:     requires: triangle bamg
1710:     args: -dm_refine_hierarchy 3 -dm_plex_box_faces 4,4 -bc_type dirichlet -petscspace_degree 1 \
1711:           -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
1712:           -snes_max_it 1 -ksp_converged_reason \
1713:           -ksp_rtol 1e-8 -pc_type mg -pc_mg_galerkin -pc_mg_adapt_interp -pc_mg_adapt_interp_coarse_space eigenvector -pc_mg_adapt_interp_n 1 \
1714:             -pc_mg_mesp_ksp_type richardson -pc_mg_mesp_ksp_richardson_self_scale -pc_mg_mesp_ksp_max_it 100 -pc_mg_mesp_pc_type none

1716: TEST*/