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 {
 24:   NEUMANN,
 25:   DIRICHLET,
 26:   NONE
 27: } BCType;
 28: typedef enum {
 29:   RUN_FULL,
 30:   RUN_EXACT,
 31:   RUN_TEST,
 32:   RUN_PERF
 33: } RunType;
 34: typedef enum {
 35:   COEFF_NONE,
 36:   COEFF_ANALYTIC,
 37:   COEFF_FIELD,
 38:   COEFF_NONLINEAR,
 39:   COEFF_BALL,
 40:   COEFF_CROSS,
 41:   COEFF_CHECKERBOARD_0,
 42:   COEFF_CHECKERBOARD_1
 43: } CoeffType;

 45: typedef struct {
 46:   RunType   runType;    /* Whether to run tests, or solve the full problem */
 47:   PetscBool jacobianMF; /* Whether to calculate the Jacobian action on the fly */
 48:   PetscBool showInitial, showSolution, restart, quiet, nonzInit;
 49:   /* Problem definition */
 50:   BCType    bcType;
 51:   CoeffType variableCoefficient;
 52:   PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx);
 53:   PetscBool fieldBC;
 54:   void (**exactFields)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]);
 55:   PetscBool bdIntegral; /* Compute the integral of the solution on the boundary */
 56:   /* Reproducing tests from SISC 40(3), pp. A1473-A1493, 2018 */
 57:   PetscInt  div;   /* Number of divisions */
 58:   PetscInt  k;     /* Parameter for checkerboard coefficient */
 59:   PetscInt *kgrid; /* Random parameter grid */
 60:   PetscBool rand;  /* Make random assignments */
 61:   /* Solver */
 62:   PC        pcmg;     /* This is needed for error monitoring */
 63:   PetscBool checkksp; /* Whether to check the KSPSolve for runType == RUN_TEST */
 64: } AppCtx;

 66: static PetscErrorCode zero(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 67: {
 68:   u[0] = 0.0;
 69:   return PETSC_SUCCESS;
 70: }

 72: static PetscErrorCode ecks(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 73: {
 74:   u[0] = x[0];
 75:   return PETSC_SUCCESS;
 76: }

 78: /*
 79:   In 2D for Dirichlet conditions, we use exact solution:

 81:     u = x^2 + y^2
 82:     f = 4

 84:   so that

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

 88:   For Neumann conditions, we have

 90:     -\nabla u \cdot -\hat y |_{y=0} =  (2y)|_{y=0} =  0 (bottom)
 91:     -\nabla u \cdot  \hat y |_{y=1} = -(2y)|_{y=1} = -2 (top)
 92:     -\nabla u \cdot -\hat x |_{x=0} =  (2x)|_{x=0} =  0 (left)
 93:     -\nabla u \cdot  \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)

 95:   Which we can express as

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

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

101:     \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
102: */
103: static PetscErrorCode quadratic_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
104: {
105:   *u = x[0] * x[0] + x[1] * x[1];
106:   return PETSC_SUCCESS;
107: }

109: static void quadratic_u_field_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
110: {
111:   uexact[0] = a[0];
112: }

114: static PetscErrorCode ball_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
115: {
116:   const PetscReal alpha   = 500.;
117:   const PetscReal radius2 = PetscSqr(0.15);
118:   const PetscReal r2      = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5);
119:   const PetscReal xi      = alpha * (radius2 - r2);

121:   *u = PetscTanhScalar(xi) + 1.0;
122:   return PETSC_SUCCESS;
123: }

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

130:   *u = PetscSinReal(alpha * xy) * (alpha * PetscAbsReal(xy) < 2 * PETSC_PI ? 1 : 0.01);
131:   return PETSC_SUCCESS;
132: }

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

139: static void f0_ball_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
140: {
141:   PetscInt        d;
142:   const PetscReal alpha = 500., radius2 = PetscSqr(0.15);
143:   PetscReal       r2, xi;

145:   for (d = 0, r2 = 0.0; d < dim; ++d) r2 += PetscSqr(x[d] - 0.5);
146:   xi    = alpha * (radius2 - r2);
147:   f0[0] = (-2.0 * dim * alpha - 8.0 * PetscSqr(alpha) * r2 * PetscTanhReal(xi)) * PetscSqr(1.0 / PetscCoshReal(xi));
148: }

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

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

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

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

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

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

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

187:     u = sin(2 pi x)
188:     f = -4 pi^2 sin(2 pi x)

190:   so that

192:     -\Delta u + f = 4 pi^2 sin(2 pi x) - 4 pi^2 sin(2 pi x) = 0
193: */
194: static PetscErrorCode xtrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
195: {
196:   *u = PetscSinReal(2.0 * PETSC_PI * x[0]);
197:   return PETSC_SUCCESS;
198: }

200: static void f0_xtrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
201: {
202:   f0[0] = -4.0 * PetscSqr(PETSC_PI) * PetscSinReal(2.0 * PETSC_PI * x[0]);
203: }

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

208:     u = sin(2 pi x) sin(2 pi y)
209:     f = -8 pi^2 sin(2 pi x)

211:   so that

213:     -\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
214: */
215: static PetscErrorCode xytrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
216: {
217:   *u = PetscSinReal(2.0 * PETSC_PI * x[0]) * PetscSinReal(2.0 * PETSC_PI * x[1]);
218:   return PETSC_SUCCESS;
219: }

221: static void f0_xytrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
222: {
223:   f0[0] = -8.0 * PetscSqr(PETSC_PI) * PetscSinReal(2.0 * PETSC_PI * x[0]);
224: }

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

229:     u  = x^2 + y^2
230:     f  = 6 (x + y)
231:     nu = (x + y)

233:   so that

235:     -\div \nu \grad u + f = -6 (x + y) + 6 (x + y) = 0
236: */
237: static PetscErrorCode nu_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
238: {
239:   *u = x[0] + x[1];
240:   return PETSC_SUCCESS;
241: }

