xref: /petsc/src/ts/tutorials/ex45.c (revision 4e278199b78715991f5c71ebbd945c1489263e6c)
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";
5 
6 #include <petscdmplex.h>
7 #include <petscds.h>
8 #include <petscts.h>
9 
10 /*
11   Heat equation:
12 
13     du/dt - \Delta u + f = 0
14 */
15 
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"};
18 
19 typedef struct {
20   SolutionType solType; /* Type of exact solution */
21 } AppCtx;
22 
23 /*
24 Exact 2D solution:
25   u = 2t + x^2 + y^2
26   F(u) = 2 - (2 + 2) + 2 = 0
27 
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;
35 
36   *u = dim*time;
37   for (d = 0; d < dim; ++d) *u += x[d]*x[d];
38   return 0;
39 }
40 
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 }
46 
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 }
54 
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
59 
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;
67 
68   *u = dim*PetscCosReal(time);
69   for (d = 0; d < dim; ++d) *u += x[d]*x[d];
70   return 0;
71 }
72 
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 }
78 
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 }
86 
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
91 
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;
99 
100   *u = dim*PetscSqr(PETSC_PI)*time;
101   for (d = 0; d < dim; ++d) *u += PetscCosReal(PETSC_PI*x[d]);
102   return 0;
103 }
104 
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 }
110 
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 }
120 
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 }
129 
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 }
138 
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 }
146 
147 static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
148 {
149   PetscInt       sol;
150   PetscErrorCode ierr;
151 
152   PetscFunctionBeginUser;
153   options->solType = SOL_QUADRATIC_LINEAR;
154 
155   ierr = PetscOptionsBegin(comm, "", "Heat Equation Options", "DMPLEX");CHKERRQ(ierr);
156   ierr = PetscOptionsEList("-sol_type", "Type of exact solution", "ex45.c", solutionTypes, NUM_SOLUTION_TYPES, solutionTypes[options->solType], &sol, NULL);CHKERRQ(ierr);
157   options->solType = (SolutionType) sol;
158   ierr = PetscOptionsEnd();CHKERRQ(ierr);
159   PetscFunctionReturn(0);
160 }
161 
162 static PetscErrorCode CreateMesh(MPI_Comm comm, DM *dm, AppCtx *ctx)
163 {
164   PetscErrorCode ierr;
165 
166   PetscFunctionBeginUser;
167   ierr = DMCreate(comm, dm);CHKERRQ(ierr);
168   ierr = DMSetType(*dm, DMPLEX);CHKERRQ(ierr);
169   ierr = DMSetFromOptions(*dm);CHKERRQ(ierr);
170   ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr);
171   PetscFunctionReturn(0);
172 }
173 
174 static PetscErrorCode SetupProblem(DM dm, AppCtx *ctx)
175 {
176   PetscDS        ds;
177   DMLabel        label;
178   const PetscInt id = 1;
179   PetscErrorCode ierr;
180 
181   PetscFunctionBeginUser;
182   ierr = DMGetLabel(dm, "marker", &label);CHKERRQ(ierr);
183   ierr = DMGetDS(dm, &ds);CHKERRQ(ierr);
184   ierr = PetscDSSetJacobian(ds, 0, 0, g0_temp, NULL, NULL, g3_temp);CHKERRQ(ierr);
185   switch (ctx->solType) {
186     case SOL_QUADRATIC_LINEAR:
187       ierr = PetscDSSetResidual(ds, 0, f0_quad_lin,  f1_temp);CHKERRQ(ierr);
188       ierr = PetscDSSetExactSolution(ds, 0, mms_quad_lin, ctx);CHKERRQ(ierr);
189       ierr = PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_quad_lin_t, ctx);CHKERRQ(ierr);
190       ierr = 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);CHKERRQ(ierr);
191       break;
192     case SOL_QUADRATIC_TRIG:
193       ierr = PetscDSSetResidual(ds, 0, f0_quad_trig, f1_temp);CHKERRQ(ierr);
194       ierr = PetscDSSetExactSolution(ds, 0, mms_quad_trig, ctx);CHKERRQ(ierr);
195       ierr = PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_quad_trig_t, ctx);CHKERRQ(ierr);
196       ierr = 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);CHKERRQ(ierr);
197       break;
198     case SOL_TRIG_LINEAR:
199       ierr = PetscDSSetResidual(ds, 0, f0_trig_lin,  f1_temp);CHKERRQ(ierr);
200       ierr = PetscDSSetExactSolution(ds, 0, mms_trig_lin, ctx);CHKERRQ(ierr);
201       ierr = PetscDSSetExactSolutionTimeDerivative(ds, 0, mms_trig_lin_t, ctx);CHKERRQ(ierr);
202       ierr = 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);CHKERRQ(ierr);
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   PetscFunctionReturn(0);
207 }
208 
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;
216   PetscErrorCode ierr;
217 
218   PetscFunctionBeginUser;
219   ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr);
220   ierr = DMPlexGetHeightStratum(dm, 0, &cStart, NULL);CHKERRQ(ierr);
221   ierr = DMPlexGetCellType(dm, cStart, &ct);CHKERRQ(ierr);
222   simplex = DMPolytopeTypeGetNumVertices(ct) == DMPolytopeTypeGetDim(ct)+1 ? PETSC_TRUE : PETSC_FALSE;
223   /* Create finite element */
224   ierr = PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, "temp_", -1, &fe);CHKERRQ(ierr);
225   ierr = PetscObjectSetName((PetscObject) fe, "temperature");CHKERRQ(ierr);
226   /* Set discretization and boundary conditions for each mesh */
227   ierr = DMSetField(dm, 0, NULL, (PetscObject) fe);CHKERRQ(ierr);
228   ierr = DMCreateDS(dm);CHKERRQ(ierr);
229   ierr = SetupProblem(dm, ctx);CHKERRQ(ierr);
230   while (cdm) {
231     ierr = DMCopyDisc(dm, cdm);CHKERRQ(ierr);
232     ierr = DMGetCoarseDM(cdm, &cdm);CHKERRQ(ierr);
233   }
234   ierr = PetscFEDestroy(&fe);CHKERRQ(ierr);
235   PetscFunctionReturn(0);
236 }
237 
238 static PetscErrorCode SetInitialConditions(TS ts, Vec u)
239 {
240   DM             dm;
241   PetscReal      t;
242   PetscErrorCode ierr;
243 
244   PetscFunctionBegin;
245   ierr = TSGetDM(ts, &dm);CHKERRQ(ierr);
246   ierr = TSGetTime(ts, &t);CHKERRQ(ierr);
247   ierr = DMComputeExactSolution(dm, t, u, NULL);CHKERRQ(ierr);
248   PetscFunctionReturn(0);
249 }
250 
251 int main(int argc, char **argv)
252 {
253   DM             dm;
254   TS             ts;
255   Vec            u;
256   AppCtx         ctx;
257   PetscErrorCode ierr;
258 
259   ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr;
260   ierr = ProcessOptions(PETSC_COMM_WORLD, &ctx);CHKERRQ(ierr);
261   ierr = CreateMesh(PETSC_COMM_WORLD, &dm, &ctx);CHKERRQ(ierr);
262   ierr = DMSetApplicationContext(dm, &ctx);CHKERRQ(ierr);
263   ierr = SetupDiscretization(dm, &ctx);CHKERRQ(ierr);
264 
265   ierr = TSCreate(PETSC_COMM_WORLD, &ts);CHKERRQ(ierr);
266   ierr = TSSetDM(ts, dm);CHKERRQ(ierr);
267   ierr = DMTSSetBoundaryLocal(dm, DMPlexTSComputeBoundary, &ctx);CHKERRQ(ierr);
268   ierr = DMTSSetIFunctionLocal(dm, DMPlexTSComputeIFunctionFEM, &ctx);CHKERRQ(ierr);
269   ierr = DMTSSetIJacobianLocal(dm, DMPlexTSComputeIJacobianFEM, &ctx);CHKERRQ(ierr);
270   ierr = TSSetExactFinalTime(ts, TS_EXACTFINALTIME_MATCHSTEP);CHKERRQ(ierr);
271   ierr = TSSetFromOptions(ts);CHKERRQ(ierr);
272   ierr = TSSetComputeInitialCondition(ts, SetInitialConditions);CHKERRQ(ierr);
273 
274   ierr = DMCreateGlobalVector(dm, &u);CHKERRQ(ierr);
275   ierr = DMTSCheckFromOptions(ts, u);CHKERRQ(ierr);
276   ierr = SetInitialConditions(ts, u);CHKERRQ(ierr);
277   ierr = PetscObjectSetName((PetscObject) u, "temperature");CHKERRQ(ierr);
278   ierr = TSSolve(ts, u);CHKERRQ(ierr);
279   ierr = DMTSCheckFromOptions(ts, u);CHKERRQ(ierr);
280 
281   ierr = VecDestroy(&u);CHKERRQ(ierr);
282   ierr = TSDestroy(&ts);CHKERRQ(ierr);
283   ierr = DMDestroy(&dm);CHKERRQ(ierr);
284   ierr = PetscFinalize();
285   return ierr;
286 }
287 
288 /*TEST
289 
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
352 
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
415 
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
424 
425 TEST*/
426