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 {SOL_QUADRATIC_LINEAR, SOL_QUADRATIC_TRIG, SOL_TRIG_LINEAR, NUM_SOLUTION_TYPES} SolutionType;
17: const char *solutionTypes[NUM_SOLUTION_TYPES+1] = {"quadratic_linear", "quadratic_trig", "trig_linear", "unknown"};
19: typedef struct {
20: SolutionType solType; /* Type of exact solution */
21: } AppCtx;
23: /*
24: Exact 2D solution:
25: u = 2t + x^2 + y^2
26: F(u) = 2 - (2 + 2) + 2 = 0
28: Exact 3D solution:
29: u = 3t + x^2 + y^2 + z^2
30: F(u) = 3 - (2 + 2 + 2) + 3 = 0
31: */
32: static PetscErrorCode mms_quad_lin(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
33: {
34: PetscInt d;
36: *u = dim*time;
37: for (d = 0; d < dim; ++d) *u += x[d]*x[d];
38: return 0;
39: }
41: static PetscErrorCode mms_quad_lin_t(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
42: {
43: *u = dim;
44: return 0;
45: }
47: static void f0_quad_lin(PetscInt dim, PetscInt Nf, PetscInt NfAux,
48: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
49: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
50: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
51: {
52: f0[0] = u_t[0] + (PetscScalar) dim;
53: }
55: /*
56: Exact 2D solution:
57: u = 2*cos(t) + x^2 + y^2
58: F(u) = -2*sint(t) - (2 + 2) + 2*sin(t) + 4 = 0
60: Exact 3D solution:
61: u = 3*cos(t) + x^2 + y^2 + z^2
62: F(u) = -3*sin(t) - (2 + 2 + 2) + 3*sin(t) + 6 = 0
63: */
64: static PetscErrorCode mms_quad_trig(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
65: {
66: PetscInt d;
68: *u = dim*PetscCosReal(time);
69: for (d = 0; d < dim; ++d) *u += x[d]*x[d];
70: return 0;
71: }
73: static PetscErrorCode mms_quad_trig_t(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
74: {
75: *u = -dim*PetscSinReal(time);
76: return 0;
77: }
79: static void f0_quad_trig(PetscInt dim, PetscInt Nf, PetscInt NfAux,
80: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
81: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
82: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
83: {
84: f0[0] = u_t[0] + dim*(PetscSinReal(t) + 2.0);
85: }
87: /*
88: Exact 2D solution:
89: u = 2\pi^2 t + cos(\pi x) + cos(\pi y)
90: 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
92: Exact 3D solution:
93: u = 3\pi^2 t + cos(\pi x) + cos(\pi y) + cos(\pi z)
94: 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
95: */
96: static PetscErrorCode mms_trig_lin(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
97: {
98: PetscInt d;
100: *u = dim*PetscSqr(PETSC_PI)*time;
101: for (d = 0; d < dim; ++d) *u += PetscCosReal(PETSC_PI*x[d]);
102: return 0;
103: }
105: static PetscErrorCode mms_trig_lin_t(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
106: {
107: *u = dim*PetscSqr(PETSC_PI);
108: return 0;
109: }
111: static void f0_trig_lin(PetscInt dim, PetscInt Nf, PetscInt NfAux,
112: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
113: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
114: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
115: {
116: PetscInt d;
117: f0[0] = u_t[0];
118: for (d = 0; d < dim; ++d) f0[0] += PetscSqr(PETSC_PI)*(PetscCosReal(PETSC_PI*x[d]) - 1.0);
119: }
121: static void f1_temp(PetscInt dim, PetscInt Nf, PetscInt NfAux,
122: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
123: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
124: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
125: {
126: PetscInt d;
127: for (d = 0; d < dim; ++d) f1[d] = u_x[d];
128: }
130: static void g3_temp(PetscInt dim, PetscInt Nf, PetscInt NfAux,
131: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
132: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
133: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
134: {
135: PetscInt d;
136: for (d = 0; d < dim; ++d) g3[d*dim+d] = 1.0;
137: }
139: static void g0_temp(PetscInt dim, PetscInt Nf, PetscInt NfAux,
140: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
141: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
142: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[])
143: {
144: g0[0] = u_tShift*1.0;
145: }
147: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
148: {
149: PetscInt sol;
153: options->solType = SOL_QUADRATIC_LINEAR;
155: PetscOptionsBegin(comm, "", "Heat Equation Options", "DMPLEX");
156: PetscOptionsEList("-sol_type", "Type of exact solution", "ex45.c", solutionTypes, NUM_SOLUTION_TYPES, solutionTypes[options->solType], &sol, NULL);
157: options->solType = (SolutionType) sol;
158: PetscOptionsEnd();
159: return(0);
160: }
162: static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm, AppCtx *ctx)
163: {
167: DMCreate(comm, dm);
168: DMSetType(*dm, DMPLEX);
169: DMSetFromOptions(*dm);
170: DMViewFromOptions(*dm, NULL, "-dm_view");
171: return(0);
172: }
174: static PetscErrorCode SetupProblem(DM dm, AppCtx *ctx)
175: {
176: PetscDS ds;
177: DMLabel label;
178: const PetscInt id = 1;
182: DMGetLabel(dm, "marker", &label);
183: DMGetDS(dm, &ds);
184: PetscDSSetJacobian(ds, 0, 0, g0_temp, NULL, NULL, g3_temp);
185: switch (ctx->solType) {
186: case SOL_QUADRATIC_LINEAR:
187: PetscDSSetResidual(ds, 0, f0_quad_lin, f1_temp);
188: PetscDSSetExactSolution(ds, 0, mms_quad_lin, ctx);
189: PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_quad_lin_t, ctx);
190: 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);
191: break;
192: case SOL_QUADRATIC_TRIG:
193: PetscDSSetResidual(ds, 0, f0_quad_trig, f1_temp);
194: PetscDSSetExactSolution(ds, 0, mms_quad_trig, ctx);
195: PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_quad_trig_t, ctx);
196: 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);
197: break;
198: case SOL_TRIG_LINEAR:
199: PetscDSSetResidual(ds, 0, f0_trig_lin, f1_temp);
200: PetscDSSetExactSolution(ds, 0, mms_trig_lin, ctx);
201: PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_trig_lin_t, ctx);
202: 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);
203: break;
204: default: SETERRQ2(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "Invalid solution type: %s (%D)", solutionTypes[PetscMin(ctx->solType, NUM_SOLUTION_TYPES)], ctx->solType);
205: }
206: return(0);
207: }
209: static PetscErrorCode SetupDiscretization(DM dm, AppCtx* ctx)
210: {
211: DM cdm = dm;
212: PetscFE fe;
213: DMPolytopeType ct;
214: PetscBool simplex;
215: PetscInt dim, cStart;
219: DMGetDimension(dm, &dim);
220: DMPlexGetHeightStratum(dm, 0, &cStart, NULL);
221: DMPlexGetCellType(dm, cStart, &ct);
222: simplex = DMPolytopeTypeGetNumVertices(ct) == DMPolytopeTypeGetDim(ct)+1 ? PETSC_TRUE : PETSC_FALSE;
223: /* Create finite element */
224: PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, "temp_", -1, &fe);
225: PetscObjectSetName((PetscObject) fe, "temperature");
226: /* Set discretization and boundary conditions for each mesh */
227: DMSetField(dm, 0, NULL, (PetscObject) fe);
228: DMCreateDS(dm);
229: SetupProblem(dm, ctx);
230: while (cdm) {
231: DMCopyDisc(dm, cdm);
232: DMGetCoarseDM(cdm, &cdm);
233: }
234: PetscFEDestroy(&fe);
235: return(0);
236: }
238: static PetscErrorCode SetInitialConditions(TS ts, Vec u)
239: {
240: DM dm;
241: PetscReal t;
245: TSGetDM(ts, &dm);
246: TSGetTime(ts, &t);
247: DMComputeExactSolution(dm, t, u, NULL);
248: return(0);
249: }
251: int main(int argc, char **argv)
252: {
253: DM dm;
254: TS ts;
255: Vec u;
256: AppCtx ctx;
259: PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr;
260: ProcessOptions(PETSC_COMM_WORLD, &ctx);
261: CreateMesh(PETSC_COMM_WORLD, &dm, &ctx);
262: DMSetApplicationContext(dm, &ctx);
263: SetupDiscretization(dm, &ctx);
265: TSCreate(PETSC_COMM_WORLD, &ts);
266: TSSetDM(ts, dm);
267: DMTSSetBoundaryLocal(dm, DMPlexTSComputeBoundary, &ctx);
268: DMTSSetIFunctionLocal(dm, DMPlexTSComputeIFunctionFEM, &ctx);
269: DMTSSetIJacobianLocal(dm, DMPlexTSComputeIJacobianFEM, &ctx);
270: TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
271: TSSetFromOptions(ts);
272: TSSetComputeInitialCondition(ts, SetInitialConditions);
274: DMCreateGlobalVector(dm, &u);
275: DMTSCheckFromOptions(ts, u);
276: SetInitialConditions(ts, u);
277: PetscObjectSetName((PetscObject) u, "temperature");
278: TSSolve(ts, u);
279: DMTSCheckFromOptions(ts, u);
281: VecDestroy(&u);
282: TSDestroy(&ts);
283: DMDestroy(&dm);
284: PetscFinalize();
285: return ierr;
286: }
288: /*TEST
290: test:
291: suffix: 2d_p1
292: requires: triangle
293: args: -sol_type quadratic_linear -dm_refine 1 -temp_petscspace_degree 1 -dmts_check .0001 \
294: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
295: test:
296: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [1.9]
297: suffix: 2d_p1_sconv
298: requires: triangle
299: args: -sol_type quadratic_linear -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
300: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
301: test:
302: # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
303: suffix: 2d_p1_tconv
304: requires: triangle
305: args: -sol_type quadratic_trig -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
306: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
307: test:
308: suffix: 2d_p2
309: requires: triangle
310: args: -sol_type quadratic_linear -dm_refine 0 -temp_petscspace_degree 2 -dmts_check .0001 \
311: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
312: test:
313: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
314: suffix: 2d_p2_sconv
315: requires: triangle
316: args: -sol_type trig_linear -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
317: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
318: test:
319: # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
320: suffix: 2d_p2_tconv
321: requires: triangle
322: args: -sol_type quadratic_trig -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
323: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
324: test:
325: suffix: 2d_q1
326: args: -sol_type quadratic_linear -dm_plex_simplex 0 -dm_refine 1 -temp_petscspace_degree 1 -dmts_check .0001 \
327: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
328: test:
329: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [1.9]
330: suffix: 2d_q1_sconv
331: args: -sol_type quadratic_linear -dm_plex_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
332: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
333: test:
334: # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
335: suffix: 2d_q1_tconv
336: args: -sol_type quadratic_trig -dm_plex_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
337: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
338: test:
339: suffix: 2d_q2
340: args: -sol_type quadratic_linear -dm_plex_simplex 0 -dm_refine 0 -temp_petscspace_degree 2 -dmts_check .0001 \
341: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
342: test:
343: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
344: suffix: 2d_q2_sconv
345: args: -sol_type trig_linear -dm_plex_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
346: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
347: test:
348: # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
349: suffix: 2d_q2_tconv
350: args: -sol_type quadratic_trig -dm_plex_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
351: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
353: test:
354: suffix: 3d_p1
355: requires: ctetgen
356: args: -sol_type quadratic_linear -dm_refine 1 -temp_petscspace_degree 1 -dmts_check .0001 \
357: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
358: test:
359: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [1.9]
360: suffix: 3d_p1_sconv
361: requires: ctetgen
362: args: -sol_type quadratic_linear -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
363: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
364: test:
365: # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
366: suffix: 3d_p1_tconv
367: requires: ctetgen
368: args: -sol_type quadratic_trig -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
369: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
370: test:
371: suffix: 3d_p2
372: requires: ctetgen
373: args: -sol_type quadratic_linear -dm_refine 0 -temp_petscspace_degree 2 -dmts_check .0001 \
374: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
375: test:
376: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
377: suffix: 3d_p2_sconv
378: requires: ctetgen
379: args: -sol_type trig_linear -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
380: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
381: test:
382: # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
383: suffix: 3d_p2_tconv
384: requires: ctetgen
385: args: -sol_type quadratic_trig -temp_petscspace_degree 2 -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: suffix: 3d_q1
389: args: -sol_type quadratic_linear -dm_plex_simplex 0 -dm_refine 1 -temp_petscspace_degree 1 -dmts_check .0001 \
390: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
391: test:
392: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [1.9]
393: suffix: 3d_q1_sconv
394: args: -sol_type quadratic_linear -dm_plex_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
395: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
396: test:
397: # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
398: suffix: 3d_q1_tconv
399: args: -sol_type quadratic_trig -dm_plex_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
400: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
401: test:
402: suffix: 3d_q2
403: args: -sol_type quadratic_linear -dm_plex_simplex 0 -dm_refine 0 -temp_petscspace_degree 2 -dmts_check .0001 \
404: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
405: test:
406: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
407: suffix: 3d_q2_sconv
408: args: -sol_type trig_linear -dm_plex_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
409: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
410: test:
411: # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
412: suffix: 3d_q2_tconv
413: args: -sol_type quadratic_trig -dm_plex_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
414: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
416: test:
417: # 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
418: suffix: egads_sphere
419: requires: egads ctetgen
420: args: -sol_type quadratic_linear \
421: -dm_plex_boundary_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/unit_sphere.egadslite -dm_plex_boundary_label marker -bd_dm_plex_scale 40 \
422: -temp_petscspace_degree 2 -dmts_check .0001 \
423: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
425: TEST*/