Actual source code: ex45.c

  1: static char help[] = "Heat Equation in 2d and 3d with finite elements.\n\
  2: We solve the heat equation in a rectangular\n\
  3: domain, using a parallel unstructured mesh (DMPLEX) to discretize it.\n\
  4: Contributed by: Julian Andrej <juan@tf.uni-kiel.de>\n\n\n";

  6: #include <petscdmplex.h>
  7: #include <petscds.h>
  8: #include <petscts.h>

 10: /*
 11:   Heat equation:

 13:     du/dt - \Delta u + f = 0
 14: */

 16: typedef enum {
 17:   SOL_QUADRATIC_LINEAR,
 18:   SOL_QUADRATIC_TRIG,
 19:   SOL_TRIG_LINEAR,
 20:   SOL_TRIG_TRIG,
 21:   NUM_SOLUTION_TYPES
 22: } SolutionType;
 23: const char *solutionTypes[NUM_SOLUTION_TYPES + 1] = {"quadratic_linear", "quadratic_trig", "trig_linear", "trig_trig", "unknown"};

 25: typedef struct {
 26:   SolutionType solType; /* Type of exact solution */
 27:   /* Solver setup */
 28:   PetscBool expTS;  /* Flag for explicit timestepping */
 29:   PetscBool lumped; /* Lump the mass matrix */
 30: } AppCtx;

 32: /*
 33: Exact 2D solution:
 34:   u    = 2t + x^2 + y^2
 35:   u_t  = 2
 36:   \Delta u = 2 + 2 = 4
 37:   f    = 2
 38:   F(u) = 2 - (2 + 2) + 2 = 0

 40: Exact 3D solution:
 41:   u = 3t + x^2 + y^2 + z^2
 42:   F(u) = 3 - (2 + 2 + 2) + 3 = 0
 43: */
 44: static PetscErrorCode mms_quad_lin(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 45: {
 46:   PetscInt d;

 48:   *u = dim * time;
 49:   for (d = 0; d < dim; ++d) *u += x[d] * x[d];
 50:   return PETSC_SUCCESS;
 51: }

 53: static PetscErrorCode mms_quad_lin_t(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 54: {
 55:   *u = dim;
 56:   return PETSC_SUCCESS;
 57: }

 59: static void f0_quad_lin_exp(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[])
 60: {
 61:   f0[0] = -(PetscScalar)dim;
 62: }
 63: static void f0_quad_lin(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[])
 64: {
 65:   PetscScalar exp[1] = {0.};
 66:   f0_quad_lin_exp(dim, Nf, NfAux, uOff, uOff_x, u, u_t, u_x, aOff, aOff_x, a, a_t, a_x, t, x, numConstants, constants, exp);
 67:   f0[0] = u_t[0] - exp[0];
 68: }

 70: /*
 71: Exact 2D solution:
 72:   u = 2*cos(t) + x^2 + y^2
 73:   F(u) = -2*sint(t) - (2 + 2) + 2*sin(t) + 4 = 0

 75: Exact 3D solution:
 76:   u = 3*cos(t) + x^2 + y^2 + z^2
 77:   F(u) = -3*sin(t) - (2 + 2 + 2) + 3*sin(t) + 6 = 0
 78: */
 79: static PetscErrorCode mms_quad_trig(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 80: {
 81:   PetscInt d;

 83:   *u = dim * PetscCosReal(time);
 84:   for (d = 0; d < dim; ++d) *u += x[d] * x[d];
 85:   return PETSC_SUCCESS;
 86: }

 88: static PetscErrorCode mms_quad_trig_t(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
 89: {
 90:   *u = -dim * PetscSinReal(time);
 91:   return PETSC_SUCCESS;
 92: }

 94: static void f0_quad_trig_exp(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[])
 95: {
 96:   f0[0] = -dim * (PetscSinReal(t) + 2.0);
 97: }
 98: static void f0_quad_trig(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[])
 99: {
100:   PetscScalar exp[1] = {0.};
101:   f0_quad_trig_exp(dim, Nf, NfAux, uOff, uOff_x, u, u_t, u_x, aOff, aOff_x, a, a_t, a_x, t, x, numConstants, constants, exp);
102:   f0[0] = u_t[0] - exp[0];
103: }

105: /*
106: Exact 2D solution:
107:   u = 2\pi^2 t + cos(\pi x) + cos(\pi y)
108:   F(u) = 2\pi^2 - \pi^2 (cos(\pi x) + cos(\pi y)) + \pi^2 (cos(\pi x) + cos(\pi y)) - 2\pi^2 = 0

110: Exact 3D solution:
111:   u = 3\pi^2 t + cos(\pi x) + cos(\pi y) + cos(\pi z)
112:   F(u) = 3\pi^2 - \pi^2 (cos(\pi x) + cos(\pi y) + cos(\pi z)) + \pi^2 (cos(\pi x) + cos(\pi y) + cos(\pi z)) - 3\pi^2 = 0
113: */
114: static PetscErrorCode mms_trig_lin(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
115: {
116:   PetscInt d;

118:   *u = dim * PetscSqr(PETSC_PI) * time;
119:   for (d = 0; d < dim; ++d) *u += PetscCosReal(PETSC_PI * x[d]);
120:   return PETSC_SUCCESS;
121: }

123: static PetscErrorCode mms_trig_lin_t(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
124: {
125:   *u = dim * PetscSqr(PETSC_PI);
126:   return PETSC_SUCCESS;
127: }

129: static void f0_trig_lin(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[])
130: {
131:   PetscInt d;
132:   f0[0] = u_t[0];
133:   for (d = 0; d < dim; ++d) f0[0] += PetscSqr(PETSC_PI) * (PetscCosReal(PETSC_PI * x[d]) - 1.0);
134: }

136: /*
137: Exact 2D solution:
138:   u    = pi^2 cos(t) + cos(\pi x) + cos(\pi y)
139:   u_t  = -pi^2 sin(t)
140:   \Delta u = -\pi^2 (cos(\pi x) + cos(\pi y))
141:   f    = pi^2 sin(t) - \pi^2 (cos(\pi x) + cos(\pi y))
142:   F(u) = -\pi^2 sin(t) + \pi^2 (cos(\pi x) + cos(\pi y)) - \pi^2 (cos(\pi x) + cos(\pi y)) + \pi^2 sin(t) = 0

144: Exact 3D solution:
145:   u    = pi^2 cos(t) + cos(\pi x) + cos(\pi y) + cos(\pi z)
146:   u_t  = -pi^2 sin(t)
147:   \Delta u = -\pi^2 (cos(\pi x) + cos(\pi y) + cos(\pi z))
148:   f    = pi^2 sin(t) - \pi^2 (cos(\pi x) + cos(\pi y) + cos(\pi z))
149:   F(u) = -\pi^2 sin(t) + \pi^2 (cos(\pi x) + cos(\pi y) + cos(\pi z)) - \pi^2 (cos(\pi x) + cos(\pi y) + cos(\pi z)) + \pi^2 sin(t) = 0
150: */
151: static PetscErrorCode mms_trig_trig(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
152: {
153:   PetscInt d;

155:   *u = PetscSqr(PETSC_PI) * PetscCosReal(time);
156:   for (d = 0; d < dim; ++d) *u += PetscCosReal(PETSC_PI * x[d]);
157:   return PETSC_SUCCESS;
158: }

160: static PetscErrorCode mms_trig_trig_t(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
161: {
162:   *u = -PetscSqr(PETSC_PI) * PetscSinReal(time);
163:   return PETSC_SUCCESS;
164: }

166: static void f0_trig_trig_exp(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[])
167: {
168:   PetscInt d;
169:   f0[0] -= PetscSqr(PETSC_PI) * PetscSinReal(t);
170:   for (d = 0; d < dim; ++d) f0[0] += PetscSqr(PETSC_PI) * PetscCosReal(PETSC_PI * x[d]);
171: }
172: static void f0_trig_trig(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[])
173: {
174:   PetscScalar exp[1] = {0.};
175:   f0_trig_trig_exp(dim, Nf, NfAux, uOff, uOff_x, u, u_t, u_x, aOff, aOff_x, a, a_t, a_x, t, x, numConstants, constants, exp);
176:   f0[0] = u_t[0] - exp[0];
177: }

179: static void f1_temp_exp(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[])
180: {
181:   for (PetscInt d = 0; d < dim; ++d) f1[d] = -u_x[d];
182: }
183: static void f1_temp(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[])
184: {
185:   for (PetscInt d = 0; d < dim; ++d) f1[d] = u_x[d];
186: }

188: static void g3_temp(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[])
189: {
190:   for (PetscInt d = 0; d < dim; ++d) g3[d * dim + d] = 1.0;
191: }

193: static void g0_temp(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 g0[])
194: {
195:   g0[0] = u_tShift * 1.0;
196: }

198: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
199: {
200:   PetscInt sol;

202:   PetscFunctionBeginUser;
203:   options->solType = SOL_QUADRATIC_LINEAR;
204:   options->expTS   = PETSC_FALSE;
205:   options->lumped  = PETSC_TRUE;

207:   PetscOptionsBegin(comm, "", "Heat Equation Options", "DMPLEX");
208:   PetscCall(PetscOptionsEList("-sol_type", "Type of exact solution", "ex45.c", solutionTypes, NUM_SOLUTION_TYPES, solutionTypes[options->solType], &sol, NULL));
209:   options->solType = (SolutionType)sol;
210:   PetscCall(PetscOptionsBool("-explicit", "Use explicit timestepping", "ex45.c", options->expTS, &options->expTS, NULL));
211:   PetscCall(PetscOptionsBool("-lumped", "Lump the mass matrix", "ex45.c", options->lumped, &options->lumped, NULL));
212:   PetscOptionsEnd();
213:   PetscFunctionReturn(PETSC_SUCCESS);
214: }

216: static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm, AppCtx *ctx)
217: {
218:   PetscFunctionBeginUser;
219:   PetscCall(DMCreate(comm, dm));
220:   PetscCall(DMSetType(*dm, DMPLEX));
221:   PetscCall(DMSetFromOptions(*dm));
222:   PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
223:   PetscFunctionReturn(PETSC_SUCCESS);
224: }

226: static PetscErrorCode SetupProblem(DM dm, AppCtx *ctx)
227: {
228:   PetscDS        ds;
229:   DMLabel        label;
230:   const PetscInt id = 1;

232:   PetscFunctionBeginUser;
233:   PetscCall(DMGetLabel(dm, "marker", &label));
234:   PetscCall(DMGetDS(dm, &ds));
235:   PetscCall(PetscDSSetJacobian(ds, 0, 0, g0_temp, NULL, NULL, g3_temp));
236:   switch (ctx->solType) {
237:   case SOL_QUADRATIC_LINEAR:
238:     PetscCall(PetscDSSetResidual(ds, 0, f0_quad_lin, f1_temp));
239:     PetscCall(PetscDSSetRHSResidual(ds, 0, f0_quad_lin_exp, f1_temp_exp));
240:     PetscCall(PetscDSSetExactSolution(ds, 0, mms_quad_lin, ctx));
241:     PetscCall(PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_quad_lin_t, ctx));
242:     PetscCall(DMAddBoundary(dm, DM_BC_ESSENTIAL, "wall", label, 1, &id, 0, 0, NULL, (void (*)(void))mms_quad_lin, (void (*)(void))mms_quad_lin_t, ctx, NULL));
243:     break;
244:   case SOL_QUADRATIC_TRIG:
245:     PetscCall(PetscDSSetResidual(ds, 0, f0_quad_trig, f1_temp));
246:     PetscCall(PetscDSSetRHSResidual(ds, 0, f0_quad_trig_exp, f1_temp_exp));
247:     PetscCall(PetscDSSetExactSolution(ds, 0, mms_quad_trig, ctx));
248:     PetscCall(PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_quad_trig_t, ctx));
249:     PetscCall(DMAddBoundary(dm, DM_BC_ESSENTIAL, "wall", label, 1, &id, 0, 0, NULL, (void (*)(void))mms_quad_trig, (void (*)(void))mms_quad_trig_t, ctx, NULL));
250:     break;
251:   case SOL_TRIG_LINEAR:
252:     PetscCall(PetscDSSetResidual(ds, 0, f0_trig_lin, f1_temp));
253:     PetscCall(PetscDSSetExactSolution(ds, 0, mms_trig_lin, ctx));
254:     PetscCall(PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_trig_lin_t, ctx));
255:     PetscCall(DMAddBoundary(dm, DM_BC_ESSENTIAL, "wall", label, 1, &id, 0, 0, NULL, (void (*)(void))mms_trig_lin, (void (*)(void))mms_trig_lin_t, ctx, NULL));
256:     break;
257:   case SOL_TRIG_TRIG:
258:     PetscCall(PetscDSSetResidual(ds, 0, f0_trig_trig, f1_temp));
259:     PetscCall(PetscDSSetRHSResidual(ds, 0, f0_trig_trig_exp, f1_temp_exp));
260:     PetscCall(PetscDSSetExactSolution(ds, 0, mms_trig_trig, ctx));
261:     PetscCall(PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_trig_trig_t, ctx));
262:     break;
263:   default:
264:     SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONG, "Invalid solution type: %s (%d)", solutionTypes[PetscMin(ctx->solType, NUM_SOLUTION_TYPES)], ctx->solType);
265:   }
266:   PetscFunctionReturn(PETSC_SUCCESS);
267: }

269: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *ctx)
270: {
271:   DM             cdm = dm;
272:   PetscFE        fe;
273:   DMPolytopeType ct;
274:   PetscInt       dim, cStart;

276:   PetscFunctionBeginUser;
277:   PetscCall(DMGetDimension(dm, &dim));
278:   PetscCall(DMPlexGetHeightStratum(dm, 0, &cStart, NULL));
279:   PetscCall(DMPlexGetCellType(dm, cStart, &ct));
280:   /* Create finite element */
281:   PetscCall(PetscFECreateByCell(PETSC_COMM_SELF, dim, 1, ct, "temp_", -1, &fe));
282:   PetscCall(PetscObjectSetName((PetscObject)fe, "temperature"));
283:   /* Set discretization and boundary conditions for each mesh */
284:   PetscCall(DMSetField(dm, 0, NULL, (PetscObject)fe));
285:   PetscCall(DMCreateDS(dm));
286:   if (ctx->expTS) {
287:     PetscDS ds;

289:     PetscCall(DMGetDS(dm, &ds));
290:     PetscCall(PetscDSSetImplicit(ds, 0, PETSC_FALSE));
291:   }
292:   PetscCall(SetupProblem(dm, ctx));
293:   while (cdm) {
294:     PetscCall(DMCopyDisc(dm, cdm));
295:     PetscCall(DMGetCoarseDM(cdm, &cdm));
296:   }
297:   PetscCall(PetscFEDestroy(&fe));
298:   PetscFunctionReturn(PETSC_SUCCESS);
299: }

301: static PetscErrorCode SetInitialConditions(TS ts, Vec u)
302: {
303:   DM        dm;
304:   PetscReal t;

306:   PetscFunctionBeginUser;
307:   PetscCall(TSGetDM(ts, &dm));
308:   PetscCall(TSGetTime(ts, &t));
309:   PetscCall(DMComputeExactSolution(dm, t, u, NULL));
310:   PetscFunctionReturn(PETSC_SUCCESS);
311: }

313: int main(int argc, char **argv)
314: {
315:   DM     dm;
316:   TS     ts;
317:   Vec    u;
318:   AppCtx ctx;

320:   PetscFunctionBeginUser;
321:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
322:   PetscCall(ProcessOptions(PETSC_COMM_WORLD, &ctx));
323:   PetscCall(CreateMesh(PETSC_COMM_WORLD, &dm, &ctx));
324:   PetscCall(DMSetApplicationContext(dm, &ctx));
325:   PetscCall(SetupDiscretization(dm, &ctx));

327:   PetscCall(TSCreate(PETSC_COMM_WORLD, &ts));
328:   PetscCall(TSSetDM(ts, dm));
329:   PetscCall(DMTSSetBoundaryLocal(dm, DMPlexTSComputeBoundary, &ctx));
330:   if (ctx.expTS) {
331:     PetscCall(DMTSSetRHSFunctionLocal(dm, DMPlexTSComputeRHSFunctionFEM, &ctx));
332:     if (ctx.lumped) PetscCall(DMTSCreateRHSMassMatrixLumped(dm));
333:     else PetscCall(DMTSCreateRHSMassMatrix(dm));
334:   } else {
335:     PetscCall(DMTSSetIFunctionLocal(dm, DMPlexTSComputeIFunctionFEM, &ctx));
336:     PetscCall(DMTSSetIJacobianLocal(dm, DMPlexTSComputeIJacobianFEM, &ctx));
337:   }
338:   PetscCall(TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP));
339:   PetscCall(TSSetFromOptions(ts));
340:   PetscCall(TSSetComputeInitialCondition(ts, SetInitialConditions));

342:   PetscCall(DMCreateGlobalVector(dm, &u));
343:   PetscCall(DMTSCheckFromOptions(ts, u));
344:   PetscCall(SetInitialConditions(ts, u));
345:   PetscCall(PetscObjectSetName((PetscObject)u, "temperature"));
346:   PetscCall(TSSolve(ts, u));
347:   PetscCall(DMTSCheckFromOptions(ts, u));
348:   if (ctx.expTS) PetscCall(DMTSDestroyRHSMassMatrix(dm));

350:   PetscCall(VecDestroy(&u));
351:   PetscCall(TSDestroy(&ts));
352:   PetscCall(DMDestroy(&dm));
353:   PetscCall(PetscFinalize());
354:   return 0;
355: }