243: static PetscErrorCode checkerboardCoeff(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
244: {
245:   AppCtx  *user = (AppCtx *)ctx;
246:   PetscInt div  = user->div;
247:   PetscInt k    = user->k;
248:   PetscInt mask = 0, ind = 0, d;

250:   PetscFunctionBeginUser;
251:   for (d = 0; d < dim; ++d) mask = (mask + (PetscInt)(x[d] * div)) % 2;
252:   if (user->kgrid) {
253:     for (d = 0; d < dim; ++d) {
254:       if (d > 0) ind *= dim;
255:       ind += (PetscInt)(x[d] * div);
256:     }
257:     k = user->kgrid[ind];
258:   }
259:   u[0] = mask ? 1.0 : PetscPowRealInt(10.0, -k);
260:   PetscFunctionReturn(PETSC_SUCCESS);
261: }

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

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

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

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

289: void g3_field_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
290: {
291:   PetscInt d;
292:   for (d = 0; d < dim; ++d) g3[d * dim + d] = a[0];
293: }

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

298:     u  = x^2 + y^2
299:     f  = 16 (x^2 + y^2)
300:     nu = 1/2 |grad u|^2

302:   so that

304:     -\div \nu \grad u + f = -16 (x^2 + y^2) + 16 (x^2 + y^2) = 0
305: */
306: void f0_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
307: {
308:   f0[0] = 16.0 * (x[0] * x[0] + x[1] * x[1]);
309: }

311: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
312: void f1_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
313: {
314:   PetscScalar nu = 0.0;
315:   PetscInt    d;
316:   for (d = 0; d < dim; ++d) nu += u_x[d] * u_x[d];
317:   for (d = 0; d < dim; ++d) f1[d] = 0.5 * nu * u_x[d];
318: }

320: /*
321:   grad (u + eps w) - grad u = eps grad w

323:   1/2 |grad (u + eps w)|^2 grad (u + eps w) - 1/2 |grad u|^2 grad u
324: = 1/2 (|grad u|^2 + 2 eps <grad u,grad w>) (grad u + eps grad w) - 1/2 |grad u|^2 grad u
325: = 1/2 (eps |grad u|^2 grad w + 2 eps <grad u,grad w> grad u)
326: = eps (1/2 |grad u|^2 grad w + grad u <grad u,grad w>)
327: */
328: void g3_analytic_nonlinear_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
329: {
330:   PetscScalar nu = 0.0;
331:   PetscInt    d, e;
332:   for (d = 0; d < dim; ++d) nu += u_x[d] * u_x[d];
333:   for (d = 0; d < dim; ++d) {
334:     g3[d * dim + d] = 0.5 * nu;
335:     for (e = 0; e < dim; ++e) g3[d * dim + e] += u_x[d] * u_x[e];
336:   }
337: }

339: /*
340:   In 3D for Dirichlet conditions we use exact solution:

342:     u = 2/3 (x^2 + y^2 + z^2)
343:     f = 4

345:   so that

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

349:   For Neumann conditions, we have

351:     -\nabla u \cdot -\hat z |_{z=0} =  (2z)|_{z=0} =  0 (bottom)
352:     -\nabla u \cdot  \hat z |_{z=1} = -(2z)|_{z=1} = -2 (top)
353:     -\nabla u \cdot -\hat y |_{y=0} =  (2y)|_{y=0} =  0 (front)
354:     -\nabla u \cdot  \hat y |_{y=1} = -(2y)|_{y=1} = -2 (back)
355:     -\nabla u \cdot -\hat x |_{x=0} =  (2x)|_{x=0} =  0 (left)
356:     -\nabla u \cdot  \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)

358:   Which we can express as

360:     \nabla u \cdot  \hat n|_\Gamma = {2 x, 2 y, 2z} \cdot \hat n = 2 (x + y + z)
361: */
362: static PetscErrorCode quadratic_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
363: {
364:   *u = 2.0 * (x[0] * x[0] + x[1] * x[1] + x[2] * x[2]) / 3.0;
365:   return PETSC_SUCCESS;
366: }

368: static PetscErrorCode ball_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
369: {
370:   const PetscReal alpha   = 500.;
371:   const PetscReal radius2 = PetscSqr(0.15);
372:   const PetscReal r2      = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5) + PetscSqr(x[2] - 0.5);
373:   const PetscReal xi      = alpha * (radius2 - r2);

375:   *u = PetscTanhScalar(xi) + 1.0;
376:   return PETSC_SUCCESS;
377: }

379: static void quadratic_u_field_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
380: {
381:   uexact[0] = a[0];
382: }

384: static PetscErrorCode cross_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
385: {
386:   const PetscReal alpha = 50 * 4;
387:   const PetscReal xyz   = (x[0] - 0.5) * (x[1] - 0.5) * (x[2] - 0.5);

389:   *u = PetscSinReal(alpha * xyz) * (alpha * PetscAbsReal(xyz) < 2 * PETSC_PI ? (alpha * PetscAbsReal(xyz) > -2 * PETSC_PI ? 1.0 : 0.01) : 0.01);
390:   return PETSC_SUCCESS;
391: }

393: static void f0_cross_u_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
394: {
395:   const PetscReal alpha = 50 * 4;
396:   const PetscReal xyz   = (x[0] - 0.5) * (x[1] - 0.5) * (x[2] - 0.5);

398:   f0[0] = PetscSinReal(alpha * xyz) * (alpha * PetscAbsReal(xyz) < 2 * PETSC_PI ? (alpha * PetscAbsReal(xyz) > -2 * PETSC_PI ? 1.0 : 0.01) : 0.01);
399: }

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

