Actual source code: ex12.c
petsc-3.8.4 2018-03-24
1: static char help[] = "Poisson Problem in 2d and 3d with simplicial finite elements.\n\
2: We solve the Poisson problem in a rectangular\n\
3: domain, using a parallel unstructured mesh (DMPLEX) to discretize it.\n\
4: This example supports discretized auxiliary fields (conductivity) as well as\n\
5: multilevel nonlinear solvers.\n\n\n";
7: #include <petscdmplex.h>
8: #include <petscsnes.h>
9: #include <petscds.h>
10: #include <petscviewerhdf5.h>
12: typedef enum {NEUMANN, DIRICHLET, NONE} BCType;
13: typedef enum {RUN_FULL, RUN_EXACT, RUN_TEST, RUN_PERF} RunType;
14: typedef enum {COEFF_NONE, COEFF_ANALYTIC, COEFF_FIELD, COEFF_NONLINEAR} CoeffType;
16: typedef struct {
17: PetscInt debug; /* The debugging level */
18: RunType runType; /* Whether to run tests, or solve the full problem */
19: PetscBool jacobianMF; /* Whether to calculate the Jacobian action on the fly */
20: PetscLogEvent createMeshEvent;
21: PetscBool showInitial, showSolution, restart, check, quiet, nonzInit;
22: /* Domain and mesh definition */
23: PetscInt dim; /* The topological mesh dimension */
24: DMBoundaryType periodicity[3]; /* The domain periodicity */
25: PetscInt cells[3]; /* The initial domain division */
26: char filename[2048]; /* The optional mesh file */
27: PetscBool interpolate; /* Generate intermediate mesh elements */
28: PetscReal refinementLimit; /* The largest allowable cell volume */
29: PetscBool viewHierarchy; /* Whether to view the hierarchy */
30: PetscBool simplex; /* Simplicial mesh */
31: /* Problem definition */
32: BCType bcType;
33: CoeffType variableCoefficient;
34: PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx);
35: PetscBool fieldBC;
36: void (**exactFields)(PetscInt, PetscInt, PetscInt,
37: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
38: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
39: PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]);
40: /* Solver */
41: PC pcmg; /* This is needed for error monitoring */
42: } AppCtx;
44: static PetscErrorCode zero(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
45: {
46: u[0] = 0.0;
47: return 0;
48: }
50: static PetscErrorCode ecks(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
51: {
52: u[0] = x[0];
53: return 0;
54: }
56: /*
57: In 2D for Dirichlet conditions, we use exact solution:
59: u = x^2 + y^2
60: f = 4
62: so that
64: -\Delta u + f = -4 + 4 = 0
66: For Neumann conditions, we have
68: -\nabla u \cdot -\hat y |_{y=0} = (2y)|_{y=0} = 0 (bottom)
69: -\nabla u \cdot \hat y |_{y=1} = -(2y)|_{y=1} = -2 (top)
70: -\nabla u \cdot -\hat x |_{x=0} = (2x)|_{x=0} = 0 (left)
71: -\nabla u \cdot \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)
73: Which we can express as
75: \nabla u \cdot \hat n|_\Gamma = {2 x, 2 y} \cdot \hat n = 2 (x + y)
76: */
77: static PetscErrorCode quadratic_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
78: {
79: *u = x[0]*x[0] + x[1]*x[1];
80: return 0;
81: }
83: static void quadratic_u_field_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
84: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
85: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
86: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
87: {
88: uexact[0] = a[0];
89: }
91: static void f0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
92: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
93: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
94: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
95: {
96: f0[0] = 4.0;
97: }
99: static void f0_bd_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
100: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
101: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
102: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
103: {
104: PetscInt d;
105: for (d = 0, f0[0] = 0.0; d < dim; ++d) f0[0] += -n[d]*2.0*x[d];
106: }
108: static void f1_bd_zero(PetscInt dim, PetscInt Nf, PetscInt NfAux,
109: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
110: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
111: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
112: {
113: PetscInt comp;
114: for (comp = 0; comp < dim; ++comp) f1[comp] = 0.0;
115: }
117: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
118: static void f1_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
119: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
120: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
121: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
122: {
123: PetscInt d;
124: for (d = 0; d < dim; ++d) f1[d] = u_x[d];
125: }
127: /* < \nabla v, \nabla u + {\nabla u}^T >
128: This just gives \nabla u, give the perdiagonal for the transpose */
129: static void g3_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
130: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
131: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
132: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
133: {
134: PetscInt d;
135: for (d = 0; d < dim; ++d) g3[d*dim+d] = 1.0;
136: }
138: /*
139: In 2D for x periodicity and y Dirichlet conditions, we use exact solution:
141: u = sin(2 pi x)
142: f = -4 pi^2 sin(2 pi x)
144: so that
146: -\Delta u + f = 4 pi^2 sin(2 pi x) - 4 pi^2 sin(2 pi x) = 0
147: */
148: static PetscErrorCode xtrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
149: {
150: *u = PetscSinReal(2.0*PETSC_PI*x[0]);
151: return 0;
152: }
154: static void f0_xtrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
155: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
156: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
157: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
158: {
159: f0[0] = -4.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
160: }
162: /*
163: In 2D for x-y periodicity, we use exact solution:
165: u = sin(2 pi x) sin(2 pi y)
166: f = -8 pi^2 sin(2 pi x)
168: so that
170: -\Delta u + f = 4 pi^2 sin(2 pi x) sin(2 pi y) + 4 pi^2 sin(2 pi x) sin(2 pi y) - 8 pi^2 sin(2 pi x) = 0
171: */
172: static PetscErrorCode xytrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
173: {
174: *u = PetscSinReal(2.0*PETSC_PI*x[0])*PetscSinReal(2.0*PETSC_PI*x[1]);
175: return 0;
176: }
178: static void f0_xytrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
179: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
180: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
181: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
182: {
183: f0[0] = -8.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
184: }
186: /*
187: In 2D for Dirichlet conditions with a variable coefficient, we use exact solution:
189: u = x^2 + y^2
190: f = 6 (x + y)
191: nu = (x + y)
193: so that
195: -\div \nu \grad u + f = -6 (x + y) + 6 (x + y) = 0
196: */
197: static PetscErrorCode nu_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
198: {
199: *u = x[0] + x[1];
200: return 0;
201: }
203: void f0_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
204: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
205: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
206: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
207: {
208: f0[0] = 6.0*(x[0] + x[1]);
209: }
211: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
212: void f1_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
213: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
214: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
215: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
216: {
217: PetscInt d;
218: for (d = 0; d < dim; ++d) f1[d] = (x[0] + x[1])*u_x[d];
219: }
221: void f1_field_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
222: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
223: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
224: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
225: {
226: PetscInt d;
227: for (d = 0; d < dim; ++d) f1[d] = a[0]*u_x[d];
228: }
230: /* < \nabla v, \nabla u + {\nabla u}^T >
231: This just gives \nabla u, give the perdiagonal for the transpose */
232: void g3_analytic_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
233: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
234: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
235: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
236: {
237: PetscInt d;
238: for (d = 0; d < dim; ++d) g3[d*dim+d] = x[0] + x[1];
239: }
241: void g3_field_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
242: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
243: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
244: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
245: {
246: PetscInt d;
247: for (d = 0; d < dim; ++d) g3[d*dim+d] = a[0];
248: }
250: /*
251: In 2D for Dirichlet conditions with a nonlinear coefficient (p-Laplacian with p = 4), we use exact solution:
253: u = x^2 + y^2
254: f = 16 (x^2 + y^2)
255: nu = 1/2 |grad u|^2
257: so that
259: -\div \nu \grad u + f = -16 (x^2 + y^2) + 16 (x^2 + y^2) = 0
260: */
261: void f0_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
262: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
263: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
264: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
265: {
266: f0[0] = 16.0*(x[0]*x[0] + x[1]*x[1]);
267: }
269: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
270: void f1_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
271: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
272: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
273: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
274: {
275: PetscScalar nu = 0.0;
276: PetscInt d;
277: for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
278: for (d = 0; d < dim; ++d) f1[d] = 0.5*nu*u_x[d];
279: }
281: /*
282: grad (u + eps w) - grad u = eps grad w
284: 1/2 |grad (u + eps w)|^2 grad (u + eps w) - 1/2 |grad u|^2 grad u
285: = 1/2 (|grad u|^2 + 2 eps <grad u,grad w>) (grad u + eps grad w) - 1/2 |grad u|^2 grad u
286: = 1/2 (eps |grad u|^2 grad w + 2 eps <grad u,grad w> grad u)
287: = eps (1/2 |grad u|^2 grad w + grad u <grad u,grad w>)
288: */
289: void g3_analytic_nonlinear_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
290: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
291: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
292: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
293: {
294: PetscScalar nu = 0.0;
295: PetscInt d, e;
296: for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
297: for (d = 0; d < dim; ++d) {
298: g3[d*dim+d] = 0.5*nu;
299: for (e = 0; e < dim; ++e) {
300: g3[d*dim+e] += u_x[d]*u_x[e];
301: }
302: }
303: }
305: /*
306: In 3D for Dirichlet conditions we use exact solution:
308: u = 2/3 (x^2 + y^2 + z^2)
309: f = 4
311: so that
313: -\Delta u + f = -2/3 * 6 + 4 = 0
315: For Neumann conditions, we have
317: -\nabla u \cdot -\hat z |_{z=0} = (2z)|_{z=0} = 0 (bottom)
318: -\nabla u \cdot \hat z |_{z=1} = -(2z)|_{z=1} = -2 (top)
319: -\nabla u \cdot -\hat y |_{y=0} = (2y)|_{y=0} = 0 (front)
320: -\nabla u \cdot \hat y |_{y=1} = -(2y)|_{y=1} = -2 (back)
321: -\nabla u \cdot -\hat x |_{x=0} = (2x)|_{x=0} = 0 (left)
322: -\nabla u \cdot \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)
324: Which we can express as
326: \nabla u \cdot \hat n|_\Gamma = {2 x, 2 y, 2z} \cdot \hat n = 2 (x + y + z)
327: */
328: static PetscErrorCode quadratic_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
329: {
330: *u = 2.0*(x[0]*x[0] + x[1]*x[1] + x[2]*x[2])/3.0;
331: return 0;
332: }
334: static void quadratic_u_field_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
335: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
336: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
337: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
338: {
339: uexact[0] = a[0];
340: }
342: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
343: {
344: const char *bcTypes[3] = {"neumann", "dirichlet", "none"};
345: const char *runTypes[4] = {"full", "exact", "test", "perf"};
346: const char *coeffTypes[4] = {"none", "analytic", "field", "nonlinear"};
347: PetscInt bd, bc, run, coeff, n;
348: PetscBool flg;
352: options->debug = 0;
353: options->runType = RUN_FULL;
354: options->dim = 2;
355: options->periodicity[0] = DM_BOUNDARY_NONE;
356: options->periodicity[1] = DM_BOUNDARY_NONE;
357: options->periodicity[2] = DM_BOUNDARY_NONE;
358: options->cells[0] = 1;
359: options->cells[1] = 1;
360: options->cells[2] = 1;
361: options->filename[0] = '\0';
362: options->interpolate = PETSC_FALSE;
363: options->refinementLimit = 0.0;
364: options->bcType = DIRICHLET;
365: options->variableCoefficient = COEFF_NONE;
366: options->fieldBC = PETSC_FALSE;
367: options->jacobianMF = PETSC_FALSE;
368: options->showInitial = PETSC_FALSE;
369: options->showSolution = PETSC_FALSE;
370: options->restart = PETSC_FALSE;
371: options->check = PETSC_FALSE;
372: options->viewHierarchy = PETSC_FALSE;
373: options->simplex = PETSC_TRUE;
374: options->quiet = PETSC_FALSE;
375: options->nonzInit = PETSC_FALSE;
377: PetscOptionsBegin(comm, "", "Poisson Problem Options", "DMPLEX");
378: PetscOptionsInt("-debug", "The debugging level", "ex12.c", options->debug, &options->debug, NULL);
379: run = options->runType;
380: PetscOptionsEList("-run_type", "The run type", "ex12.c", runTypes, 3, runTypes[options->runType], &run, NULL);
382: options->runType = (RunType) run;
384: PetscOptionsInt("-dim", "The topological mesh dimension", "ex12.c", options->dim, &options->dim, NULL);
385: bd = options->periodicity[0];
386: PetscOptionsEList("-x_periodicity", "The x-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[0]], &bd, NULL);
387: options->periodicity[0] = (DMBoundaryType) bd;
388: bd = options->periodicity[1];
389: PetscOptionsEList("-y_periodicity", "The y-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[1]], &bd, NULL);
390: options->periodicity[1] = (DMBoundaryType) bd;
391: bd = options->periodicity[2];
392: PetscOptionsEList("-z_periodicity", "The z-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[2]], &bd, NULL);
393: options->periodicity[2] = (DMBoundaryType) bd;
394: n = 3;
395: PetscOptionsIntArray("-cells", "The initial mesh division", "ex12.c", options->cells, &n, NULL);
396: PetscOptionsString("-f", "Mesh filename to read", "ex12.c", options->filename, options->filename, sizeof(options->filename), &flg);
397: PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex12.c", options->interpolate, &options->interpolate, NULL);
398: PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex12.c", options->refinementLimit, &options->refinementLimit, NULL);
399: bc = options->bcType;
400: PetscOptionsEList("-bc_type","Type of boundary condition","ex12.c",bcTypes,3,bcTypes[options->bcType],&bc,NULL);
401: options->bcType = (BCType) bc;
402: coeff = options->variableCoefficient;
403: PetscOptionsEList("-variable_coefficient","Type of variable coefficent","ex12.c",coeffTypes,4,coeffTypes[options->variableCoefficient],&coeff,NULL);
404: options->variableCoefficient = (CoeffType) coeff;
406: PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL);
407: PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL);
408: PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL);
409: PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL);
410: PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL);
411: PetscOptionsBool("-check", "Compare with default integration routines", "ex12.c", options->check, &options->check, NULL);
412: PetscOptionsBool("-dm_view_hierarchy", "View the coarsened hierarchy", "ex12.c", options->viewHierarchy, &options->viewHierarchy, NULL);
413: PetscOptionsBool("-simplex", "Simplicial (true) or tensor (false) mesh", "ex12.c", options->simplex, &options->simplex, NULL);
414: PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL);
415: PetscOptionsBool("-nonzero_initial_guess", "nonzero intial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL);
416: PetscOptionsEnd();
417: PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
418: return(0);
419: }
421: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
422: {
423: DMLabel label;
427: DMCreateLabel(dm, name);
428: DMGetLabel(dm, name, &label);
429: DMPlexMarkBoundaryFaces(dm, label);
430: DMPlexLabelComplete(dm, label);
431: return(0);
432: }
434: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
435: {
436: PetscInt dim = user->dim;
437: const char *filename = user->filename;
438: PetscBool interpolate = user->interpolate;
439: PetscReal refinementLimit = user->refinementLimit;
440: size_t len;
444: PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
445: PetscStrlen(filename, &len);
446: if (!len) {
447: if (user->simplex) {
448: DMPlexCreateBoxMesh(comm, dim, dim == 2 ? 2 : 1, interpolate, dm);
449: } else {
450: PetscInt d;
452: if (user->periodicity[0] || user->periodicity[1] || user->periodicity[2]) for (d = 0; d < dim; ++d) user->cells[d] = PetscMax(user->cells[d], 3);
453: DMPlexCreateHexBoxMesh(comm, dim, user->cells, user->periodicity[0], user->periodicity[1], user->periodicity[2], dm);
454: }
455: PetscObjectSetName((PetscObject) *dm, "Mesh");
456: } else {
457: DMPlexCreateFromFile(comm, filename, interpolate, dm);
458: DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
459: }
460: {
461: PetscPartitioner part;
462: DM refinedMesh = NULL;
463: DM distributedMesh = NULL;
465: /* Refine mesh using a volume constraint */
466: if (refinementLimit > 0.0) {
467: DMPlexSetRefinementLimit(*dm, refinementLimit);
468: DMRefine(*dm, comm, &refinedMesh);
469: if (refinedMesh) {
470: const char *name;
472: PetscObjectGetName((PetscObject) *dm, &name);
473: PetscObjectSetName((PetscObject) refinedMesh, name);
474: DMDestroy(dm);
475: *dm = refinedMesh;
476: }
477: }
478: /* Distribute mesh over processes */
479: DMPlexGetPartitioner(*dm,&part);
480: PetscPartitionerSetFromOptions(part);
481: DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
482: if (distributedMesh) {
483: DMDestroy(dm);
484: *dm = distributedMesh;
485: }
486: }
487: if (user->bcType == NEUMANN) {
488: DMLabel label;
490: DMCreateLabel(*dm, "boundary");
491: DMGetLabel(*dm, "boundary", &label);
492: DMPlexMarkBoundaryFaces(*dm, label);
493: } else if (user->bcType == DIRICHLET) {
494: PetscBool hasLabel;
496: DMHasLabel(*dm,"marker",&hasLabel);
497: if (!hasLabel) {CreateBCLabel(*dm, "marker");}
498: }
499: {
500: char convType[256];
501: PetscBool flg;
503: PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
504: PetscOptionsFList("-dm_plex_convert_type","Convert DMPlex to another format","ex12",DMList,DMPLEX,convType,256,&flg);
505: PetscOptionsEnd();
506: if (flg) {
507: DM dmConv;
509: DMConvert(*dm,convType,&dmConv);
510: if (dmConv) {
511: DMDestroy(dm);
512: *dm = dmConv;
513: }
514: }
515: }
516: DMLocalizeCoordinates(*dm); /* needed for periodic */
517: DMSetFromOptions(*dm);
518: DMViewFromOptions(*dm, NULL, "-dm_view");
519: if (user->viewHierarchy) {
520: DM cdm = *dm;
521: PetscInt i = 0;
522: char buf[256];
524: while (cdm) {
525: DMSetUp(cdm);
526: DMGetCoarseDM(cdm, &cdm);
527: ++i;
528: }
529: cdm = *dm;
530: while (cdm) {
531: PetscViewer viewer;
532: PetscBool isHDF5, isVTK;
534: --i;
535: PetscViewerCreate(comm,&viewer);
536: PetscViewerSetType(viewer,PETSCVIEWERHDF5);
537: PetscViewerSetOptionsPrefix(viewer,"hierarchy_");
538: PetscViewerSetFromOptions(viewer);
539: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&isHDF5);
540: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERVTK,&isVTK);
541: if (isHDF5) {
542: PetscSNPrintf(buf, 256, "ex12-%d.h5", i);
543: }
544: else if (isVTK) {
545: PetscSNPrintf(buf, 256, "ex12-%d.vtu", i);
546: PetscViewerPushFormat(viewer,PETSC_VIEWER_VTK_VTU);
547: }
548: else {
549: PetscSNPrintf(buf, 256, "ex12-%d", i);
550: }
551: PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
552: PetscViewerFileSetName(viewer,buf);
553: DMView(cdm, viewer);
554: PetscViewerDestroy(&viewer);
555: DMGetCoarseDM(cdm, &cdm);
556: }
557: }
558: PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
559: return(0);
560: }
562: static PetscErrorCode SetupProblem(PetscDS prob, AppCtx *user)
563: {
564: const PetscInt id = 1;
568: switch (user->variableCoefficient) {
569: case COEFF_NONE:
570: if (user->periodicity[0]) {
571: if (user->periodicity[1]) {
572: PetscDSSetResidual(prob, 0, f0_xytrig_u, f1_u);
573: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
574: } else {
575: PetscDSSetResidual(prob, 0, f0_xtrig_u, f1_u);
576: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
577: }
578: } else {
579: PetscDSSetResidual(prob, 0, f0_u, f1_u);
580: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
581: }
582: break;
583: case COEFF_ANALYTIC:
584: PetscDSSetResidual(prob, 0, f0_analytic_u, f1_analytic_u);
585: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_analytic_uu);
586: break;
587: case COEFF_FIELD:
588: PetscDSSetResidual(prob, 0, f0_analytic_u, f1_field_u);
589: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_field_uu);
590: break;
591: case COEFF_NONLINEAR:
592: PetscDSSetResidual(prob, 0, f0_analytic_nonlinear_u, f1_analytic_nonlinear_u);
593: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_analytic_nonlinear_uu);
594: break;
595: default: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid variable coefficient type %d", user->variableCoefficient);
596: }
597: switch (user->dim) {
598: case 2:
599: if (user->periodicity[0]) {
600: if (user->periodicity[1]) {
601: user->exactFuncs[0] = xytrig_u_2d;
602: } else {
603: user->exactFuncs[0] = xtrig_u_2d;
604: }
605: } else {
606: user->exactFuncs[0] = quadratic_u_2d;
607: user->exactFields[0] = quadratic_u_field_2d;
608: }
609: if (user->bcType == NEUMANN) {PetscDSSetBdResidual(prob, 0, f0_bd_u, f1_bd_zero);}
610: break;
611: case 3:
612: user->exactFuncs[0] = quadratic_u_3d;
613: user->exactFields[0] = quadratic_u_field_3d;
614: if (user->bcType == NEUMANN) {PetscDSSetBdResidual(prob, 0, f0_bd_u, f1_bd_zero);}
615: break;
616: default:
617: SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %d", user->dim);
618: }
619: PetscDSAddBoundary(prob, user->bcType == DIRICHLET ? (user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL) : DM_BC_NATURAL,
620: "wall", user->bcType == DIRICHLET ? "marker" : "boundary", 0, 0, NULL,
621: user->fieldBC ? (void (*)()) user->exactFields[0] : (void (*)()) user->exactFuncs[0], 1, &id, user);
622: return(0);
623: }
625: static PetscErrorCode SetupMaterial(DM dm, DM dmAux, AppCtx *user)
626: {
627: PetscErrorCode (*matFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {nu_2d};
628: Vec nu;
632: DMCreateLocalVector(dmAux, &nu);
633: DMProjectFunctionLocal(dmAux, 0.0, matFuncs, NULL, INSERT_ALL_VALUES, nu);
634: PetscObjectCompose((PetscObject) dm, "A", (PetscObject) nu);
635: VecDestroy(&nu);
636: return(0);
637: }
639: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
640: {
641: PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
642: Vec uexact;
643: PetscInt dim;
647: DMGetDimension(dm, &dim);
648: if (dim == 2) bcFuncs[0] = quadratic_u_2d;
649: else bcFuncs[0] = quadratic_u_3d;
650: DMCreateLocalVector(dmAux, &uexact);
651: DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
652: PetscObjectCompose((PetscObject) dm, "A", (PetscObject) uexact);
653: VecDestroy(&uexact);
654: return(0);
655: }
657: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
658: {
659: DM cdm = dm;
660: const PetscInt dim = user->dim;
661: PetscFE feAux = NULL;
662: PetscFE feCh = NULL;
663: PetscFE fe;
664: PetscDS prob, probAux = NULL, probCh = NULL;
665: PetscBool simplex = user->simplex;
669: /* Create finite element */
670: PetscFECreateDefault(dm, dim, 1, simplex, NULL, -1, &fe);
671: PetscObjectSetName((PetscObject) fe, "potential");
672: if (user->variableCoefficient == COEFF_FIELD) {
673: PetscQuadrature q;
675: PetscFECreateDefault(dm, dim, 1, simplex, "mat_", -1, &feAux);
676: PetscFEGetQuadrature(fe, &q);
677: PetscFESetQuadrature(feAux, q);
678: PetscDSCreate(PetscObjectComm((PetscObject)dm),&probAux);
679: PetscDSSetDiscretization(probAux, 0, (PetscObject) feAux);
680: } else if (user->fieldBC) {
681: PetscQuadrature q;
683: PetscFECreateDefault(dm, dim, 1, simplex, "bc_", -1, &feAux);
684: PetscFEGetQuadrature(fe, &q);
685: PetscFESetQuadrature(feAux, q);
686: PetscDSCreate(PetscObjectComm((PetscObject)dm),&probAux);
687: PetscDSSetDiscretization(probAux, 0, (PetscObject) feAux);
688: }
689: if (user->check) {
690: PetscFECreateDefault(dm, dim, 1, simplex, "ch_", -1, &feCh);
691: PetscDSCreate(PetscObjectComm((PetscObject)dm),&probCh);
692: PetscDSSetDiscretization(probCh, 0, (PetscObject) feCh);
693: }
694: /* Set discretization and boundary conditions for each mesh */
695: DMGetDS(dm, &prob);
696: PetscDSSetDiscretization(prob, 0, (PetscObject) fe);
697: SetupProblem(prob, user);
698: while (cdm) {
699: DM coordDM;
701: DMSetDS(cdm,prob);
702: DMGetCoordinateDM(cdm,&coordDM);
703: if (feAux) {
704: DM dmAux;
706: DMClone(cdm, &dmAux);
707: DMSetCoordinateDM(dmAux, coordDM);
708: DMSetDS(dmAux, probAux);
709: PetscObjectCompose((PetscObject) dm, "dmAux", (PetscObject) dmAux);
710: if (user->fieldBC) {SetupBC(cdm, dmAux, user);}
711: else {SetupMaterial(cdm, dmAux, user);}
712: DMDestroy(&dmAux);
713: }
714: if (feCh) {
715: DM dmCh;
717: DMClone(cdm, &dmCh);
718: DMSetCoordinateDM(dmCh, coordDM);
719: DMSetDS(dmCh, probCh);
720: PetscObjectCompose((PetscObject) dm, "dmCh", (PetscObject) dmCh);
721: DMDestroy(&dmCh);
722: }
723: if (user->bcType == DIRICHLET) {
724: PetscBool hasLabel;
726: DMHasLabel(cdm, "marker", &hasLabel);
727: if (!hasLabel) {CreateBCLabel(cdm, "marker");}
728: }
729: DMGetCoarseDM(cdm, &cdm);
730: }
731: PetscFEDestroy(&fe);
732: PetscFEDestroy(&feAux);
733: PetscFEDestroy(&feCh);
734: PetscDSDestroy(&probAux);
735: PetscDSDestroy(&probCh);
736: return(0);
737: }
739: #include "petsc/private/petscimpl.h"
741: /*@C
742: KSPMonitorError - Outputs the error at each iteration of an iterative solver.
744: Collective on KSP
746: Input Parameters:
747: + ksp - the KSP
748: . its - iteration number
749: . rnorm - 2-norm, preconditioned residual value (may be estimated).
750: - ctx - monitor context
752: Level: intermediate
754: .keywords: KSP, default, monitor, residual
755: .seealso: KSPMonitorSet(), KSPMonitorTrueResidualNorm(), KSPMonitorDefault()
756: @*/
757: static PetscErrorCode KSPMonitorError(KSP ksp, PetscInt its, PetscReal rnorm, void *ctx)
758: {
759: AppCtx *user = (AppCtx *) ctx;
760: DM dm;
761: Vec du, r;
762: PetscInt level = 0;
763: PetscBool hasLevel;
764: PetscViewer viewer;
765: char buf[256];
769: KSPGetDM(ksp, &dm);
770: /* Calculate solution */
771: {
772: PC pc = user->pcmg; /* The MG PC */
773: DM fdm = NULL, cdm;
774: KSP fksp, cksp;
775: Vec fu, cu;
776: PetscInt levels, l;
778: KSPBuildSolution(ksp, NULL, &du);
779: PetscObjectComposedDataGetInt((PetscObject) ksp, PetscMGLevelId, level, hasLevel);
780: PCMGGetLevels(pc, &levels);
781: PCMGGetSmoother(pc, levels-1, &fksp);
782: KSPBuildSolution(fksp, NULL, &fu);
783: for (l = levels-1; l > level; --l) {
784: Mat R;
785: Vec s;
787: PCMGGetSmoother(pc, l-1, &cksp);
788: KSPGetDM(cksp, &cdm);
789: DMGetGlobalVector(cdm, &cu);
790: PCMGGetRestriction(pc, l, &R);
791: PCMGGetRScale(pc, l, &s);
792: MatRestrict(R, fu, cu);
793: VecPointwiseMult(cu, cu, s);
794: if (l < levels-1) {DMRestoreGlobalVector(fdm, &fu);}
795: fdm = cdm;
796: fu = cu;
797: }
798: if (levels-1 > level) {
799: VecAXPY(du, 1.0, cu);
800: DMRestoreGlobalVector(cdm, &cu);
801: }
802: }
803: /* Calculate error */
804: DMGetGlobalVector(dm, &r);
805: DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
806: VecAXPY(r,-1.0,du);
807: PetscObjectSetName((PetscObject) r, "solution error");
808: /* View error */
809: PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
810: PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
811: VecView(r, viewer);
812: /* Cleanup */
813: PetscViewerDestroy(&viewer);
814: DMRestoreGlobalVector(dm, &r);
815: return(0);
816: }
818: /*@C
819: SNESMonitorError - Outputs the error at each iteration of an iterative solver.
821: Collective on SNES
823: Input Parameters:
824: + snes - the SNES
825: . its - iteration number
826: . rnorm - 2-norm of residual
827: - ctx - user context
829: Level: intermediate
831: .keywords: SNES, nonlinear, default, monitor, norm
832: .seealso: SNESMonitorDefault(), SNESMonitorSet(), SNESMonitorSolution()
833: @*/
834: static PetscErrorCode SNESMonitorError(SNES snes, PetscInt its, PetscReal rnorm, void *ctx)
835: {
836: AppCtx *user = (AppCtx *) ctx;
837: DM dm;
838: Vec u, r;
839: PetscInt level = -1;
840: PetscBool hasLevel;
841: PetscViewer viewer;
842: char buf[256];
846: SNESGetDM(snes, &dm);
847: /* Calculate error */
848: SNESGetSolution(snes, &u);
849: DMGetGlobalVector(dm, &r);
850: PetscObjectSetName((PetscObject) r, "solution error");
851: DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
852: VecAXPY(r, -1.0, u);
853: /* View error */
854: PetscObjectComposedDataGetInt((PetscObject) snes, PetscMGLevelId, level, hasLevel);
855: PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
856: PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
857: VecView(r, viewer);
858: /* Cleanup */
859: PetscViewerDestroy(&viewer);
860: DMRestoreGlobalVector(dm, &r);
861: return(0);
862: }
864: int main(int argc, char **argv)
865: {
866: DM dm; /* Problem specification */
867: SNES snes; /* nonlinear solver */
868: Vec u,r; /* solution, residual vectors */
869: Mat A,J; /* Jacobian matrix */
870: MatNullSpace nullSpace; /* May be necessary for Neumann conditions */
871: AppCtx user; /* user-defined work context */
872: JacActionCtx userJ; /* context for Jacobian MF action */
873: PetscInt its; /* iterations for convergence */
874: PetscReal error = 0.0; /* L_2 error in the solution */
875: PetscBool isFAS;
878: PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
879: ProcessOptions(PETSC_COMM_WORLD, &user);
880: SNESCreate(PETSC_COMM_WORLD, &snes);
881: CreateMesh(PETSC_COMM_WORLD, &user, &dm);
882: SNESSetDM(snes, dm);
883: DMSetApplicationContext(dm, &user);
885: PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
886: SetupDiscretization(dm, &user);
888: DMCreateGlobalVector(dm, &u);
889: PetscObjectSetName((PetscObject) u, "potential");
890: VecDuplicate(u, &r);
892: DMCreateMatrix(dm, &J);
893: if (user.jacobianMF) {
894: PetscInt M, m, N, n;
896: MatGetSize(J, &M, &N);
897: MatGetLocalSize(J, &m, &n);
898: MatCreate(PETSC_COMM_WORLD, &A);
899: MatSetSizes(A, m, n, M, N);
900: MatSetType(A, MATSHELL);
901: MatSetUp(A);
902: #if 0
903: MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction);
904: #endif
906: userJ.dm = dm;
907: userJ.J = J;
908: userJ.user = &user;
910: DMCreateLocalVector(dm, &userJ.u);
911: if (user.fieldBC) {DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);}
912: else {DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);}
913: MatShellSetContext(A, &userJ);
914: } else {
915: A = J;
916: }
917: if (user.bcType == NEUMANN) {
918: MatNullSpaceCreate(PetscObjectComm((PetscObject) dm), PETSC_TRUE, 0, NULL, &nullSpace);
919: MatSetNullSpace(A, nullSpace);
920: }
922: DMPlexSetSNESLocalFEM(dm,&user,&user,&user);
923: SNESSetJacobian(snes, A, J, NULL, NULL);
925: SNESSetFromOptions(snes);
927: if (user.fieldBC) {DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u);}
928: else {DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);}
929: if (user.restart) {
930: #if defined(PETSC_HAVE_HDF5)
931: PetscViewer viewer;
933: PetscViewerCreate(PETSC_COMM_WORLD, &viewer);
934: PetscViewerSetType(viewer, PETSCVIEWERHDF5);
935: PetscViewerFileSetMode(viewer, FILE_MODE_READ);
936: PetscViewerFileSetName(viewer, user.filename);
937: PetscViewerHDF5PushGroup(viewer, "/fields");
938: VecLoad(u, viewer);
939: PetscViewerHDF5PopGroup(viewer);
940: PetscViewerDestroy(&viewer);
941: #endif
942: }
943: if (user.showInitial) {
944: Vec lv;
945: DMGetLocalVector(dm, &lv);
946: DMGlobalToLocalBegin(dm, u, INSERT_VALUES, lv);
947: DMGlobalToLocalEnd(dm, u, INSERT_VALUES, lv);
948: DMPrintLocalVec(dm, "Local function", 1.0e-10, lv);
949: DMRestoreLocalVector(dm, &lv);
950: }
951: if (user.viewHierarchy) {
952: SNES lsnes;
953: KSP ksp;
954: PC pc;
955: PetscInt numLevels, l;
956: PetscBool isMG;
958: PetscObjectTypeCompare((PetscObject) snes, SNESFAS, &isFAS);
959: if (isFAS) {
960: SNESFASGetLevels(snes, &numLevels);
961: for (l = 0; l < numLevels; ++l) {
962: SNESFASGetCycleSNES(snes, l, &lsnes);
963: SNESMonitorSet(lsnes, SNESMonitorError, &user, NULL);
964: }
965: } else {
966: SNESGetKSP(snes, &ksp);
967: KSPGetPC(ksp, &pc);
968: PetscObjectTypeCompare((PetscObject) pc, PCMG, &isMG);
969: if (isMG) {
970: user.pcmg = pc;
971: PCMGGetLevels(pc, &numLevels);
972: for (l = 0; l < numLevels; ++l) {
973: PCMGGetSmootherDown(pc, l, &ksp);
974: KSPMonitorSet(ksp, KSPMonitorError, &user, NULL);
975: }
976: }
977: }
978: }
979: if (user.runType == RUN_FULL || user.runType == RUN_EXACT) {
980: PetscErrorCode (*initialGuess[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {zero};
982: if (user.nonzInit) initialGuess[0] = ecks;
983: if (user.runType == RUN_FULL) {
984: DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
985: }
986: if (user.debug) {
987: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
988: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
989: }
990: SNESSolve(snes, NULL, u);
991: SNESGetIterationNumber(snes, &its);
992: PetscPrintf(PETSC_COMM_WORLD, "Number of SNES iterations = %D\n", its);
993: DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
994: if (error < 1.0e-11) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < 1.0e-11\n");}
995: else {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", error);}
996: if (user.showSolution) {
997: PetscPrintf(PETSC_COMM_WORLD, "Solution\n");
998: VecChop(u, 3.0e-9);
999: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1000: }
1001: } else if (user.runType == RUN_PERF) {
1002: PetscReal res = 0.0;
1004: SNESComputeFunction(snes, u, r);
1005: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1006: VecChop(r, 1.0e-10);
1007: VecNorm(r, NORM_2, &res);
1008: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1009: } else {
1010: PetscReal res = 0.0;
1012: /* Check discretization error */
1013: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1014: if (!user.quiet) {VecView(u, PETSC_VIEWER_STDOUT_WORLD);}
1015: DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
1016: if (error < 1.0e-11) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < 1.0e-11\n");}
1017: else {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", error);}
1018: /* Check residual */
1019: SNESComputeFunction(snes, u, r);
1020: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1021: VecChop(r, 1.0e-10);
1022: if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1023: VecNorm(r, NORM_2, &res);
1024: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1025: /* Check Jacobian */
1026: {
1027: Vec b;
1029: SNESComputeJacobian(snes, u, A, A);
1030: VecDuplicate(u, &b);
1031: VecSet(r, 0.0);
1032: SNESComputeFunction(snes, r, b);
1033: MatMult(A, u, r);
1034: VecAXPY(r, 1.0, b);
1035: VecDestroy(&b);
1036: PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n");
1037: VecChop(r, 1.0e-10);
1038: if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1039: VecNorm(r, NORM_2, &res);
1040: PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", res);
1041: }
1042: }
1043: VecViewFromOptions(u, NULL, "-vec_view");
1045: if (user.bcType == NEUMANN) {MatNullSpaceDestroy(&nullSpace);}
1046: if (user.jacobianMF) {VecDestroy(&userJ.u);}
1047: if (A != J) {MatDestroy(&A);}
1048: MatDestroy(&J);
1049: VecDestroy(&u);
1050: VecDestroy(&r);
1051: SNESDestroy(&snes);
1052: DMDestroy(&dm);
1053: PetscFree2(user.exactFuncs, user.exactFields);
1054: PetscFinalize();
1055: return ierr;
1056: }
1058: /*TEST
1059: build:
1060: requires: hdf5 triangle
1061: # 2D serial P1 test 0-4
1062: test:
1063: suffix: 0
1064: requires: hdf5 triangle
1065: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 0 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1066: test:
1067: suffix: 1
1068: requires: hdf5 triangle
1069: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1070: test:
1071: suffix: 2
1072: requires: hdf5 triangle
1073: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1074: test:
1075: suffix: 3
1076: requires: hdf5 triangle
1077: args: -run_type test -refinement_limit 0.0 -bc_type neumann -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
1078: test:
1079: suffix: 4
1080: requires: hdf5 triangle
1081: args: -run_type test -refinement_limit 0.0625 -bc_type neumann -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1082: # 2D serial P2 test 5-8
1083: test:
1084: suffix: 5
1085: requires: hdf5 triangle
1086: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1087: test:
1088: suffix: 6
1089: requires: hdf5 triangle
1090: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1091: test:
1092: suffix: 7
1093: requires: hdf5 triangle
1094: args: -run_type test -refinement_limit 0.0 -bc_type neumann -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
1095: test:
1096: suffix: 8
1097: requires: hdf5 triangle
1098: args: -run_type test -refinement_limit 0.0625 -bc_type neumann -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
1099: # 3D serial P1 test 9-12
1100: test:
1101: suffix: 9
1102: requires: hdf5 ctetgen
1103: args: -run_type test -dim 3 -refinement_limit 0.0 -bc_type dirichlet -interpolate 0 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view
1104: test:
1105: suffix: 10
1106: requires: hdf5 ctetgen
1107: args: -run_type test -dim 3 -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view
1108: test:
1109: suffix: 11
1110: requires: hdf5 ctetgen
1111: args: -run_type test -dim 3 -refinement_limit 0.0125 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1 -dm_view
1112: test:
1113: suffix: 12
1114: requires: hdf5 ctetgen
1115: args: -run_type test -dim 3 -refinement_limit 0.0 -bc_type neumann -interpolate 1 -petscspace_order 1 -snes_fd -show_initial -dm_plex_print_fem 1 -dm_view
1116: # Analytic variable coefficient 13-20
1117: test:
1118: suffix: 13
1119: requires: hdf5 triangle
1120: args: -run_type test -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1121: test:
1122: suffix: 14
1123: requires: hdf5 triangle
1124: args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1125: test:
1126: suffix: 15
1127: requires: hdf5 triangle
1128: args: -run_type test -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1129: test:
1130: suffix: 16
1131: requires: hdf5 triangle
1132: args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1133: test:
1134: suffix: 17
1135: requires: hdf5 ctetgen
1136: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1137: test:
1138: suffix: 18
1139: requires: hdf5 ctetgen
1140: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient analytic -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1141: test:
1142: suffix: 19
1143: requires: hdf5 ctetgen
1144: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1145: test:
1146: suffix: 20
1147: requires: hdf5 ctetgen
1148: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient analytic -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1149: # P1 variable coefficient 21-28
1150: test:
1151: suffix: 21
1152: requires: hdf5 triangle
1153: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1154: test:
1155: suffix: 22
1156: requires: hdf5 triangle
1157: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1158: test:
1159: suffix: 23
1160: requires: hdf5 triangle
1161: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1162: test:
1163: suffix: 24
1164: requires: hdf5 triangle
1165: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1166: test:
1167: suffix: 25
1168: requires: hdf5 ctetgen
1169: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1170: test:
1171: suffix: 26
1172: requires: hdf5 ctetgen
1173: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field -interpolate 1 -petscspace_order 1 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1174: test:
1175: suffix: 27
1176: requires: hdf5 ctetgen
1177: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1178: test:
1179: suffix: 28
1180: requires: hdf5 ctetgen
1181: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field -interpolate 1 -petscspace_order 2 -mat_petscspace_order 1 -show_initial -dm_plex_print_fem 1
1182: # P0 variable coefficient 29-36
1183: test:
1184: suffix: 29
1185: requires: hdf5 triangle
1186: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1187: test:
1188: suffix: 30
1189: requires: hdf5 triangle
1190: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1191: test:
1192: suffix: 31
1193: requires: hdf5 triangle
1194: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1195: test:
1196: requires: hdf5 triangle
1197: suffix: 32
1198: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1199: test:
1200: requires: hdf5 ctetgen
1201: suffix: 33
1202: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1203: test:
1204: suffix: 34
1205: requires: hdf5 ctetgen
1206: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field -interpolate 1 -petscspace_order 1 -show_initial -dm_plex_print_fem 1
1207: test:
1208: suffix: 35
1209: requires: hdf5 ctetgen
1210: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1211: test:
1212: suffix: 36
1213: requires: hdf5 ctetgen
1214: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field -interpolate 1 -petscspace_order 2 -show_initial -dm_plex_print_fem 1
1215: # Full solve 39-44
1216: test:
1217: suffix: 39
1218: requires: hdf5 triangle
1219: args: -run_type full -refinement_limit 0.015625 -interpolate 1 -petscspace_order 2 -pc_type gamg -ksp_rtol 1.0e-10 -ksp_monitor_short -ksp_converged_reason -snes_monitor_short -snes_converged_reason
1220: test:
1221: suffix: 40
1222: requires: hdf5 triangle
1223: args: -run_type full -refinement_limit 0.015625 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 2 -pc_type svd -ksp_rtol 1.0e-10 -snes_monitor_short -snes_converged_reason
1224: test:
1225: suffix: 41
1226: requires: hdf5 triangle
1227: args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -dm_refine_hierarchy 1 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short
1228: test:
1229: suffix: 42
1230: requires: hdf5 triangle
1231: args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 3 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -dm_refine_hierarchy 2 -dm_plex_print_fem 0 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short -fas_levels_2_snes_type newtonls -fas_levels_2_pc_type svd -fas_levels_2_ksp_rtol 1.0e-10 -fas_levels_2_snes_atol 1.0e-11 -fas_levels_2_snes_monitor_short
1232: test:
1233: suffix: 43
1234: requires: hdf5 triangle
1235: nsize: 2
1236: args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -dm_refine_hierarchy 1 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short
1237: test:
1238: suffix: 44
1239: requires: hdf5 triangle
1240: nsize: 2
1241: args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 3 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -dm_refine_hierarchy 2 -dm_plex_print_fem 0 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short -fas_levels_2_snes_type newtonls -fas_levels_2_pc_type svd -fas_levels_2_ksp_rtol 1.0e-10 -fas_levels_2_snes_atol 1.0e-11 -fas_levels_2_snes_monitor_short
1242: # Restarting
1243: test:
1244: suffix: restart
1245: requires: hdf5 triangle !complex
1246: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_order 1
1247: test:
1248: args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1249: test:
1250: args: -f sol.h5 -restart
1251: # Periodicity
1252: test:
1253: suffix: periodic_0
1254: requires: hdf5 triangle
1255: args: -run_type full -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_order 1
1256: # 2D serial P1 test with field bc
1257: test:
1258: suffix: field_bc_p1_0
1259: requires: hdf5 triangle
1260: args: -run_type test -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1261: test:
1262: suffix: field_bc_p1_1
1263: requires: hdf5 triangle
1264: args: -run_type test -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1265: test:
1266: suffix: field_bc_p1_2
1267: requires: hdf5 triangle
1268: args: -run_type test -interpolate 1 -bc_type neumann -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1269: test:
1270: suffix: field_bc_p1_3
1271: requires: hdf5 triangle
1272: args: -run_type test -dm_refine 1 -interpolate 1 -bc_type neumann -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1273: # 3D serial P1 test with field bc
1274: test:
1275: suffix: field_bc_p1_4
1276: requires: hdf5 ctetgen
1277: args: -run_type test -dim 3 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1278: test:
1279: suffix: field_bc_p1_5
1280: requires: hdf5 ctetgen
1281: args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1282: test:
1283: suffix: field_bc_p1_6
1284: requires: hdf5 ctetgen
1285: args: -run_type test -dim 3 -interpolate 1 -bc_type neumann -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1286: test:
1287: suffix: field_bc_p1_7
1288: requires: hdf5 ctetgen
1289: args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type neumann -field_bc -petscspace_order 1 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1290: # 2D serial P2 test with field bc
1291: test:
1292: suffix: field_bc_p2_0
1293: requires: hdf5 triangle
1294: args: -run_type test -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1295: test:
1296: suffix: field_bc_p2_1
1297: requires: hdf5 triangle
1298: args: -run_type test -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1299: test:
1300: suffix: field_bc_p2_2
1301: requires: hdf5 triangle
1302: args: -run_type test -interpolate 1 -bc_type neumann -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1303: test:
1304: suffix: field_bc_p2_3
1305: requires: hdf5 triangle
1306: args: -run_type test -dm_refine 1 -interpolate 1 -bc_type neumann -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1307: # 3D serial P2 test with field bc
1308: test:
1309: suffix: field_bc_p2_4
1310: requires: hdf5 ctetgen
1311: args: -run_type test -dim 3 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1312: test:
1313: suffix: field_bc_p2_5
1314: requires: hdf5 ctetgen
1315: args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1316: test:
1317: suffix: field_bc_p2_6
1318: requires: hdf5 ctetgen
1319: args: -run_type test -dim 3 -interpolate 1 -bc_type neumann -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1320: test:
1321: suffix: field_bc_p2_7
1322: requires: hdf5 ctetgen
1323: args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type neumann -field_bc -petscspace_order 2 -bc_petscspace_order 2 -show_initial -dm_plex_print_fem 1
1324: # Full solve simplex: Convergence
1325: test:
1326: suffix: tet_conv_p1_r0
1327: requires: hdf5 ctetgen
1328: args: -run_type full -dim 3 -dm_refine 0 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -dm_view -snes_converged_reason -pc_type lu
1329: test:
1330: suffix: tet_conv_p1_r2
1331: requires: hdf5 ctetgen
1332: args: -run_type full -dim 3 -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -dm_view -snes_converged_reason -pc_type lu
1333: test:
1334: suffix: tet_conv_p1_r3
1335: requires: hdf5 ctetgen
1336: args: -run_type full -dim 3 -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -dm_view -snes_converged_reason -pc_type lu
1337: test:
1338: suffix: tet_conv_p2_r0
1339: requires: hdf5 ctetgen
1340: args: -run_type full -dim 3 -dm_refine 0 -bc_type dirichlet -interpolate 1 -petscspace_order 2 -dm_view -snes_converged_reason -pc_type lu
1341: test:
1342: suffix: tet_conv_p2_r2
1343: requires: hdf5 ctetgen
1344: args: -run_type full -dim 3 -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_order 2 -dm_view -snes_converged_reason -pc_type lu
1345: # Full solve simplex: BDDC
1346: test:
1347: suffix: tri_bddc
1348: requires: hdf5 triangle
1349: nsize: 5
1350: args: -run_type full -petscpartitioner_type simple -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -dm_mat_type is -pc_type bddc -snes_monitor_short -ksp_monitor_short -snes_converged_reason -ksp_converged_reason -snes_view -show_solution 0
1351: # Full solve simplex: ASM
1352: test:
1353: suffix: tri_q2q1_asm_lu
1354: requires: hdf5 triangle
1355: args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type asm -pc_asm_type restrict -pc_asm_blocks 4 -sub_pc_type lu -snes_monitor_short -ksp_monitor_short -snes_converged_reason -ksp_converged_reason -snes_view -show_solution 0
1356: test:
1357: suffix: tri_q2q1_msm_lu
1358: requires: hdf5 triangle
1359: args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type asm -pc_asm_type restrict -pc_asm_local_type multiplicative -pc_asm_blocks 4 -sub_pc_type lu -snes_monitor_short -ksp_monitor_short -snes_converged_reason -ksp_converged_reason -snes_view -show_solution 0
1360: test:
1361: suffix: tri_q2q1_asm_sor
1362: requires: hdf5 triangle
1363: args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type asm -pc_asm_type restrict -pc_asm_blocks 4 -sub_pc_type sor -snes_monitor_short -ksp_monitor_short -snes_converged_reason -ksp_converged_reason -snes_view -show_solution 0
1364: test:
1365: suffix: tri_q2q1_msm_sor
1366: requires: hdf5 triangle
1367: args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_order 1 -ksp_type gmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type asm -pc_asm_type restrict -pc_asm_local_type multiplicative -pc_asm_blocks 4 -sub_pc_type sor -snes_monitor_short -ksp_monitor_short -snes_converged_reason -ksp_converged_reason -snes_view -show_solution 0
1368: # Full solve simplex: FAS
1369: test:
1370: suffix: fas_newton_0
1371: requires: hdf5 triangle
1372: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -dm_refine_hierarchy 1 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short
1373: test:
1374: suffix: fas_newton_1
1375: requires: hdf5 triangle
1376: args: -run_type full -dm_refine_hierarchy 3 -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 3 -ksp_rtol 1.0e-10 -fas_coarse_pc_type lu -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -snes_view -fas_levels_snes_type newtonls -fas_levels_snes_linesearch_type basic -fas_levels_ksp_rtol 1.0e-10 -fas_levels_snes_monitor_short
1377: test:
1378: suffix: fas_ngs_0
1379: requires: hdf5 triangle
1380: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -dm_refine_hierarchy 1 -snes_view -fas_levels_1_snes_type ngs -fas_levels_1_snes_monitor_short
1381: test:
1382: suffix: fas_newton_coarse_0
1383: requires: hdf5 pragmatic triangle
1384: TODO: broken
1385: args: -run_type full -dm_refine 2 -dm_plex_hash_location -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -dm_coarsen_hierarchy 1 -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short
1386: test:
1387: suffix: mg_newton_coarse_0
1388: requires: hdf5 triangle pragmatic
1389: args: -run_type full -dm_refine 3 -interpolate 1 -petscspace_order 1 -snes_monitor_short -ksp_monitor_true_residual -snes_linesearch_type basic -snes_converged_reason -dm_coarsen_hierarchy 3 -dm_plex_hash_location -snes_view -dm_view -ksp_type richardson -pc_type mg -pc_mg_levels 4 -snes_atol 1.0e-8 -ksp_atol 1.0e-8 -snes_rtol 0.0 -ksp_rtol 0.0 -ksp_norm_type unpreconditioned -mg_levels_ksp_type gmres -mg_levels_pc_type ilu -mg_levels_ksp_max_it 10
1390: test:
1391: suffix: mg_newton_coarse_1
1392: requires: hdf5 triangle pragmatic
1393: args: -run_type full -dm_refine 5 -interpolate 1 -petscspace_order 1 -dm_coarsen_hierarchy 5 -dm_plex_hash_location -dm_plex_separate_marker -dm_plex_coarsen_bd_label marker -dm_plex_remesh_bd -ksp_type richardson -ksp_rtol 1.0e-12 -pc_type mg -pc_mg_levels 3 -mg_levels_ksp_max_it 2 -snes_monitor -ksp_monitor_true_residual -mg_levels_ksp_monitor_true_residual -dm_view -ksp_view
1394: test:
1395: suffix: mg_newton_coarse_2
1396: requires: hdf5 triangle pragmatic
1397: args: -run_type full -dm_refine 5 -interpolate 1 -petscspace_order 1 -dm_coarsen_hierarchy 5 -dm_plex_hash_location -dm_plex_separate_marker -dm_plex_remesh_bd -ksp_type richardson -ksp_rtol 1.0e-12 -pc_type mg -pc_mg_levels 3 -mg_levels_ksp_max_it 2 -snes_monitor -ksp_monitor_true_residual -mg_levels_ksp_monitor_true_residual -dm_view -ksp_view
1398: # Full solve tensor
1399: test:
1400: suffix: tensor_plex_2d
1401: requires: hdf5
1402: args: -run_type test -refinement_limit 0.0 -simplex 0 -bc_type dirichlet -petscspace_order 1 -dm_refine_hierarchy 2 -cells 2,2
1403: test:
1404: suffix: tensor_p4est_2d
1405: requires: hdf5 p4est
1406: args: -run_type test -refinement_limit 0.0 -simplex 0 -bc_type dirichlet -petscspace_order 1 -dm_forest_initial_refinement 2 -dm_forest_minimum_refinement 0 -dm_plex_convert_type p4est -cells 2,2
1407: test:
1408: suffix: tensor_plex_3d
1409: requires: hdf5
1410: args: -run_type test -refinement_limit 0.0 -simplex 0 -bc_type dirichlet -petscspace_order 1 -dim 3 -dm_refine_hierarchy 1 -cells 2,2,2
1411: test:
1412: suffix: tensor_p4est_3d
1413: requires: hdf5 p4est
1414: args: -run_type test -refinement_limit 0.0 -simplex 0 -bc_type dirichlet -petscspace_order 1 -dm_forest_initial_refinement 1 -dm_forest_minimum_refinement 0 -dim 3 -dm_plex_convert_type p8est -cells 2,2,2
1415: # Full solve tensor: AMR
1416: test:
1417: suffix: amr_0
1418: requires: hdf5
1419: nsize: 5
1420: args: -run_type test -petscpartitioner_type simple -refinement_limit 0.0 -simplex 0 -bc_type dirichlet -petscspace_order 1 -dm_refine 1 -cells 2,2
1421: test:
1422: suffix: amr_1
1423: requires: hdf5 p4est !complex
1424: args: -run_type test -refinement_limit 0.0 -simplex 0 -bc_type dirichlet -petscspace_order 1 -dm_plex_convert_type p4est -dm_p4est_refine_pattern center -dm_forest_maximum_refinement 5 -dm_view vtk:amr.vtu:vtk_vtu -vec_view vtk:amr.vtu:vtk_vtu:append -cells 2,2
1425: test:
1426: suffix: p4est_test_q2_conformal_serial
1427: requires: hdf5 p4est
1428: args: -run_type test -interpolate 1 -petscspace_order 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -cells 2,2
1429: test:
1430: suffix: p4est_test_q2_conformal_parallel
1431: requires: hdf5 p4est
1432: nsize: 7
1433: args: -run_type test -interpolate 1 -petscspace_order 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type simple -cells 2,2
1434: test:
1435: suffix: p4est_test_q2_nonconformal_serial
1436: requires: hdf5 p4est
1437: filter: grep -v "CG or CGNE: variant"
1438: args: -run_type test -interpolate 1 -petscspace_order 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2
1439: test:
1440: suffix: p4est_test_q2_nonconformal_parallel
1441: requires: hdf5 p4est
1442: filter: grep -v "CG or CGNE: variant"
1443: nsize: 7
1444: args: -run_type test -interpolate 1 -petscspace_order 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -cells 2,2
1445: test:
1446: suffix: p4est_exact_q2_conformal_serial
1447: requires: hdf5 p4est
1448: args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type none -fas_levels_ksp_type preonly -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -cells 2,2
1449: test:
1450: suffix: p4est_exact_q2_conformal_parallel
1451: requires: hdf5 p4est
1452: nsize: 7
1453: args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type none -fas_levels_ksp_type preonly -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -cells 2,2
1454: test:
1455: suffix: p4est_exact_q2_nonconformal_serial
1456: requires: hdf5 p4est
1457: args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type none -fas_levels_ksp_type preonly -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2
1458: test:
1459: suffix: p4est_exact_q2_nonconformal_parallel
1460: requires: hdf5 p4est
1461: nsize: 7
1462: args: -run_type exact -interpolate 1 -petscspace_order 2 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -pc_type none -ksp_type preonly -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type none -fas_levels_ksp_type preonly -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -cells 2,2
1463: test:
1464: suffix: p4est_full_q2_nonconformal_serial
1465: requires: hdf5 p4est
1466: filter: grep -v "CG or CGNE: variant"
1467: args: -run_type full -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type cg -fas_coarse_pc_type jacobi -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type jacobi -fas_levels_ksp_type cg -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2
1468: test:
1469: suffix: p4est_full_q2_nonconformal_parallel
1470: requires: p4est hdf5
1471: filter: grep -v "CG or CGNE: variant"
1472: nsize: 7
1473: args: -run_type full -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type cg -fas_coarse_pc_type jacobi -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type jacobi -fas_levels_ksp_type cg -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -cells 2,2
1474: test:
1475: suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1476: requires: hdf5 p4est
1477: filter: grep -v "CG or CGNE: variant" | sed -e "s/BDDC: Graph max count: 9223372036854775807/BDDC: Graph max count: 2147483647/g"
1478: nsize: 7
1479: args: -run_type full -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -dm_mat_type is -pc_type bddc -ksp_type cg -fas_coarse_pc_type bddc -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type bddc -fas_levels_ksp_type cg -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -cells 2,2
1480: test:
1481: suffix: p4est_full_q2_nonconformal_parallel_bddc
1482: requires: hdf5 p4est
1483: filter: grep -v "CG or CGNE: variant" | sed -e "s/BDDC: Graph max count: 9223372036854775807/BDDC: Graph max count: 2147483647/g"
1484: nsize: 7
1485: args: -run_type full -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type newtonls -dm_mat_type is -pc_type bddc -ksp_type cg -snes_monitor_short -snes_linesearch_type basic -snes_converged_reason -snes_view -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -petscpartitioner_type simple -cells 2,2
1486: test:
1487: suffix: p4est_fas_q2_conformal_serial
1488: requires: hdf5 p4est broken
1489: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type gmres -fas_coarse_pc_type svd -fas_coarse_ksp_type gmres -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type svd -fas_levels_ksp_type gmres -fas_levels_snes_monitor_short -simplex 0 -dm_refine_hierarchy 3 -cells 2,2
1490: test:
1491: suffix: p4est_fas_q2_nonconformal_serial
1492: requires: hdf5 p4est broken
1493: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -pc_type jacobi -ksp_type gmres -fas_coarse_pc_type jacobi -fas_coarse_ksp_type gmres -fas_coarse_ksp_monitor_true_residual -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -snes_view -fas_levels_snes_type newtonls -fas_levels_pc_type jacobi -fas_levels_ksp_type gmres -fas_levels_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2
1494: test:
1495: suffix: fas_newton_0_p4est
1496: requires: hdf5 p4est
1497: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_order 1 -snes_type fas -snes_fas_levels 2 -pc_type svd -ksp_rtol 1.0e-10 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -snes_linesearch_type basic -fas_coarse_snes_linesearch_type basic -snes_converged_reason -snes_view -fas_levels_1_snes_type newtonls -fas_levels_1_pc_type svd -fas_levels_1_ksp_rtol 1.0e-10 -fas_levels_1_snes_monitor_short -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2
1499: TEST*/