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: char filename[PETSC_MAX_PATH_LEN]; /* Mesh filename */
21: char bdfilename[PETSC_MAX_PATH_LEN]; /* Mesh boundary filename */
22: PetscReal scale; /* Scale factor for mesh */
23: SolutionType solType; /* Type of exact solution */
24: } AppCtx;
26: /*
27: Exact 2D solution:
28: u = 2t + x^2 + y^2
29: F(u) = 2 - (2 + 2) + 2 = 0
31: Exact 3D solution:
32: u = 3t + x^2 + y^2 + z^2
33: F(u) = 3 - (2 + 2 + 2) + 3 = 0
34: */
35: static PetscErrorCode mms_quad_lin(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
36: {
37: PetscInt d;
39: *u = dim*time;
40: for (d = 0; d < dim; ++d) *u += x[d]*x[d];
41: return 0;
42: }
44: static PetscErrorCode mms_quad_lin_t(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
45: {
46: *u = dim;
47: return 0;
48: }
50: static void f0_quad_lin(PetscInt dim, PetscInt Nf, PetscInt NfAux,
51: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
52: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
53: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
54: {
55: f0[0] = u_t[0] + (PetscScalar) dim;
56: }
58: /*
59: Exact 2D solution:
60: u = 2*cos(t) + x^2 + y^2
61: F(u) = -2*sint(t) - (2 + 2) + 2*sin(t) + 4 = 0
63: Exact 3D solution:
64: u = 3*cos(t) + x^2 + y^2 + z^2
65: F(u) = -3*sin(t) - (2 + 2 + 2) + 3*sin(t) + 6 = 0
66: */
67: static PetscErrorCode mms_quad_trig(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
68: {
69: PetscInt d;
71: *u = dim*PetscCosReal(time);
72: for (d = 0; d < dim; ++d) *u += x[d]*x[d];
73: return 0;
74: }
76: static PetscErrorCode mms_quad_trig_t(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
77: {
78: *u = -dim*PetscSinReal(time);
79: return 0;
80: }
82: static void f0_quad_trig(PetscInt dim, PetscInt Nf, PetscInt NfAux,
83: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
84: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
85: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
86: {
87: f0[0] = u_t[0] + dim*(PetscSinReal(t) + 2.0);
88: }
90: /*
91: Exact 2D solution:
92: u = 2\pi^2 t + cos(\pi x) + cos(\pi y)
93: 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
95: Exact 3D solution:
96: u = 3\pi^2 t + cos(\pi x) + cos(\pi y) + cos(\pi z)
97: 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
98: */
99: static PetscErrorCode mms_trig_lin(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
100: {
101: PetscInt d;
103: *u = dim*PetscSqr(PETSC_PI)*time;
104: for (d = 0; d < dim; ++d) *u += PetscCosReal(PETSC_PI*x[d]);
105: return 0;
106: }
108: static PetscErrorCode mms_trig_lin_t(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
109: {
110: *u = dim*PetscSqr(PETSC_PI);
111: return 0;
112: }
114: static void f0_trig_lin(PetscInt dim, PetscInt Nf, PetscInt NfAux,
115: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
116: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
117: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
118: {
119: PetscInt d;
120: f0[0] = u_t[0];
121: for (d = 0; d < dim; ++d) f0[0] += PetscSqr(PETSC_PI)*(PetscCosReal(PETSC_PI*x[d]) - 1.0);
122: }
124: static void f1_temp(PetscInt dim, PetscInt Nf, PetscInt NfAux,
125: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
126: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
127: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
128: {
129: PetscInt d;
130: for (d = 0; d < dim; ++d) f1[d] = u_x[d];
131: }
133: static void g3_temp(PetscInt dim, PetscInt Nf, PetscInt NfAux,
134: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
135: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
136: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
137: {
138: PetscInt d;
139: for (d = 0; d < dim; ++d) g3[d*dim+d] = 1.0;
140: }
142: static void g0_temp(PetscInt dim, PetscInt Nf, PetscInt NfAux,
143: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
144: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
145: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[])
146: {
147: g0[0] = u_tShift*1.0;
148: }
150: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
151: {
152: PetscInt sol;
156: options->filename[0] = '\0';
157: options->bdfilename[0] = '\0';
158: options->scale = 0.0;
159: options->solType = SOL_QUADRATIC_LINEAR;
161: PetscOptionsBegin(comm, "", "Heat Equation Options", "DMPLEX");
162: PetscOptionsString("-filename", "The mesh file", "ex45.c", options->filename, options->filename, PETSC_MAX_PATH_LEN, NULL);
163: PetscOptionsString("-bd_filename", "The mesh boundary file", "ex45.c", options->bdfilename, options->bdfilename, PETSC_MAX_PATH_LEN, NULL);
164: PetscOptionsReal("-scale", "Scale factor for the mesh", "ex45.c", options->scale, &options->scale, NULL);
165: sol = options->solType;
166: PetscOptionsEList("-sol_type", "Type of exact solution", "ex45.c", solutionTypes, NUM_SOLUTION_TYPES, solutionTypes[options->solType], &sol, NULL);
167: options->solType = (SolutionType) sol;
168: PetscOptionsEnd();
169: return(0);
170: }
172: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
173: {
174: DM plex;
175: DMLabel label;
176: PetscBool hasLabel;
180: DMHasLabel(dm, name, &hasLabel);
181: if (hasLabel) return(0);
182: DMCreateLabel(dm, name);
183: DMGetLabel(dm, name, &label);
184: DMConvert(dm, DMPLEX, &plex);
185: DMPlexMarkBoundaryFaces(plex, 1, label);
186: DMDestroy(&plex);
187: return(0);
188: }
190: static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm, AppCtx *ctx)
191: {
192: size_t len, lenbd;
196: PetscStrlen(ctx->filename, &len);
197: PetscStrlen(ctx->bdfilename, &lenbd);
198: if (lenbd) {
199: DM bdm;
201: DMPlexCreateFromFile(comm, ctx->bdfilename, PETSC_TRUE, &bdm);
202: PetscObjectSetOptionsPrefix((PetscObject) bdm, "bd_");
203: DMSetFromOptions(bdm);
204: if (ctx->scale != 0.0) {
205: Vec coordinates, coordinatesLocal;
207: DMGetCoordinates(bdm, &coordinates);
208: DMGetCoordinatesLocal(bdm, &coordinatesLocal);
209: VecScale(coordinates, ctx->scale);
210: VecScale(coordinatesLocal, ctx->scale);
211: }
212: DMViewFromOptions(bdm, NULL, "-dm_view");
213: DMPlexGenerate(bdm, NULL, PETSC_TRUE, dm);
214: DMDestroy(&bdm);
215: } else if (len) {
216: DMPlexCreateFromFile(comm, ctx->filename, PETSC_TRUE, dm);
217: } else {
218: DMPlexCreateBoxMesh(comm, 2, PETSC_TRUE, NULL, NULL, NULL, NULL, PETSC_TRUE, dm);
219: }
220: DMSetFromOptions(*dm);
221: PetscObjectSetName((PetscObject) *dm, "Mesh");
222: CreateBCLabel(*dm, "marker");
223: DMViewFromOptions(*dm, NULL, "-dm_view");
224: return(0);
225: }
227: static PetscErrorCode SetupProblem(DM dm, AppCtx *ctx)
228: {
229: PetscDS ds;
230: const PetscInt id = 1;
234: DMGetDS(dm, &ds);
235: PetscDSSetJacobian(ds, 0, 0, g0_temp, NULL, NULL, g3_temp);
236: switch (ctx->solType) {
237: case SOL_QUADRATIC_LINEAR:
238: PetscDSSetResidual(ds, 0, f0_quad_lin, f1_temp);
239: PetscDSSetExactSolution(ds, 0, mms_quad_lin, ctx);
240: PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_quad_lin_t, ctx);
241: DMAddBoundary(dm, DM_BC_ESSENTIAL, "wall", "marker", 0, 0, NULL, (void (*)(void)) mms_quad_lin, (void (*)(void)) mms_quad_lin_t, 1, &id, ctx);
242: break;
243: case SOL_QUADRATIC_TRIG:
244: PetscDSSetResidual(ds, 0, f0_quad_trig, f1_temp);
245: PetscDSSetExactSolution(ds, 0, mms_quad_trig, ctx);
246: PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_quad_trig_t, ctx);
247: DMAddBoundary(dm, DM_BC_ESSENTIAL, "wall", "marker", 0, 0, NULL, (void (*)(void)) mms_quad_trig, (void (*)(void)) mms_quad_trig_t, 1, &id, ctx);
248: break;
249: case SOL_TRIG_LINEAR:
250: PetscDSSetResidual(ds, 0, f0_trig_lin, f1_temp);
251: PetscDSSetExactSolution(ds, 0, mms_trig_lin, ctx);
252: PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_trig_lin_t, ctx);
253: DMAddBoundary(dm, DM_BC_ESSENTIAL, "wall", "marker", 0, 0, NULL, (void (*)(void)) mms_trig_lin, (void (*)(void)) mms_trig_lin_t, 1, &id, ctx);
254: break;
255: default: SETERRQ2(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "Invalid solution type: %s (%D)", solutionTypes[PetscMin(ctx->solType, NUM_SOLUTION_TYPES)], ctx->solType);
256: }
257: return(0);
258: }
260: static PetscErrorCode SetupDiscretization(DM dm, AppCtx* ctx)
261: {
262: DM cdm = dm;
263: PetscFE fe;
264: DMPolytopeType ct;
265: PetscBool simplex;
266: PetscInt dim, cStart;
270: DMGetDimension(dm, &dim);
271: DMPlexGetHeightStratum(dm, 0, &cStart, NULL);
272: DMPlexGetCellType(dm, cStart, &ct);
273: simplex = DMPolytopeTypeGetNumVertices(ct) == DMPolytopeTypeGetDim(ct)+1 ? PETSC_TRUE : PETSC_FALSE;
274: /* Create finite element */
275: PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, "temp_", -1, &fe);
276: PetscObjectSetName((PetscObject) fe, "temperature");
277: /* Set discretization and boundary conditions for each mesh */
278: DMSetField(dm, 0, NULL, (PetscObject) fe);
279: DMCreateDS(dm);
280: SetupProblem(dm, ctx);
281: while (cdm) {
282: CreateBCLabel(cdm, "marker");
283: DMCopyDisc(dm, cdm);
284: DMGetCoarseDM(cdm, &cdm);
285: }
286: PetscFEDestroy(&fe);
287: return(0);
288: }
290: static PetscErrorCode SetInitialConditions(TS ts, Vec u)
291: {
292: DM dm;
293: PetscReal t;
297: TSGetDM(ts, &dm);
298: TSGetTime(ts, &t);
299: DMComputeExactSolution(dm, t, u, NULL);
300: return(0);
301: }
303: int main(int argc, char **argv)
304: {
305: DM dm;
306: TS ts;
307: Vec u;
308: AppCtx ctx;
311: PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr;
312: ProcessOptions(PETSC_COMM_WORLD, &ctx);
313: CreateMesh(PETSC_COMM_WORLD, &dm, &ctx);
314: DMSetApplicationContext(dm, &ctx);
315: SetupDiscretization(dm, &ctx);
317: TSCreate(PETSC_COMM_WORLD, &ts);
318: TSSetDM(ts, dm);
319: DMTSSetBoundaryLocal(dm, DMPlexTSComputeBoundary, &ctx);
320: DMTSSetIFunctionLocal(dm, DMPlexTSComputeIFunctionFEM, &ctx);
321: DMTSSetIJacobianLocal(dm, DMPlexTSComputeIJacobianFEM, &ctx);
322: TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);
323: TSSetFromOptions(ts);
324: TSSetComputeInitialCondition(ts, SetInitialConditions);
326: DMCreateGlobalVector(dm, &u);
327: DMTSCheckFromOptions(ts, u);
328: SetInitialConditions(ts, u);
329: PetscObjectSetName((PetscObject) u, "temperature");
330: TSSolve(ts, u);
331: DMTSCheckFromOptions(ts, u);
333: VecDestroy(&u);
334: TSDestroy(&ts);
335: DMDestroy(&dm);
336: PetscFinalize();
337: return ierr;
338: }
340: /*TEST
342: test:
343: suffix: 2d_p1
344: requires: triangle
345: args: -sol_type quadratic_linear -dm_refine 1 -temp_petscspace_degree 1 -dmts_check .0001 \
346: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
347: test:
348: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [1.9]
349: suffix: 2d_p1_sconv
350: requires: triangle
351: args: -sol_type quadratic_linear -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
352: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
353: test:
354: # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
355: suffix: 2d_p1_tconv
356: requires: triangle
357: args: -sol_type quadratic_trig -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
358: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
359: test:
360: suffix: 2d_p2
361: requires: triangle
362: args: -sol_type quadratic_linear -dm_refine 0 -temp_petscspace_degree 2 -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: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
366: suffix: 2d_p2_sconv
367: requires: triangle
368: args: -sol_type trig_linear -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
369: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
370: test:
371: # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
372: suffix: 2d_p2_tconv
373: requires: triangle
374: args: -sol_type quadratic_trig -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
375: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
376: test:
377: suffix: 2d_q1
378: args: -sol_type quadratic_linear -dm_plex_box_simplex 0 -dm_refine 1 -temp_petscspace_degree 1 -dmts_check .0001 \
379: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
380: test:
381: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [1.9]
382: suffix: 2d_q1_sconv
383: args: -sol_type quadratic_linear -dm_plex_box_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
384: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
385: test:
386: # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
387: suffix: 2d_q1_tconv
388: args: -sol_type quadratic_trig -dm_plex_box_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
389: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
390: test:
391: suffix: 2d_q2
392: args: -sol_type quadratic_linear -dm_plex_box_simplex 0 -dm_refine 0 -temp_petscspace_degree 2 -dmts_check .0001 \
393: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
394: test:
395: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
396: suffix: 2d_q2_sconv
397: args: -sol_type trig_linear -dm_plex_box_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
398: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
399: test:
400: # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
401: suffix: 2d_q2_tconv
402: args: -sol_type quadratic_trig -dm_plex_box_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
403: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
405: test:
406: suffix: 3d_p1
407: requires: ctetgen
408: args: -sol_type quadratic_linear -dm_refine 1 -temp_petscspace_degree 1 -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: [1.9]
412: suffix: 3d_p1_sconv
413: requires: ctetgen
414: args: -sol_type quadratic_linear -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
415: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
416: test:
417: # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
418: suffix: 3d_p1_tconv
419: requires: ctetgen
420: args: -sol_type quadratic_trig -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
421: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
422: test:
423: suffix: 3d_p2
424: requires: ctetgen
425: args: -sol_type quadratic_linear -dm_refine 0 -temp_petscspace_degree 2 -dmts_check .0001 \
426: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
427: test:
428: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
429: suffix: 3d_p2_sconv
430: requires: ctetgen
431: args: -sol_type trig_linear -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
432: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
433: test:
434: # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
435: suffix: 3d_p2_tconv
436: requires: ctetgen
437: args: -sol_type quadratic_trig -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
438: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
439: test:
440: suffix: 3d_q1
441: args: -sol_type quadratic_linear -dm_plex_box_simplex 0 -dm_refine 1 -temp_petscspace_degree 1 -dmts_check .0001 \
442: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
443: test:
444: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [1.9]
445: suffix: 3d_q1_sconv
446: args: -sol_type quadratic_linear -dm_plex_box_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
447: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00001 -snes_error_if_not_converged -pc_type lu
448: test:
449: # -dm_refine 4 -convest_num_refine 3 get L_2 convergence rate: [1.2]
450: suffix: 3d_q1_tconv
451: args: -sol_type quadratic_trig -dm_plex_box_simplex 0 -temp_petscspace_degree 1 -ts_convergence_estimate -convest_num_refine 1 \
452: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
453: test:
454: suffix: 3d_q2
455: args: -sol_type quadratic_linear -dm_plex_box_simplex 0 -dm_refine 0 -temp_petscspace_degree 2 -dmts_check .0001 \
456: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
457: test:
458: # -dm_refine 2 -convest_num_refine 3 get L_2 convergence rate: [2.9]
459: suffix: 3d_q2_sconv
460: args: -sol_type trig_linear -dm_plex_box_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -ts_convergence_temporal 0 -convest_num_refine 1 \
461: -ts_type beuler -ts_max_steps 1 -ts_dt 0.00000001 -snes_error_if_not_converged -pc_type lu
462: test:
463: # -dm_refine 3 -convest_num_refine 3 get L_2 convergence rate: [1.0]
464: suffix: 3d_q2_tconv
465: args: -sol_type quadratic_trig -dm_plex_box_simplex 0 -temp_petscspace_degree 2 -ts_convergence_estimate -convest_num_refine 1 \
466: -ts_type beuler -ts_max_steps 4 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
468: test:
469: # 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
470: suffix: egads_sphere
471: requires: egads ctetgen
472: args: -sol_type quadratic_linear -bd_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/unit_sphere.egadslite -scale 40 \
473: -temp_petscspace_degree 2 -dmts_check .0001 \
474: -ts_type beuler -ts_max_steps 5 -ts_dt 0.1 -snes_error_if_not_converged -pc_type lu
476: TEST*/