406: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
407: {
408:   const char *bcTypes[3]    = {"neumann", "dirichlet", "none"};
409:   const char *runTypes[4]   = {"full", "exact", "test", "perf"};
410:   const char *coeffTypes[8] = {"none", "analytic", "field", "nonlinear", "ball", "cross", "checkerboard_0", "checkerboard_1"};
411:   PetscInt    bc, run, coeff;

413:   PetscFunctionBeginUser;
414:   options->runType             = RUN_FULL;
415:   options->bcType              = DIRICHLET;
416:   options->variableCoefficient = COEFF_NONE;
417:   options->fieldBC             = PETSC_FALSE;
418:   options->jacobianMF          = PETSC_FALSE;
419:   options->showInitial         = PETSC_FALSE;
420:   options->showSolution        = PETSC_FALSE;
421:   options->restart             = PETSC_FALSE;
422:   options->quiet               = PETSC_FALSE;
423:   options->nonzInit            = PETSC_FALSE;
424:   options->bdIntegral          = PETSC_FALSE;
425:   options->checkksp            = PETSC_FALSE;
426:   options->div                 = 4;
427:   options->k                   = 1;
428:   options->kgrid               = NULL;
429:   options->rand                = PETSC_FALSE;

431:   PetscOptionsBegin(comm, "", "Poisson Problem Options", "DMPLEX");
432:   run = options->runType;
433:   PetscCall(PetscOptionsEList("-run_type", "The run type", "ex12.c", runTypes, 4, runTypes[options->runType], &run, NULL));
434:   options->runType = (RunType)run;
435:   bc               = options->bcType;
436:   PetscCall(PetscOptionsEList("-bc_type", "Type of boundary condition", "ex12.c", bcTypes, 3, bcTypes[options->bcType], &bc, NULL));
437:   options->bcType = (BCType)bc;
438:   coeff           = options->variableCoefficient;
439:   PetscCall(PetscOptionsEList("-variable_coefficient", "Type of variable coefficient", "ex12.c", coeffTypes, 8, coeffTypes[options->variableCoefficient], &coeff, NULL));
440:   options->variableCoefficient = (CoeffType)coeff;

442:   PetscCall(PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL));
443:   PetscCall(PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL));
444:   PetscCall(PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL));
445:   PetscCall(PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL));
446:   PetscCall(PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL));
447:   PetscCall(PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL));
448:   PetscCall(PetscOptionsBool("-nonzero_initial_guess", "nonzero initial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL));
449:   PetscCall(PetscOptionsBool("-bd_integral", "Compute the integral of the solution on the boundary", "ex12.c", options->bdIntegral, &options->bdIntegral, NULL));
450:   if (options->runType == RUN_TEST) PetscCall(PetscOptionsBool("-run_test_check_ksp", "Check solution of KSP", "ex12.c", options->checkksp, &options->checkksp, NULL));
451:   PetscCall(PetscOptionsInt("-div", "The number of division for the checkerboard coefficient", "ex12.c", options->div, &options->div, NULL));
452:   PetscCall(PetscOptionsInt("-k", "The exponent for the checkerboard coefficient", "ex12.c", options->k, &options->k, NULL));
453:   PetscCall(PetscOptionsBool("-k_random", "Assign random k values to checkerboard", "ex12.c", options->rand, &options->rand, NULL));
454:   PetscOptionsEnd();
455:   PetscFunctionReturn(PETSC_SUCCESS);
456: }

458: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
459: {
460:   DM      plex;
461:   DMLabel label;

463:   PetscFunctionBeginUser;
464:   PetscCall(DMCreateLabel(dm, name));
465:   PetscCall(DMGetLabel(dm, name, &label));
466:   PetscCall(DMConvert(dm, DMPLEX, &plex));
467:   PetscCall(DMPlexMarkBoundaryFaces(plex, 1, label));
468:   PetscCall(DMDestroy(&plex));
469:   PetscFunctionReturn(PETSC_SUCCESS);
470: }

472: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
473: {
474:   PetscFunctionBeginUser;
475:   PetscCall(DMCreate(comm, dm));
476:   PetscCall(DMSetType(*dm, DMPLEX));
477:   PetscCall(DMSetFromOptions(*dm));
478:   {
479:     char      convType[256];
480:     PetscBool flg;

482:     PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
483:     PetscCall(PetscOptionsFList("-dm_plex_convert_type", "Convert DMPlex to another format", "ex12", DMList, DMPLEX, convType, 256, &flg));
484:     PetscOptionsEnd();
485:     if (flg) {
486:       DM dmConv;

488:       PetscCall(DMConvert(*dm, convType, &dmConv));
489:       if (dmConv) {
490:         PetscCall(DMDestroy(dm));
491:         *dm = dmConv;
492:       }
493:       PetscCall(DMSetFromOptions(*dm));
494:       PetscCall(DMSetUp(*dm));
495:     }
496:   }
497:   PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
498:   if (user->rand) {
499:     PetscRandom r;
500:     PetscReal   val;
501:     PetscInt    dim, N, i;

503:     PetscCall(DMGetDimension(*dm, &dim));
504:     N = PetscPowInt(user->div, dim);
505:     PetscCall(PetscMalloc1(N, &user->kgrid));
506:     PetscCall(PetscRandomCreate(PETSC_COMM_SELF, &r));
507:     PetscCall(PetscRandomSetFromOptions(r));
508:     PetscCall(PetscRandomSetInterval(r, 0.0, user->k));
509:     PetscCall(PetscRandomSetSeed(r, 1973));
510:     PetscCall(PetscRandomSeed(r));
511:     for (i = 0; i < N; ++i) {
512:       PetscCall(PetscRandomGetValueReal(r, &val));
513:       user->kgrid[i] = 1 + (PetscInt)val;
514:     }
515:     PetscCall(PetscRandomDestroy(&r));
516:   }
517:   PetscFunctionReturn(PETSC_SUCCESS);
518: }

520: static PetscErrorCode SetupProblem(DM dm, AppCtx *user)
521: {
522:   PetscDS          ds;
523:   DMLabel          label;
524:   PetscWeakForm    wf;
525:   const PetscReal *L;
526:   const PetscInt   id = 1;
527:   PetscInt         bd, dim;

529:   PetscFunctionBeginUser;
530:   PetscCall(DMGetDS(dm, &ds));
531:   PetscCall(DMGetDimension(dm, &dim));
532:   PetscCall(DMGetPeriodicity(dm, NULL, NULL, &L));
533:   switch (user->variableCoefficient) {
534:   case COEFF_NONE:
535:     if (L && L[0]) {
536:       if (L && L[1]) {
537:         PetscCall(PetscDSSetResidual(ds, 0, f0_xytrig_u, f1_u));
538:         PetscCall(PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu));
539:       } else {
540:         PetscCall(PetscDSSetResidual(ds, 0, f0_xtrig_u, f1_u));
541:         PetscCall(PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu));
542:       }
543:     } else {
544:       PetscCall(PetscDSSetResidual(ds, 0, f0_u, f1_u));
545:       PetscCall(PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu));
546:     }
547:     break;
548:   case COEFF_ANALYTIC:
549:     PetscCall(PetscDSSetResidual(ds, 0, f0_analytic_u, f1_analytic_u));
550:     PetscCall(PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_analytic_uu));
551:     break;
552:   case COEFF_FIELD:
553:     PetscCall(PetscDSSetResidual(ds, 0, f0_analytic_u, f1_field_u));
554:     PetscCall(PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_field_uu));
555:     break;
556:   case COEFF_NONLINEAR:
557:     PetscCall(PetscDSSetResidual(ds, 0, f0_analytic_nonlinear_u, f1_analytic_nonlinear_u));
558:     PetscCall(PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_analytic_nonlinear_uu));
559:     break;
560:   case COEFF_BALL:
561:     PetscCall(PetscDSSetResidual(ds, 0, f0_ball_u, f1_u));
562:     PetscCall(PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu));
563:     break;
564:   case COEFF_CROSS:
565:     switch (dim) {
566:     case 2:
567:       PetscCall(PetscDSSetResidual(ds, 0, f0_cross_u_2d, f1_u));
568:       break;
569:     case 3:
570:       PetscCall(PetscDSSetResidual(ds, 0, f0_cross_u_3d, f1_u));
571:       break;
572:     default:
573:       SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %" PetscInt_FMT, dim);
574:     }
575:     PetscCall(PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu));
576:     break;
577:   case COEFF_CHECKERBOARD_0:
578:     PetscCall(PetscDSSetResidual(ds, 0, f0_checkerboard_0_u, f1_field_u));
579:     PetscCall(PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_field_uu));
580:     break;
581:   default:
582:     SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid variable coefficient type %d", user->variableCoefficient);
583:   }
584:   switch (dim) {
585:   case 2:
586:     switch (user->variableCoefficient) {
587:     case COEFF_BALL:
588:       user->exactFuncs[0] = ball_u_2d;
589:       break;
590:     case COEFF_CROSS:
591:       user->exactFuncs[0] = cross_u_2d;
592:       break;
593:     case COEFF_CHECKERBOARD_0:
594:       user->exactFuncs[0] = zero;
595:       break;
596:     default:
597:       if (L && L[0]) {
598:         if (L && L[1]) {
599:           user->exactFuncs[0] = xytrig_u_2d;
600:         } else {
601:           user->exactFuncs[0] = xtrig_u_2d;
602:         }
603:       } else {
604:         user->exactFuncs[0]  = quadratic_u_2d;
605:         user->exactFields[0] = quadratic_u_field_2d;
606:       }
607:     }
608:     if (user->bcType == NEUMANN) {
609:       PetscCall(DMGetLabel(dm, "boundary", &label));
610:       PetscCall(DMAddBoundary(dm, DM_BC_NATURAL, "wall", label, 1, &id, 0, 0, NULL, NULL, NULL, user, &bd));
611:       PetscCall(PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL));
612:       PetscCall(PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, 0, f0_bd_u, 0, NULL));
613:     }
614:     break;
615:   case 3:
616:     switch (user->variableCoefficient) {
617:     case COEFF_BALL:
618:       user->exactFuncs[0] = ball_u_3d;
619:       break;
620:     case COEFF_CROSS:
621:       user->exactFuncs[0] = cross_u_3d;
622:       break;
623:     default:
624:       user->exactFuncs[0]  = quadratic_u_3d;
625:       user->exactFields[0] = quadratic_u_field_3d;
626:     }
627:     if (user->bcType == NEUMANN) {
628:       PetscCall(DMGetLabel(dm, "boundary", &label));
629:       PetscCall(DMAddBoundary(dm, DM_BC_NATURAL, "wall", label, 1, &id, 0, 0, NULL, NULL, NULL, user, &bd));
630:       PetscCall(PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL));
631:       PetscCall(PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, 0, f0_bd_u, 0, NULL));
632:     }
633:     break;
634:   default:
635:     SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %" PetscInt_FMT, dim);
636:   }
637:   /* Setup constants */
638:   switch (user->variableCoefficient) {
639:   case COEFF_CHECKERBOARD_0: {
640:     PetscScalar constants[2];

642:     constants[0] = user->div;
643:     constants[1] = user->k;
644:     PetscCall(PetscDSSetConstants(ds, 2, constants));
645:   } break;
646:   default:
647:     break;
648:   }
649:   PetscCall(PetscDSSetExactSolution(ds, 0, user->exactFuncs[0], user));
650:   /* Setup Boundary Conditions */
651:   if (user->bcType == DIRICHLET) {
652:     PetscCall(DMGetLabel(dm, "marker", &label));
653:     if (!label) {
654:       /* Right now, p4est cannot create labels immediately */
655:       PetscCall(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));
656:     } else {
657:       PetscCall(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));
658:     }
659:   }
660:   PetscFunctionReturn(PETSC_SUCCESS);
661: }

