Actual source code: ex36.c
1: static char help[] = "First example in homogenization book\n\n";
3: #include <petscsnes.h>
4: #include <petscdmplex.h>
5: #include <petscds.h>
6: #include <petscconvest.h>
7: #include <petscbag.h>
9: /*
10: To control the refinement use -dm_plex_box_faces <n> or -dm_refine <k>, or both
12: To see the exact and computed solutions
14: -compare_view draw -draw_size 500,500 -draw_pause -1
16: To see the delay in convergence of the discretization use
18: -snes_convergence_estimate -convest_num_refine 7 -convest_monitor
20: and to see the proper rate use
22: -dm_refine 5 -snes_convergence_estimate -convest_num_refine 2 -convest_monitor
23: */
25: typedef enum {
26: MOD_CONSTANT,
27: MOD_OSCILLATORY,
28: MOD_TANH,
29: NUM_MOD_TYPES
30: } ModType;
31: const char *modTypes[NUM_MOD_TYPES + 1] = {"constant", "oscillatory", "tanh", "unknown"};
33: /* Constants */
34: enum {
35: EPSILON,
36: NUM_CONSTANTS
37: };
39: typedef struct {
40: PetscReal epsilon; /* Wavelength of fine scale oscillation */
41: } Parameter;
43: typedef struct {
44: PetscBag bag; /* Holds problem parameters */
45: ModType modType; /* Model type */
46: } AppCtx;
48: static PetscErrorCode trig_homogeneous_u(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
49: {
50: PetscInt d;
51: *u = 1.0;
52: for (d = 0; d < dim; ++d) *u *= PetscSinReal(2.0 * PETSC_PI * x[d]);
53: return 0;
54: }
56: static PetscErrorCode oscillatory_u(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
57: {
58: Parameter *param = (Parameter *)ctx;
59: const PetscReal eps = param->epsilon;
61: u[0] = x[0] - x[0] * x[0] + (eps / (2. * PETSC_PI)) * (0.5 - x[0]) * PetscSinReal(2. * PETSC_PI * x[0] / eps) + PetscSqr(eps / (2. * PETSC_PI)) * (1. - PetscCosReal(2. * PETSC_PI * x[0] / eps));
62: return 0;
63: }
65: static void f0_trig_homogeneous_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[])
66: {
67: PetscInt d;
68: for (d = 0; d < dim; ++d) {
69: PetscScalar v = 1.;
70: for (PetscInt e = 0; e < dim; e++) {
71: if (e == d) {
72: v *= -4.0 * PetscSqr(PETSC_PI) * PetscSinReal(2.0 * PETSC_PI * x[d]);
73: } else {
74: v *= PetscSinReal(2.0 * PETSC_PI * x[d]);
75: }
76: }
77: f0[0] += v;
78: }
79: }
81: 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[])
82: {
83: PetscInt d;
84: for (d = 0; d < dim; ++d) f1[d] = u_x[d];
85: }
87: 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[])
88: {
89: PetscInt d;
90: for (d = 0; d < dim; ++d) g3[d * dim + d] = 1.0;
91: }
93: static void f0_oscillatory_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[])
94: {
95: f0[0] = -1.;
96: }
98: static void f1_oscillatory_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[])
99: {
100: const PetscReal eps = PetscRealPart(constants[EPSILON]);
102: f1[0] = u_x[0] / (2. + PetscCosReal(2. * PETSC_PI * x[0] / eps));
103: }
105: static void g3_oscillatory_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[])
106: {
107: const PetscReal eps = PetscRealPart(constants[EPSILON]);
109: g3[0] = 1. / (2. + PetscCosReal(2. * PETSC_PI * x[0] / eps));
110: }
112: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
113: {
114: PetscInt mod;
117: options->modType = MOD_CONSTANT;
118: PetscOptionsBegin(comm, "", "Homogenization Problem Options", "DMPLEX");
119: mod = options->modType;
120: PetscOptionsEList("-mod_type", "The model type", "ex36.c", modTypes, NUM_MOD_TYPES, modTypes[options->modType], &mod, NULL);
121: options->modType = (ModType)mod;
122: PetscOptionsEnd();
123: return 0;
124: }
126: static PetscErrorCode SetupParameters(MPI_Comm comm, AppCtx *user)
127: {
128: PetscBag bag;
129: Parameter *p;
132: PetscBagCreate(comm, sizeof(Parameter), &user->bag);
133: PetscBagGetData(user->bag, (void **)&p);
134: PetscBagSetName(user->bag, "par", "Homogenization parameters");
135: bag = user->bag;
136: PetscBagRegisterReal(bag, &p->epsilon, 1.0, "epsilon", "Wavelength of fine scale oscillation");
137: return 0;
138: }
140: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
141: {
143: DMCreate(comm, dm);
144: DMSetType(*dm, DMPLEX);
145: DMSetFromOptions(*dm);
146: DMSetApplicationContext(*dm, user);
147: DMViewFromOptions(*dm, NULL, "-dm_view");
148: return 0;
149: }
151: static PetscErrorCode SetupPrimalProblem(DM dm, AppCtx *user)
152: {
153: PetscDS ds;
154: DMLabel label;
155: PetscSimplePointFunc ex;
156: const PetscInt id = 1;
157: void *ctx;
160: DMGetDS(dm, &ds);
161: PetscBagGetData(user->bag, (void **)&ctx);
162: switch (user->modType) {
163: case MOD_CONSTANT:
164: PetscDSSetResidual(ds, 0, f0_trig_homogeneous_u, f1_u);
165: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
166: DMGetLabel(dm, "marker", &label);
167: ex = trig_homogeneous_u;
168: DMAddBoundary(dm, DM_BC_ESSENTIAL, "wall", label, 1, &id, 0, 0, NULL, (void (*)(void))ex, NULL, ctx, NULL);
169: break;
170: case MOD_OSCILLATORY:
171: PetscDSSetResidual(ds, 0, f0_oscillatory_u, f1_oscillatory_u);
172: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_oscillatory_uu);
173: DMGetLabel(dm, "marker", &label);
174: ex = oscillatory_u;
175: DMAddBoundary(dm, DM_BC_ESSENTIAL, "wall", label, 1, &id, 0, 0, NULL, (void (*)(void))ex, NULL, ctx, NULL);
176: break;
177: default:
178: SETERRQ(PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_WRONG, "Unsupported model type: %s (%d)", modTypes[PetscMin(user->modType, NUM_MOD_TYPES)], user->modType);
179: }
180: PetscDSSetExactSolution(ds, 0, ex, ctx);
181: /* Setup constants */
182: {
183: Parameter *param;
184: PetscScalar constants[NUM_CONSTANTS];
186: PetscBagGetData(user->bag, (void **)¶m);
188: constants[EPSILON] = param->epsilon;
189: PetscDSSetConstants(ds, NUM_CONSTANTS, constants);
190: }
191: return 0;
192: }
194: static PetscErrorCode SetupDiscretization(DM dm, const char name[], PetscErrorCode (*setup)(DM, AppCtx *), AppCtx *user)
195: {
196: DM cdm = dm;
197: PetscFE fe;
198: PetscBool simplex;
199: PetscInt dim;
200: char prefix[PETSC_MAX_PATH_LEN];
203: DMGetDimension(dm, &dim);
204: DMPlexIsSimplex(dm, &simplex);
205: /* Create finite element */
206: PetscSNPrintf(prefix, PETSC_MAX_PATH_LEN, "%s_", name);
207: PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, name ? prefix : NULL, -1, &fe);
208: PetscObjectSetName((PetscObject)fe, name);
209: /* Set discretization and boundary conditions for each mesh */
210: DMSetField(dm, 0, NULL, (PetscObject)fe);
211: DMCreateDS(dm);
212: (*setup)(dm, user);
213: while (cdm) {
214: DMCopyDisc(dm, cdm);
215: DMGetCoarseDM(cdm, &cdm);
216: }
217: PetscFEDestroy(&fe);
218: return 0;
219: }
221: static PetscErrorCode CompareView(Vec u)
222: {
223: DM dm;
224: Vec v[2], lv[2], exact;
225: PetscOptions options;
226: PetscViewer viewer;
227: PetscViewerFormat format;
228: PetscBool flg;
229: PetscInt i;
230: const char *name, *prefix;
233: VecGetDM(u, &dm);
234: PetscObjectGetOptions((PetscObject)dm, &options);
235: PetscObjectGetOptionsPrefix((PetscObject)dm, &prefix);
236: PetscOptionsGetViewer(PetscObjectComm((PetscObject)dm), options, prefix, "-compare_view", &viewer, &format, &flg);
237: if (flg) {
238: DMGetGlobalVector(dm, &exact);
239: DMComputeExactSolution(dm, 0.0, exact, NULL);
240: v[0] = u;
241: v[1] = exact;
242: for (i = 0; i < 2; ++i) {
243: DMGetLocalVector(dm, &lv[i]);
244: PetscObjectGetName((PetscObject)v[i], &name);
245: PetscObjectSetName((PetscObject)lv[i], name);
246: DMGlobalToLocalBegin(dm, v[i], INSERT_VALUES, lv[i]);
247: DMGlobalToLocalEnd(dm, v[i], INSERT_VALUES, lv[i]);
248: DMPlexInsertBoundaryValues(dm, PETSC_TRUE, lv[i], 0., NULL, NULL, NULL);
249: }
250: DMPlexVecView1D(dm, 2, lv, viewer);
251: for (i = 0; i < 2; ++i) DMRestoreLocalVector(dm, &lv[i]);
252: DMRestoreGlobalVector(dm, &exact);
253: PetscViewerDestroy(&viewer);
254: }
255: return 0;
256: }
258: typedef struct {
259: Mat Mcoarse; /* Mass matrix on the coarse space */
260: Mat Mfine; /* Mass matrix on the fine space */
261: Mat Ifine; /* Interpolator from coarse to fine */
262: Vec Iscale; /* Scaling vector for restriction */
263: KSP kspCoarse; /* Solver for the coarse mass matrix */
264: Vec tmpfine; /* Temporary vector in the fine space */
265: Vec tmpcoarse; /* Temporary vector in the coarse space */
266: } ProjStruct;
268: static PetscErrorCode DestroyCoarseProjection(Mat Pi)
269: {
270: ProjStruct *ctx;
272: MatShellGetContext(Pi, (void **)&ctx);
273: MatDestroy(&ctx->Mcoarse);
274: MatDestroy(&ctx->Mfine);
275: MatDestroy(&ctx->Ifine);
276: VecDestroy(&ctx->Iscale);
277: KSPDestroy(&ctx->kspCoarse);
278: VecDestroy(&ctx->tmpcoarse);
279: VecDestroy(&ctx->tmpfine);
280: PetscFree(ctx);
281: MatShellSetContext(Pi, NULL);
282: return 0;
283: }
285: static PetscErrorCode CoarseProjection(Mat Pi, Vec x, Vec y)
286: {
287: ProjStruct *ctx;
289: MatShellGetContext(Pi, (void **)&ctx);
290: MatMult(ctx->Mfine, x, ctx->tmpfine);
291: PetscObjectSetName((PetscObject)ctx->tmpfine, "Fine DG RHS");
292: PetscObjectSetOptionsPrefix((PetscObject)ctx->tmpfine, "fine_dg_");
293: VecViewFromOptions(ctx->tmpfine, NULL, "-rhs_view");
294: MatMultTranspose(ctx->Ifine, ctx->tmpfine, ctx->tmpcoarse);
295: PetscObjectSetName((PetscObject)ctx->tmpcoarse, "Coarse DG RHS");
296: PetscObjectSetOptionsPrefix((PetscObject)ctx->tmpcoarse, "coarse_dg_");
297: VecViewFromOptions(ctx->tmpcoarse, NULL, "-rhs_view");
298: VecPointwiseMult(ctx->tmpcoarse, ctx->Iscale, ctx->tmpcoarse);
299: VecViewFromOptions(ctx->tmpcoarse, NULL, "-rhs_view");
300: KSPSolve(ctx->kspCoarse, ctx->tmpcoarse, y);
301: return 0;
302: }
304: static PetscErrorCode CreateCoarseProjection(DM dmc, DM dmf, Mat *Pi)
305: {
306: ProjStruct *ctx;
307: PetscInt m, n, M, N;
309: PetscMalloc1(1, &ctx);
310: DMCreateGlobalVector(dmc, &ctx->tmpcoarse);
311: DMCreateGlobalVector(dmf, &ctx->tmpfine);
312: VecGetLocalSize(ctx->tmpcoarse, &m);
313: VecGetSize(ctx->tmpcoarse, &M);
314: VecGetLocalSize(ctx->tmpfine, &n);
315: VecGetSize(ctx->tmpfine, &N);
316: DMCreateMassMatrix(dmc, dmc, &ctx->Mcoarse);
317: DMCreateMassMatrix(dmf, dmf, &ctx->Mfine);
318: DMCreateInterpolation(dmc, dmf, &ctx->Ifine, &ctx->Iscale);
319: KSPCreate(PetscObjectComm((PetscObject)dmc), &ctx->kspCoarse);
320: PetscObjectSetOptionsPrefix((PetscObject)ctx->kspCoarse, "coarse_");
321: KSPSetOperators(ctx->kspCoarse, ctx->Mcoarse, ctx->Mcoarse);
322: KSPSetFromOptions(ctx->kspCoarse);
323: MatCreateShell(PetscObjectComm((PetscObject)dmc), m, n, M, N, ctx, Pi);
324: MatShellSetOperation(*Pi, MATOP_DESTROY, (void (*)(void))DestroyCoarseProjection);
325: MatShellSetOperation(*Pi, MATOP_MULT, (void (*)(void))CoarseProjection);
326: return 0;
327: }
329: typedef struct {
330: Mat Ifdg; /* Embed the fine space into the DG version */
331: Mat Pi; /* The L_2 stable projection to the DG coarse space */
332: Vec tmpc; /* A temporary vector in the DG coarse space */
333: Vec tmpf; /* A temporary vector in the DG fine space */
334: } QuasiInterp;
336: static PetscErrorCode DestroyQuasiInterpolator(Mat P)
337: {
338: QuasiInterp *ctx;
340: MatShellGetContext(P, (void **)&ctx);
341: MatDestroy(&ctx->Ifdg);
342: MatDestroy(&ctx->Pi);
343: VecDestroy(&ctx->tmpc);
344: VecDestroy(&ctx->tmpf);
345: PetscFree(ctx);
346: MatShellSetContext(P, NULL);
347: return 0;
348: }
350: static PetscErrorCode QuasiInterpolate(Mat P, Vec x, Vec y)
351: {
352: QuasiInterp *ctx;
353: DM dmcdg, dmc;
354: Vec ly;
356: MatShellGetContext(P, (void **)&ctx);
357: MatMult(ctx->Ifdg, x, ctx->tmpf);
359: PetscObjectSetName((PetscObject)ctx->tmpf, "Fine DG Potential");
360: PetscObjectSetOptionsPrefix((PetscObject)ctx->tmpf, "fine_dg_");
361: VecViewFromOptions(ctx->tmpf, NULL, "-vec_view");
362: MatMult(ctx->Pi, x, ctx->tmpc);
364: PetscObjectSetName((PetscObject)ctx->tmpc, "Coarse DG Potential");
365: PetscObjectSetOptionsPrefix((PetscObject)ctx->tmpc, "coarse_dg_");
366: VecViewFromOptions(ctx->tmpc, NULL, "-vec_view");
367: VecGetDM(ctx->tmpc, &dmcdg);
369: VecGetDM(y, &dmc);
370: DMGetLocalVector(dmc, &ly);
371: DMPlexComputeClementInterpolant(dmcdg, ctx->tmpc, ly);
372: DMLocalToGlobalBegin(dmc, ly, INSERT_VALUES, y);
373: DMLocalToGlobalEnd(dmc, ly, INSERT_VALUES, y);
374: DMRestoreLocalVector(dmc, &ly);
375: return 0;
376: }
378: static PetscErrorCode CreateQuasiInterpolator(DM dmc, DM dmf, Mat *P)
379: {
380: QuasiInterp *ctx;
381: DM dmcdg, dmfdg;
382: PetscFE fe;
383: Vec x, y;
384: DMPolytopeType ct;
385: PetscInt dim, cStart, m, n, M, N;
387: PetscCalloc1(1, &ctx);
388: DMGetGlobalVector(dmc, &x);
389: DMGetGlobalVector(dmf, &y);
390: VecGetLocalSize(x, &m);
391: VecGetSize(x, &M);
392: VecGetLocalSize(y, &n);
393: VecGetSize(y, &N);
394: DMRestoreGlobalVector(dmc, &x);
395: DMRestoreGlobalVector(dmf, &y);
397: DMClone(dmf, &dmfdg);
398: DMGetDimension(dmfdg, &dim);
399: DMPlexGetHeightStratum(dmfdg, 0, &cStart, NULL);
400: DMPlexGetCellType(dmfdg, cStart, &ct);
401: PetscFECreateByCell(PETSC_COMM_SELF, dim, 1, ct, "fine_dg_", PETSC_DETERMINE, &fe);
402: DMSetField(dmfdg, 0, NULL, (PetscObject)fe);
403: PetscFEDestroy(&fe);
404: DMCreateDS(dmfdg);
405: DMCreateInterpolation(dmf, dmfdg, &ctx->Ifdg, NULL);
406: DMCreateGlobalVector(dmfdg, &ctx->tmpf);
407: DMDestroy(&dmfdg);
409: DMClone(dmc, &dmcdg);
410: DMGetDimension(dmcdg, &dim);
411: DMPlexGetHeightStratum(dmcdg, 0, &cStart, NULL);
412: DMPlexGetCellType(dmcdg, cStart, &ct);
413: PetscFECreateByCell(PETSC_COMM_SELF, dim, 1, ct, "coarse_dg_", PETSC_DETERMINE, &fe);
414: DMSetField(dmcdg, 0, NULL, (PetscObject)fe);
415: PetscFEDestroy(&fe);
416: DMCreateDS(dmcdg);
418: CreateCoarseProjection(dmcdg, dmf, &ctx->Pi);
419: DMCreateGlobalVector(dmcdg, &ctx->tmpc);
420: DMDestroy(&dmcdg);
422: MatCreateShell(PetscObjectComm((PetscObject)dmc), m, n, M, N, ctx, P);
423: MatShellSetOperation(*P, MATOP_DESTROY, (void (*)(void))DestroyQuasiInterpolator);
424: MatShellSetOperation(*P, MATOP_MULT, (void (*)(void))QuasiInterpolate);
425: return 0;
426: }
428: static PetscErrorCode CoarseTest(DM dm, Vec u, AppCtx *user)
429: {
430: DM dmc;
431: Mat P; /* The quasi-interpolator to the coarse space */
432: Vec uc;
434: if (user->modType == MOD_CONSTANT) return 0;
435: DMCreate(PetscObjectComm((PetscObject)dm), &dmc);
436: DMSetType(dmc, DMPLEX);
437: PetscObjectSetOptionsPrefix((PetscObject)dmc, "coarse_");
438: DMSetApplicationContext(dmc, user);
439: DMSetFromOptions(dmc);
440: DMViewFromOptions(dmc, NULL, "-dm_view");
442: SetupDiscretization(dmc, "potential", SetupPrimalProblem, user);
443: DMCreateGlobalVector(dmc, &uc);
444: PetscObjectSetName((PetscObject)uc, "potential");
445: PetscObjectSetOptionsPrefix((PetscObject)uc, "coarse_");
447: CreateQuasiInterpolator(dmc, dm, &P);
448: #if 1
449: MatMult(P, u, uc);
450: #else
451: {
452: Mat In;
453: Vec sc;
455: DMCreateInterpolation(dmc, dm, &In, &sc);
456: MatMultTranspose(In, u, uc);
457: VecPointwiseMult(uc, sc, uc);
458: MatDestroy(&In);
459: VecDestroy(&sc);
460: }
461: #endif
462: CompareView(uc);
464: MatDestroy(&P);
465: VecDestroy(&uc);
466: DMDestroy(&dmc);
467: return 0;
468: }
470: int main(int argc, char **argv)
471: {
472: DM dm; /* Problem specification */
473: SNES snes; /* Nonlinear solver */
474: Vec u; /* Solutions */
475: AppCtx user; /* User-defined work context */
478: PetscInitialize(&argc, &argv, NULL, help);
479: ProcessOptions(PETSC_COMM_WORLD, &user);
480: SetupParameters(PETSC_COMM_WORLD, &user);
481: /* Primal system */
482: SNESCreate(PETSC_COMM_WORLD, &snes);
483: CreateMesh(PETSC_COMM_WORLD, &user, &dm);
484: SNESSetDM(snes, dm);
485: SetupDiscretization(dm, "potential", SetupPrimalProblem, &user);
486: DMCreateGlobalVector(dm, &u);
487: VecSet(u, 0.0);
488: PetscObjectSetName((PetscObject)u, "potential");
489: DMPlexSetSNESLocalFEM(dm, &user, &user, &user);
490: SNESSetFromOptions(snes);
491: DMSNESCheckFromOptions(snes, u);
492: SNESSolve(snes, NULL, u);
493: SNESGetSolution(snes, &u);
494: VecViewFromOptions(u, NULL, "-potential_view");
495: CompareView(u);
496: /* Looking at a coarse problem */
497: CoarseTest(dm, u, &user);
498: /* Cleanup */
499: VecDestroy(&u);
500: SNESDestroy(&snes);
501: DMDestroy(&dm);
502: PetscBagDestroy(&user.bag);
503: PetscFinalize();
504: return 0;
505: }
507: /*TEST
509: test:
510: suffix: 1d_p1_constant
511: args: -dm_plex_dim 1 -dm_plex_box_faces 4 -potential_petscspace_degree 1 -dmsnes_check
513: test:
514: suffix: 1d_p1_constant_conv
515: args: -dm_plex_dim 1 -dm_plex_box_faces 4 -potential_petscspace_degree 1 \
516: -snes_convergence_estimate -convest_num_refine 2
518: test:
519: suffix: 1d_p1_oscillatory
520: args: -mod_type oscillatory -epsilon 0.03125 \
521: -dm_plex_dim 1 -dm_plex_box_faces 4 -potential_petscspace_degree 1 -dm_refine 2 -dmsnes_check \
522: -coarse_dm_plex_dim 1 -coarse_dm_plex_box_faces 4 -coarse_dm_plex_hash_location \
523: -fine_dg_petscspace_degree 1 -fine_dg_petscdualspace_lagrange_continuity 0 \
524: -coarse_dg_petscspace_degree 1 -coarse_dg_petscdualspace_lagrange_continuity 0
526: TEST*/