357: /*TEST

359:   test:
360:     suffix: 2d_p1
361:     requires: triangle
362:     args: -sol_type quadratic_linear -dm_refine 1 -temp_petscspace_degree 1 -dmts_check .0001 \
363:           -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
364:   test:
365:     suffix: 2d_p1_exp
366:     requires: triangle
367:     args: -sol_type quadratic_linear -dm_refine 1 -temp_petscspace_degree 1 -explicit \
368:           -ts_type euler -ts_max_steps 4 -ts_dt 1e-3 -ts_monitor_error
369:   test:
370:     # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [1.9]
371:     suffix: 2d_p1_sconv
372:     requires: triangle
373:     args: -sol_type quadratic_linear -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
374:           -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
375:   test:
376:     # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.1]
377:     suffix: 2d_p1_sconv_2
378:     requires: triangle
379:     args: -sol_type trig_trig -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
380:           -ts_type beuler -ts_max_steps 1 -ts_dt 1e-6 -snes_error_if_not_converged -pc_type lu
381:   test:
382:     # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
383:     suffix: 2d_p1_tconv
384:     requires: triangle
385:     args: -sol_type quadratic_trig -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
386:           -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
387:   test:
388:     # -dm_refine 6 -convest_num_refine 3 get L_2 convergence rate: [1.0]
389:     suffix: 2d_p1_tconv_2
390:     requires: triangle
391:     args: -sol_type trig_trig -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
392:           -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
393:   test:
394:     # The L_2 convergence rate cannot be seen since stability of the explicit integrator requires that is be more accurate than the grid
395:     suffix: 2d_p1_exp_tconv_2
396:     requires: triangle
397:     args: -sol_type trig_trig -temp_petscspace_degree 1 -explicit -ts_convergence_estimate -convest_num_refine 1 \
398:           -ts_type euler -ts_max_steps 4 -ts_dt 1e-4 -lumped 0 -mass_pc_type lu
399:   test:
400:     # The L_2 convergence rate cannot be seen since stability of the explicit integrator requires that is be more accurate than the grid
401:     suffix: 2d_p1_exp_tconv_2_lump
402:     requires: triangle
403:     args: -sol_type trig_trig -temp_petscspace_degree 1 -explicit -ts_convergence_estimate -convest_num_refine 1 \
404:           -ts_type euler -ts_max_steps 4 -ts_dt 1e-4
405:   test:
406:     suffix: 2d_p2
407:     requires: triangle
408:     args: -sol_type quadratic_linear -dm_refine 0 -temp_petscspace_degree 2 -dmts_check .0001 \
409:           -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
410:   test:
411:     # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
412:     suffix: 2d_p2_sconv
413:     requires: triangle
414:     args: -sol_type trig_linear -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
415:           -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
416:   test:
417:     # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [3.1]
418:     suffix: 2d_p2_sconv_2
419:     requires: triangle
420:     args: -sol_type trig_trig -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
421:           -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
422:   test:
423:     # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
424:     suffix: 2d_p2_tconv
425:     requires: triangle
426:     args: -sol_type quadratic_trig -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
427:           -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
428:   test:
429:     # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
430:     suffix: 2d_p2_tconv_2
431:     requires: triangle
432:     args: -sol_type trig_trig -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
433:           -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
434:   test:
435:     suffix: 2d_q1
436:     args: -sol_type quadratic_linear -dm_plex_simplex 0 -dm_refine 1 -temp_petscspace_degree 1 -dmts_check .0001 \
437:           -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
438:   test:
439:     # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [1.9]
440:     suffix: 2d_q1_sconv
441:     args: -sol_type quadratic_linear -dm_plex_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
442:           -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
443:   test:
444:     # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
445:     suffix: 2d_q1_tconv
446:     args: -sol_type quadratic_trig -dm_plex_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
447:           -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
448:   test:
449:     suffix: 2d_q2
450:     args: -sol_type quadratic_linear -dm_plex_simplex 0 -dm_refine 0 -temp_petscspace_degree 2 -dmts_check .0001 \
451:           -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
452:   test:
453:     # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
454:     suffix: 2d_q2_sconv
455:     args: -sol_type trig_linear -dm_plex_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
456:           -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
457:   test:
458:     # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
459:     suffix: 2d_q2_tconv
460:     args: -sol_type quadratic_trig -dm_plex_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
461:           -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu

463:   test:
464:     suffix: 3d_p1
465:     requires: ctetgen
466:     args: -sol_type quadratic_linear -dm_refine 1 -temp_petscspace_degree 1 -dmts_check .0001 \
467:           -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
468:   test:
469:     # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [1.9]
470:     suffix: 3d_p1_sconv
471:     requires: ctetgen
472:     args: -sol_type quadratic_linear -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
473:           -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
474:   test:
475:     # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
476:     suffix: 3d_p1_tconv
477:     requires: ctetgen
478:     args: -sol_type quadratic_trig -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
479:           -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
480:   test:
481:     suffix: 3d_p2
482:     requires: ctetgen
483:     args: -sol_type quadratic_linear -dm_refine 0 -temp_petscspace_degree 2 -dmts_check .0001 \
484:           -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
485:   test:
486:     # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
487:     suffix: 3d_p2_sconv
488:     requires: ctetgen
489:     args: -sol_type trig_linear -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
490:           -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
491:   test:
492:     # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
493:     suffix: 3d_p2_tconv
494:     requires: ctetgen
495:     args: -sol_type quadratic_trig -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
496:           -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
497:   test:
498:     suffix: 3d_q1
499:     args: -sol_type quadratic_linear -dm_plex_simplex 0 -dm_refine 1 -temp_petscspace_degree 1 -dmts_check .0001 \
500:           -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
501:   test:
502:     # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [1.9]
503:     suffix: 3d_q1_sconv
504:     args: -sol_type quadratic_linear -dm_plex_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
505:           -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
506:   test:
507:     # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
508:     suffix: 3d_q1_tconv
509:     args: -sol_type quadratic_trig -dm_plex_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
510:           -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
511:   test:
512:     suffix: 3d_q2
513:     args: -sol_type quadratic_linear -dm_plex_simplex 0 -dm_refine 0 -temp_petscspace_degree 2 -dmts_check .0001 \
514:           -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
515:   test:
516:     # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
517:     suffix: 3d_q2_sconv
518:     args: -sol_type trig_linear -dm_plex_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
519:           -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
520:   test:
521:     # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
522:     suffix: 3d_q2_tconv
523:     args: -sol_type quadratic_trig -dm_plex_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
524:           -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu

526:   test:
527:     # For a nice picture, -bd_dm_refine 2 -dm_refine 1 -dm_view hdf5:${PETSC_DIR}/sol.h5 -ts_monitor_solution hdf5:${PETSC_DIR}/sol.h5::append
528:     suffix: egads_sphere
529:     requires: egads ctetgen
530:     args: -sol_type quadratic_linear \
531:           -dm_plex_boundary_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/sphere_example.egadslite -dm_plex_boundary_label marker \
532:           -temp_petscspace_degree 2 -dmts_check .0001 \
533:           -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu

535: TEST*/