663: static PetscErrorCode SetupMaterial(DM dm, DM dmAux, AppCtx *user)
664: {
665:   PetscErrorCode (*matFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {nu_2d};
666:   void *ctx[1];
667:   Vec   nu;

669:   PetscFunctionBegin;
670:   ctx[0] = user;
671:   if (user->variableCoefficient == COEFF_CHECKERBOARD_0) matFuncs[0] = checkerboardCoeff;
672:   PetscCall(DMCreateLocalVector(dmAux, &nu));
673:   PetscCall(PetscObjectSetName((PetscObject)nu, "Coefficient"));
674:   PetscCall(DMProjectFunctionLocal(dmAux, 0.0, matFuncs, ctx, INSERT_ALL_VALUES, nu));
675:   PetscCall(DMSetAuxiliaryVec(dm, NULL, 0, 0, nu));
676:   PetscCall(VecDestroy(&nu));
677:   PetscFunctionReturn(PETSC_SUCCESS);
678: }

680: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
681: {
682:   PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
683:   Vec      uexact;
684:   PetscInt dim;

686:   PetscFunctionBegin;
687:   PetscCall(DMGetDimension(dm, &dim));
688:   if (dim == 2) bcFuncs[0] = quadratic_u_2d;
689:   else bcFuncs[0] = quadratic_u_3d;
690:   PetscCall(DMCreateLocalVector(dmAux, &uexact));
691:   PetscCall(DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact));
692:   PetscCall(DMSetAuxiliaryVec(dm, NULL, 0, 0, uexact));
693:   PetscCall(VecDestroy(&uexact));
694:   PetscFunctionReturn(PETSC_SUCCESS);
695: }

697: static PetscErrorCode SetupAuxDM(DM dm, PetscFE feAux, AppCtx *user)
698: {
699:   DM dmAux, coordDM;

701:   PetscFunctionBegin;
702:   /* MUST call DMGetCoordinateDM() in order to get p4est setup if present */
703:   PetscCall(DMGetCoordinateDM(dm, &coordDM));
704:   if (!feAux) PetscFunctionReturn(PETSC_SUCCESS);
705:   PetscCall(DMClone(dm, &dmAux));
706:   PetscCall(DMSetCoordinateDM(dmAux, coordDM));
707:   PetscCall(DMSetField(dmAux, 0, NULL, (PetscObject)feAux));
708:   PetscCall(DMCreateDS(dmAux));
709:   if (user->fieldBC) PetscCall(SetupBC(dm, dmAux, user));
710:   else PetscCall(SetupMaterial(dm, dmAux, user));
711:   PetscCall(DMDestroy(&dmAux));
712:   PetscFunctionReturn(PETSC_SUCCESS);
713: }

715: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
716: {
717:   DM        plex, cdm = dm;
718:   PetscFE   fe, feAux = NULL;
719:   PetscBool simplex;
720:   PetscInt  dim;

722:   PetscFunctionBeginUser;
723:   PetscCall(DMGetDimension(dm, &dim));
724:   PetscCall(DMConvert(dm, DMPLEX, &plex));
725:   PetscCall(DMPlexIsSimplex(plex, &simplex));
726:   PetscCall(DMDestroy(&plex));
727:   PetscCall(PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, NULL, -1, &fe));
728:   PetscCall(PetscObjectSetName((PetscObject)fe, "potential"));
729:   if (user->variableCoefficient == COEFF_FIELD || user->variableCoefficient == COEFF_CHECKERBOARD_0) {
730:     PetscCall(PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, "mat_", -1, &feAux));
731:     PetscCall(PetscObjectSetName((PetscObject)feAux, "coefficient"));
732:     PetscCall(PetscFECopyQuadrature(fe, feAux));
733:   } else if (user->fieldBC) {
734:     PetscCall(PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, "bc_", -1, &feAux));
735:     PetscCall(PetscFECopyQuadrature(fe, feAux));
736:   }
737:   /* Set discretization and boundary conditions for each mesh */
738:   PetscCall(DMSetField(dm, 0, NULL, (PetscObject)fe));
739:   PetscCall(DMCreateDS(dm));
740:   PetscCall(SetupProblem(dm, user));
741:   while (cdm) {
742:     PetscCall(SetupAuxDM(cdm, feAux, user));
743:     if (user->bcType == DIRICHLET) {
744:       PetscBool hasLabel;

746:       PetscCall(DMHasLabel(cdm, "marker", &hasLabel));
747:       if (!hasLabel) PetscCall(CreateBCLabel(cdm, "marker"));
748:     }
749:     PetscCall(DMCopyDisc(dm, cdm));
750:     PetscCall(DMGetCoarseDM(cdm, &cdm));
751:   }
752:   PetscCall(PetscFEDestroy(&fe));
753:   PetscCall(PetscFEDestroy(&feAux));
754:   PetscFunctionReturn(PETSC_SUCCESS);
755: }

757: int main(int argc, char **argv)
758: {
759:   DM           dm;          /* Problem specification */
760:   SNES         snes;        /* nonlinear solver */
761:   Vec          u;           /* solution vector */
762:   Mat          A, J;        /* Jacobian matrix */
763:   MatNullSpace nullSpace;   /* May be necessary for Neumann conditions */
764:   AppCtx       user;        /* user-defined work context */
765:   JacActionCtx userJ;       /* context for Jacobian MF action */
766:   PetscReal    error = 0.0; /* L_2 error in the solution */

768:   PetscFunctionBeginUser;
769:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
770:   PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
771:   PetscCall(SNESCreate(PETSC_COMM_WORLD, &snes));
772:   PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm));
773:   PetscCall(SNESSetDM(snes, dm));
774:   PetscCall(DMSetApplicationContext(dm, &user));

776:   PetscCall(PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields));
777:   PetscCall(SetupDiscretization(dm, &user));

779:   PetscCall(DMCreateGlobalVector(dm, &u));
780:   PetscCall(PetscObjectSetName((PetscObject)u, "potential"));

782:   PetscCall(DMCreateMatrix(dm, &J));
783:   if (user.jacobianMF) {
784:     PetscInt M, m, N, n;

786:     PetscCall(MatGetSize(J, &M, &N));
787:     PetscCall(MatGetLocalSize(J, &m, &n));
788:     PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
789:     PetscCall(MatSetSizes(A, m, n, M, N));
790:     PetscCall(MatSetType(A, MATSHELL));
791:     PetscCall(MatSetUp(A));
792: #if 0
793:     PetscCall(MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction));
794: #endif

796:     userJ.dm   = dm;
797:     userJ.J    = J;
798:     userJ.user = &user;

800:     PetscCall(DMCreateLocalVector(dm, &userJ.u));
801:     if (user.fieldBC) PetscCall(DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u));
802:     else PetscCall(DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u));
803:     PetscCall(MatShellSetContext(A, &userJ));
804:   } else {
805:     A = J;
806:   }

808:   nullSpace = NULL;
809:   if (user.bcType != DIRICHLET) {
810:     PetscCall(MatNullSpaceCreate(PetscObjectComm((PetscObject)dm), PETSC_TRUE, 0, NULL, &nullSpace));
811:     PetscCall(MatSetNullSpace(A, nullSpace));
812:   }

814:   PetscCall(DMPlexSetSNESLocalFEM(dm, PETSC_FALSE, &user));
815:   PetscCall(SNESSetJacobian(snes, A, J, NULL, NULL));

817:   PetscCall(SNESSetFromOptions(snes));

819:   if (user.fieldBC) PetscCall(DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u));
820:   else PetscCall(DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u));
821:   if (user.restart) {
822: #if defined(PETSC_HAVE_HDF5)
823:     PetscViewer viewer;
824:     char        filename[PETSC_MAX_PATH_LEN];

826:     PetscCall(PetscOptionsGetString(NULL, NULL, "-dm_plex_filename", filename, sizeof(filename), NULL));
827:     PetscCall(PetscViewerCreate(PETSC_COMM_WORLD, &viewer));
828:     PetscCall(PetscViewerSetType(viewer, PETSCVIEWERHDF5));
829:     PetscCall(PetscViewerFileSetMode(viewer, FILE_MODE_READ));
830:     PetscCall(PetscViewerFileSetName(viewer, filename));
831:     PetscCall(PetscViewerHDF5PushGroup(viewer, "/fields"));
832:     PetscCall(VecLoad(u, viewer));
833:     PetscCall(PetscViewerHDF5PopGroup(viewer));
834:     PetscCall(PetscViewerDestroy(&viewer));
835: #endif
836:   }
837:   if (user.showInitial) {
838:     Vec lv;
839:     PetscCall(DMGetLocalVector(dm, &lv));
840:     PetscCall(DMGlobalToLocalBegin(dm, u, INSERT_VALUES, lv));
841:     PetscCall(DMGlobalToLocalEnd(dm, u, INSERT_VALUES, lv));
842:     PetscCall(DMPrintLocalVec(dm, "Local function", 1.0e-10, lv));
843:     PetscCall(DMRestoreLocalVector(dm, &lv));
844:   }
845:   if (user.runType == RUN_FULL || user.runType == RUN_EXACT) {
846:     PetscErrorCode (*initialGuess[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {zero};

848:     if (user.nonzInit) initialGuess[0] = ecks;
849:     if (user.runType == RUN_FULL) PetscCall(DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u));
850:     PetscCall(VecViewFromOptions(u, NULL, "-guess_vec_view"));
851:     PetscCall(SNESSolve(snes, NULL, u));
852:     PetscCall(SNESGetSolution(snes, &u));
853:     PetscCall(SNESGetDM(snes, &dm));

855:     if (user.showSolution) {
856:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Solution\n"));
857:       PetscCall(VecFilter(u, 3.0e-9));
858:       PetscCall(VecView(u, PETSC_VIEWER_STDOUT_WORLD));
859:     }
860:   } else if (user.runType == RUN_PERF) {
861:     Vec       r;
862:     PetscReal res = 0.0;

864:     PetscCall(SNESGetFunction(snes, &r, NULL, NULL));
865:     PetscCall(SNESComputeFunction(snes, u, r));
866:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n"));
867:     PetscCall(VecFilter(r, 1.0e-10));
868:     PetscCall(VecNorm(r, NORM_2, &res));
869:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res));
870:   } else {
871:     Vec       r;
872:     PetscReal res = 0.0, tol = 1.0e-11;

874:     /* Check discretization error */
875:     PetscCall(SNESGetFunction(snes, &r, NULL, NULL));
876:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n"));
877:     if (!user.quiet) PetscCall(VecView(u, PETSC_VIEWER_STDOUT_WORLD));
878:     PetscCall(DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error));
879:     if (error < tol) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < %2.1e\n", (double)tol));
880:     else PetscCall(PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", (double)error));
881:     /* Check residual */
882:     PetscCall(SNESComputeFunction(snes, u, r));
883:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n"));
884:     PetscCall(VecFilter(r, 1.0e-10));
885:     if (!user.quiet) PetscCall(VecView(r, PETSC_VIEWER_STDOUT_WORLD));
886:     PetscCall(VecNorm(r, NORM_2, &res));
887:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res));
888:     /* Check Jacobian */
889:     {
890:       Vec b;

892:       PetscCall(SNESComputeJacobian(snes, u, A, A));
893:       PetscCall(VecDuplicate(u, &b));
894:       PetscCall(VecSet(r, 0.0));
895:       PetscCall(SNESComputeFunction(snes, r, b));
896:       PetscCall(MatMult(A, u, r));
897:       PetscCall(VecAXPY(r, 1.0, b));
898:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n"));
899:       PetscCall(VecFilter(r, 1.0e-10));
900:       if (!user.quiet) PetscCall(VecView(r, PETSC_VIEWER_STDOUT_WORLD));
901:       PetscCall(VecNorm(r, NORM_2, &res));
902:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", (double)res));
903:       /* check solver */
904:       if (user.checkksp) {
905:         KSP ksp;

907:         if (nullSpace) PetscCall(MatNullSpaceRemove(nullSpace, u));
908:         PetscCall(SNESComputeJacobian(snes, u, A, J));
909:         PetscCall(MatMult(A, u, b));
910:         PetscCall(SNESGetKSP(snes, &ksp));
911:         PetscCall(KSPSetOperators(ksp, A, J));
912:         PetscCall(KSPSolve(ksp, b, r));
913:         PetscCall(VecAXPY(r, -1.0, u));
914:         PetscCall(VecNorm(r, NORM_2, &res));
915:         PetscCall(PetscPrintf(PETSC_COMM_WORLD, "KSP Error: %g\n", (double)res));
916:       }
917:       PetscCall(VecDestroy(&b));
918:     }
919:   }
920:   PetscCall(VecViewFromOptions(u, NULL, "-vec_view"));
921:   {
922:     Vec nu;

924:     PetscCall(DMGetAuxiliaryVec(dm, NULL, 0, 0, &nu));
925:     if (nu) PetscCall(VecViewFromOptions(nu, NULL, "-coeff_view"));
926:   }

928:   if (user.bdIntegral) {
929:     DMLabel          label;
930:     PetscBdPointFunc func[1] = {bd_integral_2d};
931:     PetscInt         id      = 1;
932:     PetscScalar      bdInt   = 0.0;
933:     PetscReal        exact   = 3.3333333333;

935:     PetscCall(DMGetLabel(dm, "marker", &label));
936:     PetscCall(DMPlexComputeBdIntegral(dm, u, label, 1, &id, func, &bdInt, NULL));
937:     PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Solution boundary integral: %.4g\n", (double)PetscAbsScalar(bdInt)));
938:     PetscCheck(PetscAbsReal(PetscAbsScalar(bdInt) - exact) <= PETSC_SQRT_MACHINE_EPSILON, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Invalid boundary integral %g != %g", (double)PetscAbsScalar(bdInt), (double)exact);
939:   }

941:   PetscCall(MatNullSpaceDestroy(&nullSpace));
942:   if (user.jacobianMF) PetscCall(VecDestroy(&userJ.u));
943:   if (A != J) PetscCall(MatDestroy(&A));
944:   PetscCall(MatDestroy(&J));
945:   PetscCall(VecDestroy(&u));
946:   PetscCall(SNESDestroy(&snes));
947:   PetscCall(DMDestroy(&dm));
948:   PetscCall(PetscFree2(user.exactFuncs, user.exactFields));
949:   PetscCall(PetscFree(user.kgrid));
950:   PetscCall(PetscFinalize());
951:   return 0;
952: }

954: /*TEST
955:   # 2D serial P1 test 0-4
956:   test:
957:     suffix: 2d_p1_0
958:     requires: triangle
959:     args: -run_type test -bc_type dirichlet -dm_plex_interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

961:   test:
962:     suffix: 2d_p1_1
963:     requires: triangle
964:     args: -run_type test -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cdm_dm_plex_coordinate_dim {{2 3}}

966:   test:
967:     suffix: 2d_p1_1b
968:     requires: triangle
969:     args: -run_type test -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_refine 3 -dm_coord_space 0 \
970:           -dm_plex_option_phases proj_ -cdm_proj_dm_plex_coordinate_dim 3 -proj_dm_coord_space \
971:           -proj_dm_coord_remap -proj_dm_coord_map sinusoid -proj_dm_coord_map_params 0.1,1.,1.

973:   test:
974:     suffix: 2d_p1_2
975:     requires: triangle
976:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

978:   test:
979:     suffix: 2d_p1_neumann_0
980:     requires: triangle
981:     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

983:   test:
984:     suffix: 2d_p1_neumann_1
985:     requires: triangle
986:     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

988:   # 2D serial P2 test 5-8
989:   test:
990:     suffix: 2d_p2_0
991:     requires: triangle
992:     args: -run_type test -bc_type dirichlet -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

994:   test:
995:     suffix: 2d_p2_1
996:     requires: triangle
997:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -bc_type dirichlet -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

999:   test:
1000:     suffix: 2d_p2_neumann_0
1001:     requires: triangle
1002:     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

1004:   test:
1005:     suffix: 2d_p2_neumann_1
1006:     requires: triangle
1007:     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

1009:   test:
1010:     suffix: bd_int_0
1011:     requires: triangle
1012:     args: -run_type test -bc_type dirichlet -petscspace_degree 2 -bd_integral -dm_view -quiet

1014:   test:
1015:     suffix: bd_int_1
1016:     requires: triangle
1017:     args: -run_type test -dm_refine 2 -bc_type dirichlet -petscspace_degree 2 -bd_integral -dm_view -quiet

1019:   # 3D serial P1 test 9-12
1020:   test:
1021:     suffix: 3d_p1_0
1022:     requires: ctetgen
1023:     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

1025:   test:
1026:     suffix: 3d_p1_1
1027:     requires: ctetgen
1028:     args: -run_type test -dm_plex_dim 3 -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view

1030:   test:
1031:     suffix: 3d_p1_2
1032:     requires: ctetgen
1033:     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

1035:   test:
1036:     suffix: 3d_p1_neumann_0
1037:     requires: ctetgen
1038:     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

1040:   # Analytic variable coefficient 13-20
1041:   test:
1042:     suffix: 13
1043:     requires: triangle
1044:     args: -run_type test -variable_coefficient analytic -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1045:   test:
1046:     suffix: 14
1047:     requires: triangle
1048:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient analytic -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1049:   test:
1050:     suffix: 15
1051:     requires: triangle
1052:     args: -run_type test -variable_coefficient analytic -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1053:   test:
1054:     suffix: 16
1055:     requires: triangle
1056:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient analytic -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1057:   test:
1058:     suffix: 17
1059:     requires: ctetgen
1060:     args: -run_type test -dm_plex_dim 3 -variable_coefficient analytic -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1062:   test:
1063:     suffix: 18
1064:     requires: ctetgen
1065:     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

1067:   test:
1068:     suffix: 19
1069:     requires: ctetgen
1070:     args: -run_type test -dm_plex_dim 3 -variable_coefficient analytic -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1072:   test:
1073:     suffix: 20
1074:     requires: ctetgen
1075:     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

1077:   # P1 variable coefficient 21-28
1078:   test:
1079:     suffix: 21
1080:     requires: triangle
1081:     args: -run_type test -variable_coefficient field -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1083:   test:
1084:     suffix: 22
1085:     requires: triangle
1086:     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

1088:   test:
1089:     suffix: 23
1090:     requires: triangle
1091:     args: -run_type test -variable_coefficient field -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1093:   test:
1094:     suffix: 24
1095:     requires: triangle
1096:     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

1098:   test:
1099:     suffix: 25
1100:     requires: ctetgen
1101:     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

1103:   test:
1104:     suffix: 26
1105:     requires: ctetgen
1106:     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

1108:   test:
1109:     suffix: 27
1110:     requires: ctetgen
1111:     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

1113:   test:
1114:     suffix: 28
1115:     requires: ctetgen
1116:     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

1118:   # P0 variable coefficient 29-36
1119:   test:
1120:     suffix: 29
1121:     requires: triangle
1122:     args: -run_type test -variable_coefficient field -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1124:   test:
1125:     suffix: 30
1126:     requires: triangle
1127:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient field -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1129:   test:
1130:     suffix: 31
1131:     requires: triangle
1132:     args: -run_type test -variable_coefficient field -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1134:   test:
1135:     requires: triangle
1136:     suffix: 32
1137:     args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient field -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1139:   test:
1140:     requires: ctetgen
1141:     suffix: 33
1142:     args: -run_type test -dm_plex_dim 3 -variable_coefficient field -petscspace_degree 1 -show_initial -dm_plex_print_fem 1

1144:   test:
1145:     suffix: 34
1146:     requires: ctetgen
1147:     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

1149:   test:
1150:     suffix: 35
1151:     requires: ctetgen
1152:     args: -run_type test -dm_plex_dim 3 -variable_coefficient field -petscspace_degree 2 -show_initial -dm_plex_print_fem 1

1154:   test:
1155:     suffix: 36
1156:     requires: ctetgen
1157:     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

1159:   # Full solve 39-44
1160:   test:
1161:     suffix: 39
1162:     requires: triangle !single
1163:     args: -run_type full -dm_refine_volume_limit_pre 0.015625 -petscspace_degree 2 -pc_type gamg -pc_gamg_esteig_ksp_type cg -pc_gamg_esteig_ksp_max_it 10 -snes_rtol 1.0e-6 -ksp_rtol 1.0e-7 -ksp_monitor -ksp_converged_reason -snes_monitor_short -snes_converged_reason ::ascii_info_detail
1164:   test:
1165:     suffix: 40
1166:     requires: triangle !single
1167:     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
1168:   test:
1169:     suffix: 41
1170:     requires: triangle !single
1171:     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
1172:   test:
1173:     suffix: 42
1174:     requires: triangle !single
1175:     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
1176:   test:
1177:     suffix: 43
1178:     requires: triangle !single
1179:     nsize: 2
1180:     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

1182:   test:
1183:     suffix: 44
1184:     requires: triangle !single
1185:     nsize: 2
1186:     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 -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

1188:   # These tests use a loose tolerance just to exercise the PtAP operations for MATIS and multiple PCBDDC setup calls inside PCMG
1189:   testset:
1190:     requires: triangle !single
1191:     nsize: 3
1192:     args: -run_type full -petscspace_degree 1 -dm_mat_type is -pc_type mg -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
1193:     test:
1194:       suffix: gmg_bddc
1195:       filter: sed -e "s/CONVERGED_FNORM_RELATIVE iterations 3/CONVERGED_FNORM_RELATIVE iterations 4/g"
1196:       args: -mg_levels_pc_type jacobi
1197:     test:
1198:       filter: sed -e "s/iterations [0-4]/iterations 4/g"
1199:       suffix: gmg_bddc_lev
1200:       args: -mg_levels_pc_type bddc

1202:   # VTU viewer with empty processes
1203:   test:
1204:     requires: !complex
1205:     suffix: vtu_empty
1206:     args: -quiet -run_type test -dm_plex_simplex 0 -dm_plex_box_faces 2,2 -vec_view vtk:test.vtu:vtk_vtu -petscspace_degree 1 -petscpartitioner_type simple

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 -dm_plex_name box -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 -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 -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_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: -mat_is_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: -mat_is_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:   # These two tests are broken because DMPlexComputeInjectorFEM() only works for regularly refined meshes
1390:   test:
1391:     suffix: fas_newton_coarse_0
1392:     requires: pragmatic triangle
1393:     TODO: broken
1394:     args: -run_type full -variable_coefficient nonlinear -petscspace_degree 1 \
1395:           -dm_refine 2 -dm_coarsen_hierarchy 1 -dm_plex_hash_location -dm_adaptor pragmatic \
1396:           -snes_type fas -snes_fas_levels 2 -snes_converged_reason ::ascii_info_detail -snes_monitor_short -snes_view \
1397:             -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -fas_coarse_snes_linesearch_type basic \
1398:             -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

1400:   test:
1401:     suffix: mg_newton_coarse_0
1402:     requires: triangle pragmatic
1403:     TODO: broken
1404:     args: -run_type full -petscspace_degree 1 \
1405:           -dm_refine 3 -dm_coarsen_hierarchy 3 -dm_plex_hash_location -dm_adaptor pragmatic \
1406:           -snes_atol 1.0e-8 -snes_rtol 0.0 -snes_monitor_short -snes_converged_reason ::ascii_info_detail -snes_view \
1407:             -ksp_type richardson -ksp_atol 1.0e-8 -ksp_rtol 0.0 -ksp_norm_type unpreconditioned -ksp_monitor_true_residual \
1408:               -pc_type mg -pc_mg_levels 4 \
1409:               -mg_levels_ksp_type gmres -mg_levels_pc_type ilu -mg_levels_ksp_max_it 10

1411:   # Test cgns writer for ranks with no elements
1412:   test:
1413:     suffix: cgns
1414:     nsize: 5
1415:     requires: cgns
1416:     args: -quiet -run_type test -dm_plex_simplex 0 -petscspace_degree 1 -dm_plex_box_faces 2,2 -vec_view cgns:test.cgns -dm_refine 0 -petscpartitioner_type simple

1418:   # Full solve tensor
1419:   test:
1420:     suffix: tensor_plex_2d
1421:     args: -run_type test -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_refine_hierarchy 2

1423:   test:
1424:     suffix: tensor_p4est_2d
1425:     requires: p4est
1426:     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

1428:   test:
1429:     suffix: tensor_plex_3d
1430:     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

1432:   test:
1433:     suffix: tensor_p4est_3d
1434:     requires: p4est
1435:     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

1437:   test:
1438:     suffix: p4est_test_q2_conformal_serial
1439:     requires: p4est
1440:     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

1442:   test:
1443:     suffix: p4est_test_q2_conformal_parallel
1444:     requires: p4est
1445:     nsize: 7
1446:     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 -petscpartitioner_type simple

1448:   test:
1449:     suffix: p4est_test_q2_conformal_parallel_parmetis
1450:     requires: parmetis p4est
1451:     nsize: 4
1452:     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 -petscpartitioner_type parmetis

1454:   test:
1455:     suffix: p4est_test_q2_nonconformal_serial
1456:     requires: p4est
1457:     filter: grep -v "CG or CGNE: variant"
1458:     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

1460:   test:
1461:     suffix: p4est_test_q2_nonconformal_parallel
1462:     requires: p4est
1463:     filter: grep -v "CG or CGNE: variant"
1464:     nsize: 7
1465:     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 -petscpartitioner_type simple

1467:   test:
1468:     suffix: p4est_test_q2_nonconformal_parallel_parmetis
1469:     requires: parmetis p4est
1470:     nsize: 4
1471:     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 -petscpartitioner_type parmetis

1473:   test:
1474:     suffix: p4est_exact_q2_conformal_serial
1475:     requires: p4est !single !complex !__float128
1476:     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

1478:   test:
1479:     suffix: p4est_exact_q2_conformal_parallel
1480:     requires: p4est !single !complex !__float128
1481:     nsize: 4
1482:     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

1484:   test:
1485:     suffix: p4est_exact_q2_conformal_parallel_parmetis
1486:     requires: parmetis p4est !single
1487:     nsize: 4
1488:     args: -run_type exact -petscspace_degree 2 -fas_levels_snes_linesearch_type basic -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_snes_converged_reason -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

1490:   test:
1491:     suffix: p4est_exact_q2_nonconformal_serial
1492:     requires: p4est
1493:     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

1495:   test:
1496:     suffix: p4est_exact_q2_nonconformal_parallel
1497:     requires: p4est
1498:     nsize: 7
1499:     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 -petscpartitioner_type simple

1501:   test:
1502:     suffix: p4est_exact_q2_nonconformal_parallel_parmetis
1503:     requires: parmetis p4est
1504:     nsize: 4
1505:     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 -petscpartitioner_type parmetis

1507:   test:
1508:     suffix: p4est_full_q2_nonconformal_serial
1509:     requires: p4est !single
1510:     filter: grep -v "variant HERMITIAN"
1511:     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

1513:   test:
1514:     suffix: p4est_full_q2_nonconformal_parallel
1515:     requires: p4est !single
1516:     filter: grep -v "variant HERMITIAN"
1517:     nsize: 7
1518:     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 -petscpartitioner_type simple

1520:   test:
1521:     suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1522:     requires: p4est !single
1523:     filter: grep -v "variant HERMITIAN"
1524:     nsize: 7
1525:     args: -run_type full -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

1527:   test:
1528:     suffix: p4est_full_q2_nonconformal_parallel_bddc
1529:     requires: p4est !single
1530:     filter: grep -v "variant HERMITIAN"
1531:     nsize: 7
1532:     args: -run_type full -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

1534:   test:
1535:     TODO: broken
1536:     suffix: p4est_fas_q2_conformal_serial
1537:     requires: p4est !complex !__float128
1538:     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

1540:   test:
1541:     TODO: broken
1542:     suffix: p4est_fas_q2_nonconformal_serial
1543:     requires: p4est
1544:     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

1546:   test:
1547:     suffix: fas_newton_0_p4est
1548:     requires: p4est !single !__float128
1549:     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

1551:   # Full solve simplicial AMR
1552:   test:
1553:     suffix: tri_p1_adapt_init_pragmatic
1554:     requires: pragmatic
1555:     args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_initial 1 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor pragmatic

1557:   test:
1558:     suffix: tri_p2_adapt_init_pragmatic
1559:     requires: pragmatic
1560:     args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_initial 1 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor pragmatic

1562:   test:
1563:     suffix: tri_p1_adapt_init_mmg
1564:     requires: mmg
1565:     args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_initial 1 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor mmg

1567:   test:
1568:     suffix: tri_p2_adapt_init_mmg
1569:     requires: mmg
1570:     args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_initial 1 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor mmg

1572:   test:
1573:     suffix: tri_p1_adapt_seq_pragmatic
1574:     requires: pragmatic
1575:     args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_sequence 2 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor pragmatic

1577:   test:
1578:     suffix: tri_p2_adapt_seq_pragmatic
1579:     requires: pragmatic
1580:     args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_sequence 2 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor pragmatic

1582:   test:
1583:     suffix: tri_p1_adapt_seq_mmg
1584:     requires: mmg
1585:     args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_sequence 2 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor mmg

1587:   test:
1588:     suffix: tri_p2_adapt_seq_mmg
1589:     requires: mmg
1590:     args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_sequence 2 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor mmg

1592:   test:
1593:     suffix: tri_p1_adapt_analytic_pragmatic
1594:     requires: pragmatic
1595:     args: -run_type exact -dm_refine 3 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -dm_plex_metric_h_min 0.0001 -dm_plex_metric_h_max 0.05 -dm_adaptor pragmatic

1597:   test:
1598:     suffix: tri_p2_adapt_analytic_pragmatic
1599:     requires: pragmatic
1600:     args: -run_type exact -dm_refine 3 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -dm_plex_metric_h_min 0.0001 -dm_plex_metric_h_max 0.05 -dm_adaptor pragmatic

1602:   test:
1603:     suffix: tri_p1_adapt_analytic_mmg
1604:     requires: mmg
1605:     args: -run_type exact -dm_refine 3 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -dm_plex_metric_h_max 0.5 -dm_adaptor mmg

1607:   test:
1608:     suffix: tri_p2_adapt_analytic_mmg
1609:     requires: mmg
1610:     args: -run_type exact -dm_refine 3 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -dm_plex_metric_h_max 0.5 -dm_adaptor mmg

1612:   test:
1613:     suffix: tri_p1_adapt_uniform_pragmatic
1614:     requires: pragmatic tetgen
1615:     nsize: 2
1616:     args: -run_type full -dm_plex_box_faces 8,8,8 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 3 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor pragmatic
1617:     timeoutfactor: 2

1619:   test:
1620:     suffix: tri_p2_adapt_uniform_pragmatic
1621:     requires: pragmatic tetgen
1622:     nsize: 2
1623:     args: -run_type full -dm_plex_box_faces 8,8,8 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 1 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor pragmatic
1624:     timeoutfactor: 1

1626:   test:
1627:     suffix: tri_p1_adapt_uniform_mmg
1628:     requires: mmg tetgen
1629:     args: -run_type full -dm_plex_box_faces 4,4,4 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 3 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor mmg
1630:     timeoutfactor: 2

1632:   test:
1633:     suffix: tri_p2_adapt_uniform_mmg
1634:     requires: mmg tetgen
1635:     TODO: broken
1636:     args: -run_type full -dm_plex_box_faces 4,4,4 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 1 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor mmg
1637:     timeoutfactor: 1

1639:   test:
1640:     suffix: tri_p1_adapt_uniform_parmmg
1641:     requires: parmmg tetgen
1642:     nsize: 2
1643:     args: -run_type full -dm_plex_box_faces 8,8,8 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 3 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor parmmg
1644:     timeoutfactor: 2

1646:   test:
1647:     suffix: tri_p2_adapt_uniform_parmmg
1648:     requires: parmmg tetgen
1649:     nsize: 2
1650:     args: -run_type full -dm_plex_box_faces 8,8,8 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 1 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor parmmg
1651:     timeoutfactor: 1

1653:   # Full solve tensor AMR
1654:   test:
1655:     suffix: quad_q1_adapt_0
1656:     requires: p4est
1657:     args: -run_type exact -dm_plex_simplex 0 -dm_plex_convert_type p4est -bc_type dirichlet -petscspace_degree 1 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -dm_forest_initial_refinement 4 -snes_adapt_initial 1 -dm_view
1658:     filter: grep -v DM_

1660:   test:
1661:     suffix: amr_0
1662:     nsize: 5
1663:     args: -run_type test -petscpartitioner_type simple -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_refine 1

1665:   test:
1666:     suffix: amr_1
1667:     requires: p4est !complex
1668:     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

1670:   test:
1671:     suffix: p4est_solve_bddc
1672:     requires: p4est !complex
1673:     args: -run_type full -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
1674:     nsize: 4

1676:   test:
1677:     suffix: p4est_solve_fas
1678:     requires: p4est
1679:     args: -run_type full -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
1680:     nsize: 4
1681:     TODO: identical machine two runs produce slightly different solver trackers

1683:   test:
1684:     suffix: p4est_convergence_test_1
1685:     requires: p4est
1686:     args: -quiet -run_type test -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
1687:     nsize: 4

1689:   # Serial tests with GLVis visualization
1690:   test:
1691:     suffix: glvis_2d_tet_p1
1692:     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
1693:   test:
1694:     suffix: glvis_2d_tet_p2
1695:     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
1696:   test:
1697:     suffix: glvis_2d_hex_p1
1698:     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
1699:   test:
1700:     suffix: glvis_2d_hex_p2
1701:     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
1702:   test:
1703:     suffix: glvis_2d_hex_p2_p4est
1704:     requires: p4est
1705:     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
1706:   test:
1707:     suffix: glvis_2d_tet_p0
1708:     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 -pc_type jacobi
1709:   test:
1710:     suffix: glvis_2d_hex_p0
1711:     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 -pc_type jacobi

1713:   # PCHPDDM tests
1714:   testset:
1715:     nsize: 4
1716:     requires: hpddm slepc !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1717:     args: -run_type test -run_test_check_ksp -quiet -petscspace_degree 1 -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
1718:     test:
1719:       suffix: quad_singular_hpddm
1720:       args: -dm_plex_box_faces 6,7
1721:     test:
1722:       requires: p4est
1723:       suffix: p4est_singular_2d_hpddm
1724:       args: -dm_plex_convert_type p4est -dm_forest_minimum_refinement 1 -dm_forest_initial_refinement 3 -dm_forest_maximum_refinement 3
1725:     test:
1726:       requires: p4est
1727:       suffix: p4est_nc_singular_2d_hpddm
1728:       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
1729:   testset:
1730:     nsize: 4
1731:     requires: hpddm slepc triangle !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1732:     args: -run_type full -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
1733:     test:
1734:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1735:       suffix: tri_hpddm_reuse_baij
1736:     test:
1737:       requires: !complex
1738:       suffix: tri_hpddm_reuse
1739:   testset:
1740:     nsize: 4
1741:     requires: hpddm slepc !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1742:     args: -run_type full -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
1743:     test:
1744:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1745:       suffix: quad_hpddm_reuse_baij
1746:     test:
1747:       requires: !complex
1748:       suffix: quad_hpddm_reuse
1749:   testset:
1750:     nsize: 4
1751:     requires: hpddm slepc !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1752:     args: -run_type full -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
1753:     test:
1754:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1755:       suffix: quad_hpddm_reuse_threshold_baij
1756:     test:
1757:       requires: !complex
1758:       suffix: quad_hpddm_reuse_threshold
1759:   testset:
1760:     nsize: 4
1761:     requires: hpddm slepc parmetis !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1762:     filter: sed -e "s/linear solver iterations=17/linear solver iterations=16/g"
1763:     args: -run_type full -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 ball -dm_plex_gmsh_periodic 0 -fp_trap 0
1764:     test:
1765:       args: -pc_hpddm_coarse_mat_type baij -options_left no
1766:       filter: grep -v "      total: nonzeros=" | grep -v "      rows=" | sed -e "s/total number of linear solver iterations=[1-2][4-7]/total number of linear solver iterations=16/g"
1767:       suffix: tri_parmetis_hpddm_baij
1768:     test:
1769:       filter: grep -v "      total: nonzeros=" | grep -v "      rows=" | sed -e "s/total number of linear solver iterations=[1-2][4-7]/total number of linear solver iterations=16/g"
1770:       requires: !complex
1771:       suffix: tri_parmetis_hpddm

1773:   # 2D serial P1 tests for adaptive MG
1774:   test:
1775:     suffix: 2d_p1_adaptmg_0
1776:     requires: triangle
1777:     args: -petscpartitioner_type simple -dm_refine_hierarchy 3 -dm_plex_box_faces 4,4 -bc_type dirichlet -petscspace_degree 1 \
1778:           -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
1779:           -snes_max_it 1 -ksp_converged_reason \
1780:           -ksp_rtol 1e-8 -pc_type mg
1781:   test:
1782:     suffix: 2d_p1_adaptmg_1
1783:     TODO: broken
1784:     requires: triangle bamg
1785:     args: -petscpartitioner_type simple -dm_refine_hierarchy 3 -dm_plex_box_faces 4,4 -bc_type dirichlet -petscspace_degree 1 \
1786:           -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
1787:           -snes_max_it 1 -ksp_converged_reason \
1788:           -ksp_rtol 1e-8 -pc_type mg -pc_mg_galerkin -pc_mg_adapt_interp_coarse_space eigenvector -pc_mg_adapt_interp_n 1 \
1789:             -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
1790:   test:
1791:     suffix: 2d_p1_adaptmg_gdsw
1792:     requires: triangle
1793:     nsize: 4
1794:     args: -petscpartitioner_type simple -dm_refine 3 -dm_plex_box_faces 4,4 -bc_type dirichlet -petscspace_degree 1 \
1795:           -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
1796:           -snes_max_it 1 -ksp_converged_reason \
1797:           -ksp_rtol 1e-8 -pc_type mg -pc_mg_galerkin -pc_mg_adapt_interp_coarse_space gdsw -pc_mg_levels 2 -mg_levels_pc_type asm -dm_mat_type {{aij is}}

1799:   test:
1800:     suffix: 2d_p1_adaptmg_agdsw
1801:     requires: triangle mumps
1802:     nsize: 4
1803:     args: -petscpartitioner_type simple -dm_refine 3 -dm_plex_box_faces 4,4 -bc_type dirichlet -petscspace_degree 1 \
1804:           -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
1805:           -snes_max_it 1 -ksp_converged_reason \
1806:           -ksp_rtol 1e-8 -pc_type mg -pc_mg_galerkin -pc_mg_adapt_interp_coarse_space gdsw -pc_mg_levels 2 -mg_levels_pc_type asm -dm_mat_type is -mg_levels_gdsw_tolerance 0.1 -mg_levels_gdsw_pseudo_pc_type qr

1808:   test:
1809:     suffix: p4est_2d_asm
1810:     requires: p4est
1811:     nsize: 4
1812:     args: -run_type test -run_test_check_ksp -quiet -petscspace_degree 1 -petscpartitioner_type simple -bc_type none -dm_plex_simplex 0 \
1813:           -pc_type asm -ksp_converged_reason -dm_plex_convert_type p4est -dm_forest_minimum_refinement 1 -dm_forest_initial_refinement 3 -dm_forest_maximum_refinement 5 \
1814:           -pc_asm_dm_subdomains -dm_p4est_refine_pattern hash -dm_plex_dd_overlap 1 -sub_pc_type lu

1816: TEST*/