Actual source code: ex12.c
petsc-3.10.5 2019-03-28
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: /*
8: A visualization of the adaptation can be accomplished using:
10: -dm_adapt_view hdf5:$PWD/adapt.h5 -sol_adapt_view hdf5:$PWD/adapt.h5::append -dm_adapt_pre_view hdf5:$PWD/orig.h5 -sol_adapt_pre_view hdf5:$PWD/orig.h5::append
12: Information on refinement:
14: -info -info_exclude null,sys,vec,is,mat,ksp,snes,ts
15: */
17: #include <petscdmplex.h>
18: #include <petscdmadaptor.h>
19: #include <petscsnes.h>
20: #include <petscds.h>
21: #include <petscviewerhdf5.h>
23: typedef enum {NEUMANN, DIRICHLET, NONE} BCType;
24: typedef enum {RUN_FULL, RUN_EXACT, RUN_TEST, RUN_PERF} RunType;
25: typedef enum {COEFF_NONE, COEFF_ANALYTIC, COEFF_FIELD, COEFF_NONLINEAR, COEFF_CIRCLE, COEFF_CROSS} CoeffType;
27: typedef struct {
28: PetscInt debug; /* The debugging level */
29: RunType runType; /* Whether to run tests, or solve the full problem */
30: PetscBool jacobianMF; /* Whether to calculate the Jacobian action on the fly */
31: PetscLogEvent createMeshEvent;
32: PetscBool showInitial, showSolution, restart, check, quiet, nonzInit;
33: /* Domain and mesh definition */
34: PetscInt dim; /* The topological mesh dimension */
35: DMBoundaryType periodicity[3]; /* The domain periodicity */
36: PetscInt cells[3]; /* The initial domain division */
37: char filename[2048]; /* The optional mesh file */
38: PetscBool interpolate; /* Generate intermediate mesh elements */
39: PetscReal refinementLimit; /* The largest allowable cell volume */
40: PetscBool viewHierarchy; /* Whether to view the hierarchy */
41: PetscBool simplex; /* Simplicial mesh */
42: /* Problem definition */
43: BCType bcType;
44: CoeffType variableCoefficient;
45: PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx);
46: PetscBool fieldBC;
47: void (**exactFields)(PetscInt, PetscInt, PetscInt,
48: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
49: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
50: PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]);
51: PetscBool bdIntegral; /* Compute the integral of the solution on the boundary */
52: /* Solver */
53: PC pcmg; /* This is needed for error monitoring */
54: } AppCtx;
56: static PetscErrorCode zero(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
57: {
58: u[0] = 0.0;
59: return 0;
60: }
62: static PetscErrorCode ecks(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
63: {
64: u[0] = x[0];
65: return 0;
66: }
68: /*
69: In 2D for Dirichlet conditions, we use exact solution:
71: u = x^2 + y^2
72: f = 4
74: so that
76: -\Delta u + f = -4 + 4 = 0
78: For Neumann conditions, we have
80: -\nabla u \cdot -\hat y |_{y=0} = (2y)|_{y=0} = 0 (bottom)
81: -\nabla u \cdot \hat y |_{y=1} = -(2y)|_{y=1} = -2 (top)
82: -\nabla u \cdot -\hat x |_{x=0} = (2x)|_{x=0} = 0 (left)
83: -\nabla u \cdot \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)
85: Which we can express as
87: \nabla u \cdot \hat n|_\Gamma = {2 x, 2 y} \cdot \hat n = 2 (x + y)
89: The boundary integral of this solution is (assuming we are not orienting the edges)
91: \int^1_0 x^2 dx + \int^1_0 (1 + y^2) dy + \int^1_0 (x^2 + 1) dx + \int^1_0 y^2 dy = 1/3 + 4/3 + 4/3 + 1/3 = 3 1/3
92: */
93: static PetscErrorCode quadratic_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
94: {
95: *u = x[0]*x[0] + x[1]*x[1];
96: return 0;
97: }
99: static void quadratic_u_field_2d(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[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
103: {
104: uexact[0] = a[0];
105: }
107: static PetscErrorCode circle_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
108: {
109: const PetscReal alpha = 500.;
110: const PetscReal radius2 = PetscSqr(0.15);
111: const PetscReal r2 = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5);
112: const PetscReal xi = alpha*(radius2 - r2);
114: *u = PetscTanhScalar(xi) + 1.0;
115: return 0;
116: }
118: static PetscErrorCode cross_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
119: {
120: const PetscReal alpha = 50*4;
121: const PetscReal xy = (x[0]-0.5)*(x[1]-0.5);
123: *u = PetscSinScalar(alpha*xy) * (alpha*PetscAbsScalar(xy) < 2*PETSC_PI ? 1 : 0.01);
124: return 0;
125: }
127: static void f0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
128: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
129: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
130: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
131: {
132: f0[0] = 4.0;
133: }
135: static void f0_circle_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
136: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
137: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
138: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
139: {
140: const PetscReal alpha = 500.;
141: const PetscReal radius2 = PetscSqr(0.15);
142: const PetscReal r2 = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5);
143: const PetscReal xi = alpha*(radius2 - r2);
145: f0[0] = (-4.0*alpha - 8.0*PetscSqr(alpha)*r2*PetscTanhReal(xi)) * PetscSqr(1.0/PetscCoshReal(xi));
146: }
148: static void f0_cross_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
149: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
150: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
151: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
152: {
153: const PetscReal alpha = 50*4;
154: const PetscReal xy = (x[0]-0.5)*(x[1]-0.5);
156: f0[0] = PetscSinScalar(alpha*xy) * (alpha*PetscAbsScalar(xy) < 2*PETSC_PI ? 1 : 0.01);
157: }
159: static void f0_bd_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
160: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
161: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
162: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
163: {
164: PetscInt d;
165: for (d = 0, f0[0] = 0.0; d < dim; ++d) f0[0] += -n[d]*2.0*x[d];
166: }
168: static void f1_bd_zero(PetscInt dim, PetscInt Nf, PetscInt NfAux,
169: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
170: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
171: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
172: {
173: PetscInt comp;
174: for (comp = 0; comp < dim; ++comp) f1[comp] = 0.0;
175: }
177: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
178: static void f1_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 f1[])
182: {
183: PetscInt d;
184: for (d = 0; d < dim; ++d) f1[d] = u_x[d];
185: }
187: /* < \nabla v, \nabla u + {\nabla u}^T >
188: This just gives \nabla u, give the perdiagonal for the transpose */
189: static void g3_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
190: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
191: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
192: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
193: {
194: PetscInt d;
195: for (d = 0; d < dim; ++d) g3[d*dim+d] = 1.0;
196: }
198: /*
199: In 2D for x periodicity and y Dirichlet conditions, we use exact solution:
201: u = sin(2 pi x)
202: f = -4 pi^2 sin(2 pi x)
204: so that
206: -\Delta u + f = 4 pi^2 sin(2 pi x) - 4 pi^2 sin(2 pi x) = 0
207: */
208: static PetscErrorCode xtrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
209: {
210: *u = PetscSinReal(2.0*PETSC_PI*x[0]);
211: return 0;
212: }
214: static void f0_xtrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
215: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
216: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
217: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
218: {
219: f0[0] = -4.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
220: }
222: /*
223: In 2D for x-y periodicity, we use exact solution:
225: u = sin(2 pi x) sin(2 pi y)
226: f = -8 pi^2 sin(2 pi x)
228: so that
230: -\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
231: */
232: static PetscErrorCode xytrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
233: {
234: *u = PetscSinReal(2.0*PETSC_PI*x[0])*PetscSinReal(2.0*PETSC_PI*x[1]);
235: return 0;
236: }
238: static void f0_xytrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
239: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
240: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
241: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
242: {
243: f0[0] = -8.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
244: }
246: /*
247: In 2D for Dirichlet conditions with a variable coefficient, we use exact solution:
249: u = x^2 + y^2
250: f = 6 (x + y)
251: nu = (x + y)
253: so that
255: -\div \nu \grad u + f = -6 (x + y) + 6 (x + y) = 0
256: */
257: static PetscErrorCode nu_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
258: {
259: *u = x[0] + x[1];
260: return 0;
261: }
263: void f0_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
264: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
265: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
266: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
267: {
268: f0[0] = 6.0*(x[0] + x[1]);
269: }
271: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
272: void f1_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
273: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
274: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
275: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
276: {
277: PetscInt d;
278: for (d = 0; d < dim; ++d) f1[d] = (x[0] + x[1])*u_x[d];
279: }
281: void f1_field_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
282: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
283: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
284: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
285: {
286: PetscInt d;
287: for (d = 0; d < dim; ++d) f1[d] = a[0]*u_x[d];
288: }
290: /* < \nabla v, \nabla u + {\nabla u}^T >
291: This just gives \nabla u, give the perdiagonal for the transpose */
292: void g3_analytic_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
293: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
294: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
295: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
296: {
297: PetscInt d;
298: for (d = 0; d < dim; ++d) g3[d*dim+d] = x[0] + x[1];
299: }
301: void g3_field_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
302: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
303: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
304: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
305: {
306: PetscInt d;
307: for (d = 0; d < dim; ++d) g3[d*dim+d] = a[0];
308: }
310: /*
311: In 2D for Dirichlet conditions with a nonlinear coefficient (p-Laplacian with p = 4), we use exact solution:
313: u = x^2 + y^2
314: f = 16 (x^2 + y^2)
315: nu = 1/2 |grad u|^2
317: so that
319: -\div \nu \grad u + f = -16 (x^2 + y^2) + 16 (x^2 + y^2) = 0
320: */
321: void f0_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
322: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
323: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
324: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
325: {
326: f0[0] = 16.0*(x[0]*x[0] + x[1]*x[1]);
327: }
329: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
330: void f1_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
331: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
332: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
333: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
334: {
335: PetscScalar nu = 0.0;
336: PetscInt d;
337: for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
338: for (d = 0; d < dim; ++d) f1[d] = 0.5*nu*u_x[d];
339: }
341: /*
342: grad (u + eps w) - grad u = eps grad w
344: 1/2 |grad (u + eps w)|^2 grad (u + eps w) - 1/2 |grad u|^2 grad u
345: = 1/2 (|grad u|^2 + 2 eps <grad u,grad w>) (grad u + eps grad w) - 1/2 |grad u|^2 grad u
346: = 1/2 (eps |grad u|^2 grad w + 2 eps <grad u,grad w> grad u)
347: = eps (1/2 |grad u|^2 grad w + grad u <grad u,grad w>)
348: */
349: void g3_analytic_nonlinear_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
350: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
351: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
352: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
353: {
354: PetscScalar nu = 0.0;
355: PetscInt d, e;
356: for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
357: for (d = 0; d < dim; ++d) {
358: g3[d*dim+d] = 0.5*nu;
359: for (e = 0; e < dim; ++e) {
360: g3[d*dim+e] += u_x[d]*u_x[e];
361: }
362: }
363: }
365: /*
366: In 3D for Dirichlet conditions we use exact solution:
368: u = 2/3 (x^2 + y^2 + z^2)
369: f = 4
371: so that
373: -\Delta u + f = -2/3 * 6 + 4 = 0
375: For Neumann conditions, we have
377: -\nabla u \cdot -\hat z |_{z=0} = (2z)|_{z=0} = 0 (bottom)
378: -\nabla u \cdot \hat z |_{z=1} = -(2z)|_{z=1} = -2 (top)
379: -\nabla u \cdot -\hat y |_{y=0} = (2y)|_{y=0} = 0 (front)
380: -\nabla u \cdot \hat y |_{y=1} = -(2y)|_{y=1} = -2 (back)
381: -\nabla u \cdot -\hat x |_{x=0} = (2x)|_{x=0} = 0 (left)
382: -\nabla u \cdot \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)
384: Which we can express as
386: \nabla u \cdot \hat n|_\Gamma = {2 x, 2 y, 2z} \cdot \hat n = 2 (x + y + z)
387: */
388: static PetscErrorCode quadratic_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
389: {
390: *u = 2.0*(x[0]*x[0] + x[1]*x[1] + x[2]*x[2])/3.0;
391: return 0;
392: }
394: static void quadratic_u_field_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
395: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
396: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
397: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
398: {
399: uexact[0] = a[0];
400: }
402: static void bd_integral_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
403: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
404: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
405: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar *uint)
406: {
407: uint[0] = u[0];
408: }
410: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
411: {
412: const char *bcTypes[3] = {"neumann", "dirichlet", "none"};
413: const char *runTypes[4] = {"full", "exact", "test", "perf"};
414: const char *coeffTypes[6] = {"none", "analytic", "field", "nonlinear", "circle", "cross"};
415: PetscInt bd, bc, run, coeff, n;
416: PetscBool flg;
420: options->debug = 0;
421: options->runType = RUN_FULL;
422: options->dim = 2;
423: options->periodicity[0] = DM_BOUNDARY_NONE;
424: options->periodicity[1] = DM_BOUNDARY_NONE;
425: options->periodicity[2] = DM_BOUNDARY_NONE;
426: options->cells[0] = 2;
427: options->cells[1] = 2;
428: options->cells[2] = 2;
429: options->filename[0] = '\0';
430: options->interpolate = PETSC_FALSE;
431: options->refinementLimit = 0.0;
432: options->bcType = DIRICHLET;
433: options->variableCoefficient = COEFF_NONE;
434: options->fieldBC = PETSC_FALSE;
435: options->jacobianMF = PETSC_FALSE;
436: options->showInitial = PETSC_FALSE;
437: options->showSolution = PETSC_FALSE;
438: options->restart = PETSC_FALSE;
439: options->check = PETSC_FALSE;
440: options->viewHierarchy = PETSC_FALSE;
441: options->simplex = PETSC_TRUE;
442: options->quiet = PETSC_FALSE;
443: options->nonzInit = PETSC_FALSE;
444: options->bdIntegral = PETSC_FALSE;
446: PetscOptionsBegin(comm, "", "Poisson Problem Options", "DMPLEX");
447: PetscOptionsInt("-debug", "The debugging level", "ex12.c", options->debug, &options->debug, NULL);
448: run = options->runType;
449: PetscOptionsEList("-run_type", "The run type", "ex12.c", runTypes, 4, runTypes[options->runType], &run, NULL);
451: options->runType = (RunType) run;
453: PetscOptionsInt("-dim", "The topological mesh dimension", "ex12.c", options->dim, &options->dim, NULL);
454: bd = options->periodicity[0];
455: PetscOptionsEList("-x_periodicity", "The x-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[0]], &bd, NULL);
456: options->periodicity[0] = (DMBoundaryType) bd;
457: bd = options->periodicity[1];
458: PetscOptionsEList("-y_periodicity", "The y-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[1]], &bd, NULL);
459: options->periodicity[1] = (DMBoundaryType) bd;
460: bd = options->periodicity[2];
461: PetscOptionsEList("-z_periodicity", "The z-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[2]], &bd, NULL);
462: options->periodicity[2] = (DMBoundaryType) bd;
463: n = 3;
464: PetscOptionsIntArray("-cells", "The initial mesh division", "ex12.c", options->cells, &n, NULL);
465: PetscOptionsString("-f", "Mesh filename to read", "ex12.c", options->filename, options->filename, sizeof(options->filename), &flg);
466: PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex12.c", options->interpolate, &options->interpolate, NULL);
467: PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex12.c", options->refinementLimit, &options->refinementLimit, NULL);
468: bc = options->bcType;
469: PetscOptionsEList("-bc_type","Type of boundary condition","ex12.c",bcTypes,3,bcTypes[options->bcType],&bc,NULL);
470: options->bcType = (BCType) bc;
471: coeff = options->variableCoefficient;
472: PetscOptionsEList("-variable_coefficient","Type of variable coefficent","ex12.c",coeffTypes,6,coeffTypes[options->variableCoefficient],&coeff,NULL);
473: options->variableCoefficient = (CoeffType) coeff;
475: PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL);
476: PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL);
477: PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL);
478: PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL);
479: PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL);
480: PetscOptionsBool("-check", "Compare with default integration routines", "ex12.c", options->check, &options->check, NULL);
481: PetscOptionsBool("-dm_view_hierarchy", "View the coarsened hierarchy", "ex12.c", options->viewHierarchy, &options->viewHierarchy, NULL);
482: PetscOptionsBool("-simplex", "Simplicial (true) or tensor (false) mesh", "ex12.c", options->simplex, &options->simplex, NULL);
483: PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL);
484: PetscOptionsBool("-nonzero_initial_guess", "nonzero intial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL);
485: PetscOptionsBool("-bd_integral", "Compute the integral of the solution on the boundary", "ex12.c", options->bdIntegral, &options->bdIntegral, NULL);
486: PetscOptionsEnd();
487: PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
488: return(0);
489: }
491: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
492: {
493: DMLabel label;
497: DMCreateLabel(dm, name);
498: DMGetLabel(dm, name, &label);
499: DMPlexMarkBoundaryFaces(dm, 1, label);
500: DMPlexLabelComplete(dm, label);
501: return(0);
502: }
504: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
505: {
506: PetscInt dim = user->dim;
507: const char *filename = user->filename;
508: PetscBool interpolate = user->interpolate;
509: PetscReal refinementLimit = user->refinementLimit;
510: size_t len;
514: PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
515: PetscStrlen(filename, &len);
516: if (!len) {
517: PetscInt d;
519: if (user->periodicity[0] || user->periodicity[1] || user->periodicity[2]) for (d = 0; d < dim; ++d) user->cells[d] = PetscMax(user->cells[d], 3);
520: DMPlexCreateBoxMesh(comm, dim, user->simplex, user->cells, NULL, NULL, user->periodicity, interpolate, dm);
521: PetscObjectSetName((PetscObject) *dm, "Mesh");
522: } else {
523: DMPlexCreateFromFile(comm, filename, interpolate, dm);
524: DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
525: }
526: {
527: PetscPartitioner part;
528: DM refinedMesh = NULL;
529: DM distributedMesh = NULL;
531: /* Refine mesh using a volume constraint */
532: if (refinementLimit > 0.0) {
533: DMPlexSetRefinementLimit(*dm, refinementLimit);
534: DMRefine(*dm, comm, &refinedMesh);
535: if (refinedMesh) {
536: const char *name;
538: PetscObjectGetName((PetscObject) *dm, &name);
539: PetscObjectSetName((PetscObject) refinedMesh, name);
540: DMDestroy(dm);
541: *dm = refinedMesh;
542: }
543: }
544: /* Distribute mesh over processes */
545: DMPlexGetPartitioner(*dm,&part);
546: PetscPartitionerSetFromOptions(part);
547: DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
548: if (distributedMesh) {
549: DMDestroy(dm);
550: *dm = distributedMesh;
551: }
552: }
553: if (user->bcType == NEUMANN) {
554: DMLabel label;
556: DMCreateLabel(*dm, "boundary");
557: DMGetLabel(*dm, "boundary", &label);
558: DMPlexMarkBoundaryFaces(*dm, 1, label);
559: } else if (user->bcType == DIRICHLET) {
560: PetscBool hasLabel;
562: DMHasLabel(*dm,"marker",&hasLabel);
563: if (!hasLabel) {CreateBCLabel(*dm, "marker");}
564: }
565: {
566: char convType[256];
567: PetscBool flg;
569: PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
570: PetscOptionsFList("-dm_plex_convert_type","Convert DMPlex to another format","ex12",DMList,DMPLEX,convType,256,&flg);
571: PetscOptionsEnd();
572: if (flg) {
573: DM dmConv;
575: DMConvert(*dm,convType,&dmConv);
576: if (dmConv) {
577: DMDestroy(dm);
578: *dm = dmConv;
579: }
580: }
581: }
582: DMLocalizeCoordinates(*dm); /* needed for periodic */
583: DMSetFromOptions(*dm);
584: DMViewFromOptions(*dm, NULL, "-dm_view");
585: if (user->viewHierarchy) {
586: DM cdm = *dm;
587: PetscInt i = 0;
588: char buf[256];
590: while (cdm) {
591: DMSetUp(cdm);
592: DMGetCoarseDM(cdm, &cdm);
593: ++i;
594: }
595: cdm = *dm;
596: while (cdm) {
597: PetscViewer viewer;
598: PetscBool isHDF5, isVTK;
600: --i;
601: PetscViewerCreate(comm,&viewer);
602: PetscViewerSetType(viewer,PETSCVIEWERHDF5);
603: PetscViewerSetOptionsPrefix(viewer,"hierarchy_");
604: PetscViewerSetFromOptions(viewer);
605: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&isHDF5);
606: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERVTK,&isVTK);
607: if (isHDF5) {
608: PetscSNPrintf(buf, 256, "ex12-%d.h5", i);
609: } else if (isVTK) {
610: PetscSNPrintf(buf, 256, "ex12-%d.vtu", i);
611: PetscViewerPushFormat(viewer,PETSC_VIEWER_VTK_VTU);
612: } else {
613: PetscSNPrintf(buf, 256, "ex12-%d", i);
614: }
615: PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
616: PetscViewerFileSetName(viewer,buf);
617: DMView(cdm, viewer);
618: PetscViewerDestroy(&viewer);
619: DMGetCoarseDM(cdm, &cdm);
620: }
621: }
622: PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
623: return(0);
624: }
626: static PetscErrorCode SetupProblem(PetscDS prob, AppCtx *user)
627: {
628: const PetscInt id = 1;
632: switch (user->variableCoefficient) {
633: case COEFF_NONE:
634: if (user->periodicity[0]) {
635: if (user->periodicity[1]) {
636: PetscDSSetResidual(prob, 0, f0_xytrig_u, f1_u);
637: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
638: } else {
639: PetscDSSetResidual(prob, 0, f0_xtrig_u, f1_u);
640: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
641: }
642: } else {
643: PetscDSSetResidual(prob, 0, f0_u, f1_u);
644: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
645: }
646: break;
647: case COEFF_ANALYTIC:
648: PetscDSSetResidual(prob, 0, f0_analytic_u, f1_analytic_u);
649: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_analytic_uu);
650: break;
651: case COEFF_FIELD:
652: PetscDSSetResidual(prob, 0, f0_analytic_u, f1_field_u);
653: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_field_uu);
654: break;
655: case COEFF_NONLINEAR:
656: PetscDSSetResidual(prob, 0, f0_analytic_nonlinear_u, f1_analytic_nonlinear_u);
657: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_analytic_nonlinear_uu);
658: break;
659: case COEFF_CIRCLE:
660: PetscDSSetResidual(prob, 0, f0_circle_u, f1_u);
661: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
662: break;
663: case COEFF_CROSS:
664: PetscDSSetResidual(prob, 0, f0_cross_u, f1_u);
665: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
666: break;
667: default: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid variable coefficient type %d", user->variableCoefficient);
668: }
669: switch (user->dim) {
670: case 2:
671: switch (user->variableCoefficient) {
672: case COEFF_CIRCLE:
673: user->exactFuncs[0] = circle_u_2d;break;
674: case COEFF_CROSS:
675: user->exactFuncs[0] = cross_u_2d;break;
676: default:
677: if (user->periodicity[0]) {
678: if (user->periodicity[1]) {
679: user->exactFuncs[0] = xytrig_u_2d;
680: } else {
681: user->exactFuncs[0] = xtrig_u_2d;
682: }
683: } else {
684: user->exactFuncs[0] = quadratic_u_2d;
685: user->exactFields[0] = quadratic_u_field_2d;
686: }
687: }
688: if (user->bcType == NEUMANN) {PetscDSSetBdResidual(prob, 0, f0_bd_u, f1_bd_zero);}
689: break;
690: case 3:
691: user->exactFuncs[0] = quadratic_u_3d;
692: user->exactFields[0] = quadratic_u_field_3d;
693: if (user->bcType == NEUMANN) {PetscDSSetBdResidual(prob, 0, f0_bd_u, f1_bd_zero);}
694: break;
695: default:
696: SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %d", user->dim);
697: }
698: PetscDSAddBoundary(prob, user->bcType == DIRICHLET ? (user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL) : DM_BC_NATURAL,
699: "wall", user->bcType == DIRICHLET ? "marker" : "boundary", 0, 0, NULL,
700: user->fieldBC ? (void (*)(void)) user->exactFields[0] : (void (*)(void)) user->exactFuncs[0], 1, &id, user);
701: PetscDSSetExactSolution(prob, 0, user->exactFuncs[0]);
702: PetscDSSetFromOptions(prob);
703: return(0);
704: }
706: static PetscErrorCode SetupMaterial(DM dm, DM dmAux, AppCtx *user)
707: {
708: PetscErrorCode (*matFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {nu_2d};
709: Vec nu;
713: DMCreateLocalVector(dmAux, &nu);
714: DMProjectFunctionLocal(dmAux, 0.0, matFuncs, NULL, INSERT_ALL_VALUES, nu);
715: PetscObjectCompose((PetscObject) dm, "A", (PetscObject) nu);
716: VecDestroy(&nu);
717: return(0);
718: }
720: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
721: {
722: PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
723: Vec uexact;
724: PetscInt dim;
728: DMGetDimension(dm, &dim);
729: if (dim == 2) bcFuncs[0] = quadratic_u_2d;
730: else bcFuncs[0] = quadratic_u_3d;
731: DMCreateLocalVector(dmAux, &uexact);
732: DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
733: PetscObjectCompose((PetscObject) dm, "A", (PetscObject) uexact);
734: VecDestroy(&uexact);
735: return(0);
736: }
738: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
739: {
740: DM cdm = dm;
741: const PetscInt dim = user->dim;
742: PetscFE feAux = NULL;
743: PetscFE feCh = NULL;
744: PetscFE fe;
745: PetscDS prob, probAux = NULL, probCh = NULL;
746: PetscBool simplex = user->simplex;
747: MPI_Comm comm;
751: /* Create finite element */
752: PetscObjectGetComm((PetscObject) dm, &comm);
753: PetscFECreateDefault(comm, dim, 1, simplex, NULL, -1, &fe);
754: PetscObjectSetName((PetscObject) fe, "potential");
755: if (user->variableCoefficient == COEFF_FIELD) {
756: PetscQuadrature q;
758: PetscFECreateDefault(comm, dim, 1, simplex, "mat_", -1, &feAux);
759: PetscFEGetQuadrature(fe, &q);
760: PetscFESetQuadrature(feAux, q);
761: PetscDSCreate(PetscObjectComm((PetscObject)dm),&probAux);
762: PetscDSSetDiscretization(probAux, 0, (PetscObject) feAux);
763: } else if (user->fieldBC) {
764: PetscQuadrature q;
766: PetscFECreateDefault(comm, dim, 1, simplex, "bc_", -1, &feAux);
767: PetscFEGetQuadrature(fe, &q);
768: PetscFESetQuadrature(feAux, q);
769: PetscDSCreate(PetscObjectComm((PetscObject)dm),&probAux);
770: PetscDSSetDiscretization(probAux, 0, (PetscObject) feAux);
771: }
772: if (user->check) {
773: PetscFECreateDefault(comm, dim, 1, simplex, "ch_", -1, &feCh);
774: PetscDSCreate(PetscObjectComm((PetscObject)dm),&probCh);
775: PetscDSSetDiscretization(probCh, 0, (PetscObject) feCh);
776: }
777: /* Set discretization and boundary conditions for each mesh */
778: DMGetDS(dm, &prob);
779: PetscDSSetDiscretization(prob, 0, (PetscObject) fe);
780: SetupProblem(prob, user);
781: while (cdm) {
782: DM coordDM;
784: DMSetDS(cdm,prob);
785: DMGetCoordinateDM(cdm,&coordDM);
786: if (feAux) {
787: DM dmAux;
789: DMClone(cdm, &dmAux);
790: DMSetCoordinateDM(dmAux, coordDM);
791: DMSetDS(dmAux, probAux);
792: PetscObjectCompose((PetscObject) dm, "dmAux", (PetscObject) dmAux);
793: if (user->fieldBC) {SetupBC(cdm, dmAux, user);}
794: else {SetupMaterial(cdm, dmAux, user);}
795: DMDestroy(&dmAux);
796: }
797: if (feCh) {
798: DM dmCh;
800: DMClone(cdm, &dmCh);
801: DMSetCoordinateDM(dmCh, coordDM);
802: DMSetDS(dmCh, probCh);
803: PetscObjectCompose((PetscObject) dm, "dmCh", (PetscObject) dmCh);
804: DMDestroy(&dmCh);
805: }
806: if (user->bcType == DIRICHLET) {
807: PetscBool hasLabel;
809: DMHasLabel(cdm, "marker", &hasLabel);
810: if (!hasLabel) {CreateBCLabel(cdm, "marker");}
811: }
812: DMGetCoarseDM(cdm, &cdm);
813: }
814: PetscFEDestroy(&fe);
815: PetscFEDestroy(&feAux);
816: PetscFEDestroy(&feCh);
817: PetscDSDestroy(&probAux);
818: PetscDSDestroy(&probCh);
819: return(0);
820: }
822: #include "petsc/private/petscimpl.h"
824: /*@C
825: KSPMonitorError - Outputs the error at each iteration of an iterative solver.
827: Collective on KSP
829: Input Parameters:
830: + ksp - the KSP
831: . its - iteration number
832: . rnorm - 2-norm, preconditioned residual value (may be estimated).
833: - ctx - monitor context
835: Level: intermediate
837: .keywords: KSP, default, monitor, residual
838: .seealso: KSPMonitorSet(), KSPMonitorTrueResidualNorm(), KSPMonitorDefault()
839: @*/
840: static PetscErrorCode KSPMonitorError(KSP ksp, PetscInt its, PetscReal rnorm, void *ctx)
841: {
842: AppCtx *user = (AppCtx *) ctx;
843: DM dm;
844: Vec du = NULL, r;
845: PetscInt level = 0;
846: PetscBool hasLevel;
847: #if defined(PETSC_HAVE_HDF5)
848: PetscViewer viewer;
849: char buf[256];
850: #endif
854: KSPGetDM(ksp, &dm);
855: /* Calculate solution */
856: {
857: PC pc = user->pcmg; /* The MG PC */
858: DM fdm = NULL, cdm = NULL;
859: KSP fksp, cksp;
860: Vec fu, cu = NULL;
861: PetscInt levels, l;
863: KSPBuildSolution(ksp, NULL, &du);
864: PetscObjectComposedDataGetInt((PetscObject) ksp, PetscMGLevelId, level, hasLevel);
865: PCMGGetLevels(pc, &levels);
866: PCMGGetSmoother(pc, levels-1, &fksp);
867: KSPBuildSolution(fksp, NULL, &fu);
868: for (l = levels-1; l > level; --l) {
869: Mat R;
870: Vec s;
872: PCMGGetSmoother(pc, l-1, &cksp);
873: KSPGetDM(cksp, &cdm);
874: DMGetGlobalVector(cdm, &cu);
875: PCMGGetRestriction(pc, l, &R);
876: PCMGGetRScale(pc, l, &s);
877: MatRestrict(R, fu, cu);
878: VecPointwiseMult(cu, cu, s);
879: if (l < levels-1) {DMRestoreGlobalVector(fdm, &fu);}
880: fdm = cdm;
881: fu = cu;
882: }
883: if (levels-1 > level) {
884: VecAXPY(du, 1.0, cu);
885: DMRestoreGlobalVector(cdm, &cu);
886: }
887: }
888: /* Calculate error */
889: DMGetGlobalVector(dm, &r);
890: DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
891: VecAXPY(r,-1.0,du);
892: PetscObjectSetName((PetscObject) r, "solution error");
893: /* View error */
894: #if defined(PETSC_HAVE_HDF5)
895: PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
896: PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
897: VecView(r, viewer);
898: PetscViewerDestroy(&viewer);
899: #endif
900: DMRestoreGlobalVector(dm, &r);
901: return(0);
902: }
904: /*@C
905: SNESMonitorError - Outputs the error at each iteration of an iterative solver.
907: Collective on SNES
909: Input Parameters:
910: + snes - the SNES
911: . its - iteration number
912: . rnorm - 2-norm of residual
913: - ctx - user context
915: Level: intermediate
917: .keywords: SNES, nonlinear, default, monitor, norm
918: .seealso: SNESMonitorDefault(), SNESMonitorSet(), SNESMonitorSolution()
919: @*/
920: static PetscErrorCode SNESMonitorError(SNES snes, PetscInt its, PetscReal rnorm, void *ctx)
921: {
922: AppCtx *user = (AppCtx *) ctx;
923: DM dm;
924: Vec u, r;
925: PetscInt level = -1;
926: PetscBool hasLevel;
927: #if defined(PETSC_HAVE_HDF5)
928: PetscViewer viewer;
929: #endif
930: char buf[256];
934: SNESGetDM(snes, &dm);
935: /* Calculate error */
936: SNESGetSolution(snes, &u);
937: DMGetGlobalVector(dm, &r);
938: PetscObjectSetName((PetscObject) r, "solution error");
939: DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
940: VecAXPY(r, -1.0, u);
941: /* View error */
942: PetscObjectComposedDataGetInt((PetscObject) snes, PetscMGLevelId, level, hasLevel);
943: PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
944: #if defined(PETSC_HAVE_HDF5)
945: PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
946: VecView(r, viewer);
947: PetscViewerDestroy(&viewer);
948: /* Cleanup */
949: DMRestoreGlobalVector(dm, &r);
950: return(0);
951: #else
952: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"You need to configure with --download-hdf5");
953: #endif
954: }
956: int main(int argc, char **argv)
957: {
958: DM dm; /* Problem specification */
959: SNES snes; /* nonlinear solver */
960: Vec u; /* solution vector */
961: Mat A,J; /* Jacobian matrix */
962: MatNullSpace nullSpace; /* May be necessary for Neumann conditions */
963: AppCtx user; /* user-defined work context */
964: JacActionCtx userJ; /* context for Jacobian MF action */
965: PetscReal error = 0.0; /* L_2 error in the solution */
966: PetscBool isFAS;
969: PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
970: ProcessOptions(PETSC_COMM_WORLD, &user);
971: SNESCreate(PETSC_COMM_WORLD, &snes);
972: CreateMesh(PETSC_COMM_WORLD, &user, &dm);
973: SNESSetDM(snes, dm);
974: DMSetApplicationContext(dm, &user);
976: PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
977: SetupDiscretization(dm, &user);
979: DMCreateGlobalVector(dm, &u);
980: PetscObjectSetName((PetscObject) u, "potential");
982: DMCreateMatrix(dm, &J);
983: if (user.jacobianMF) {
984: PetscInt M, m, N, n;
986: MatGetSize(J, &M, &N);
987: MatGetLocalSize(J, &m, &n);
988: MatCreate(PETSC_COMM_WORLD, &A);
989: MatSetSizes(A, m, n, M, N);
990: MatSetType(A, MATSHELL);
991: MatSetUp(A);
992: #if 0
993: MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction);
994: #endif
996: userJ.dm = dm;
997: userJ.J = J;
998: userJ.user = &user;
1000: DMCreateLocalVector(dm, &userJ.u);
1001: if (user.fieldBC) {DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);}
1002: else {DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);}
1003: MatShellSetContext(A, &userJ);
1004: } else {
1005: A = J;
1006: }
1007: if (user.bcType == NEUMANN) {
1008: MatNullSpaceCreate(PetscObjectComm((PetscObject) dm), PETSC_TRUE, 0, NULL, &nullSpace);
1009: MatSetNullSpace(A, nullSpace);
1010: }
1012: DMPlexSetSNESLocalFEM(dm,&user,&user,&user);
1013: SNESSetJacobian(snes, A, J, NULL, NULL);
1015: SNESSetFromOptions(snes);
1017: if (user.fieldBC) {DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u);}
1018: else {DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);}
1019: if (user.restart) {
1020: #if defined(PETSC_HAVE_HDF5)
1021: PetscViewer viewer;
1023: PetscViewerCreate(PETSC_COMM_WORLD, &viewer);
1024: PetscViewerSetType(viewer, PETSCVIEWERHDF5);
1025: PetscViewerFileSetMode(viewer, FILE_MODE_READ);
1026: PetscViewerFileSetName(viewer, user.filename);
1027: PetscViewerHDF5PushGroup(viewer, "/fields");
1028: VecLoad(u, viewer);
1029: PetscViewerHDF5PopGroup(viewer);
1030: PetscViewerDestroy(&viewer);
1031: #endif
1032: }
1033: if (user.showInitial) {
1034: Vec lv;
1035: DMGetLocalVector(dm, &lv);
1036: DMGlobalToLocalBegin(dm, u, INSERT_VALUES, lv);
1037: DMGlobalToLocalEnd(dm, u, INSERT_VALUES, lv);
1038: DMPrintLocalVec(dm, "Local function", 1.0e-10, lv);
1039: DMRestoreLocalVector(dm, &lv);
1040: }
1041: if (user.viewHierarchy) {
1042: SNES lsnes;
1043: KSP ksp;
1044: PC pc;
1045: PetscInt numLevels, l;
1046: PetscBool isMG;
1048: PetscObjectTypeCompare((PetscObject) snes, SNESFAS, &isFAS);
1049: if (isFAS) {
1050: SNESFASGetLevels(snes, &numLevels);
1051: for (l = 0; l < numLevels; ++l) {
1052: SNESFASGetCycleSNES(snes, l, &lsnes);
1053: SNESMonitorSet(lsnes, SNESMonitorError, &user, NULL);
1054: }
1055: } else {
1056: SNESGetKSP(snes, &ksp);
1057: KSPGetPC(ksp, &pc);
1058: PetscObjectTypeCompare((PetscObject) pc, PCMG, &isMG);
1059: if (isMG) {
1060: user.pcmg = pc;
1061: PCMGGetLevels(pc, &numLevels);
1062: for (l = 0; l < numLevels; ++l) {
1063: PCMGGetSmootherDown(pc, l, &ksp);
1064: KSPMonitorSet(ksp, KSPMonitorError, &user, NULL);
1065: }
1066: }
1067: }
1068: }
1069: if (user.runType == RUN_FULL || user.runType == RUN_EXACT) {
1070: PetscErrorCode (*initialGuess[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {zero};
1072: if (user.nonzInit) initialGuess[0] = ecks;
1073: if (user.runType == RUN_FULL) {
1074: DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
1075: }
1076: if (user.debug) {
1077: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1078: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1079: }
1080: SNESSolve(snes, NULL, u);
1081: SNESGetSolution(snes, &u);
1082: SNESGetDM(snes, &dm);
1084: if (user.showSolution) {
1085: PetscPrintf(PETSC_COMM_WORLD, "Solution\n");
1086: VecChop(u, 3.0e-9);
1087: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1088: }
1089: VecViewFromOptions(u, NULL, "-vec_view");
1090: } else if (user.runType == RUN_PERF) {
1091: Vec r;
1092: PetscReal res = 0.0;
1094: SNESGetFunction(snes, &r, NULL, NULL);
1095: SNESComputeFunction(snes, u, r);
1096: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1097: VecChop(r, 1.0e-10);
1098: VecNorm(r, NORM_2, &res);
1099: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1100: } else {
1101: Vec r;
1102: PetscReal res = 0.0, tol = 1.0e-11;
1104: /* Check discretization error */
1105: SNESGetFunction(snes, &r, NULL, NULL);
1106: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1107: if (!user.quiet) {VecView(u, PETSC_VIEWER_STDOUT_WORLD);}
1108: DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
1109: if (error < tol) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < %2.1e\n", tol);}
1110: else {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", error);}
1111: /* Check residual */
1112: SNESComputeFunction(snes, u, r);
1113: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1114: VecChop(r, 1.0e-10);
1115: if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1116: VecNorm(r, NORM_2, &res);
1117: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1118: /* Check Jacobian */
1119: {
1120: Vec b;
1122: SNESComputeJacobian(snes, u, A, A);
1123: VecDuplicate(u, &b);
1124: VecSet(r, 0.0);
1125: SNESComputeFunction(snes, r, b);
1126: MatMult(A, u, r);
1127: VecAXPY(r, 1.0, b);
1128: VecDestroy(&b);
1129: PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n");
1130: VecChop(r, 1.0e-10);
1131: if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1132: VecNorm(r, NORM_2, &res);
1133: PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", res);
1134: }
1135: }
1136: VecViewFromOptions(u, NULL, "-vec_view");
1138: if (user.bdIntegral) {
1139: DMLabel label;
1140: PetscInt id = 1;
1141: PetscScalar bdInt = 0.0;
1142: PetscReal exact = 3.3333333333;
1144: DMGetLabel(dm, "marker", &label);
1145: DMPlexComputeBdIntegral(dm, u, label, 1, &id, bd_integral_2d, &bdInt, NULL);
1146: PetscPrintf(PETSC_COMM_WORLD, "Solution boundary integral: %.4g\n", (double) PetscAbsScalar(bdInt));
1147: if (PetscAbsReal(PetscAbsScalar(bdInt) - exact) > PETSC_SQRT_MACHINE_EPSILON) SETERRQ2(PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Invalid boundary integral %g != %g", bdInt, exact);
1148: }
1150: if (user.bcType == NEUMANN) {MatNullSpaceDestroy(&nullSpace);}
1151: if (user.jacobianMF) {VecDestroy(&userJ.u);}
1152: if (A != J) {MatDestroy(&A);}
1153: MatDestroy(&J);
1154: VecDestroy(&u);
1155: SNESDestroy(&snes);
1156: DMDestroy(&dm);
1157: PetscFree2(user.exactFuncs, user.exactFields);
1158: PetscFinalize();
1159: return ierr;
1160: }
1162: /*TEST
1163: # 2D serial P1 test 0-4
1164: test:
1165: suffix: 0
1166: requires: triangle
1167: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1169: test:
1170: suffix: 1
1171: requires: triangle
1172: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1174: test:
1175: suffix: 2
1176: requires: triangle
1177: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1179: test:
1180: suffix: 3
1181: requires: triangle
1182: args: -run_type test -refinement_limit 0.0 -bc_type neumann -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
1184: test:
1185: suffix: 4
1186: requires: triangle
1187: args: -run_type test -refinement_limit 0.0625 -bc_type neumann -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1189: # 2D serial P2 test 5-8
1190: test:
1191: suffix: 5
1192: requires: triangle
1193: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1195: test:
1196: suffix: 6
1197: requires: triangle
1198: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1200: test:
1201: suffix: 7
1202: requires: triangle
1203: args: -run_type test -refinement_limit 0.0 -bc_type neumann -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
1205: test:
1206: suffix: 8
1207: requires: triangle
1208: args: -run_type test -refinement_limit 0.0625 -bc_type neumann -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
1210: test:
1211: suffix: bd_int_0
1212: requires: triangle
1213: args: -run_type test -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet
1215: test:
1216: suffix: bd_int_1
1217: requires: triangle
1218: args: -run_type test -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet
1220: # 3D serial P1 test 9-12
1221: test:
1222: suffix: 9
1223: requires: ctetgen
1224: args: -run_type test -dim 3 -refinement_limit 0.0 -bc_type dirichlet -interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1
1226: test:
1227: suffix: 10
1228: requires: ctetgen
1229: args: -run_type test -dim 3 -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1
1231: test:
1232: suffix: 11
1233: requires: ctetgen
1234: args: -run_type test -dim 3 -refinement_limit 0.0125 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1
1236: test:
1237: suffix: 12
1238: requires: ctetgen
1239: args: -run_type test -dim 3 -refinement_limit 0.0 -bc_type neumann -interpolate 1 -petscspace_degree 1 -snes_fd -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1
1241: # Analytic variable coefficient 13-20
1242: test:
1243: suffix: 13
1244: requires: triangle
1245: args: -run_type test -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1246: test:
1247: suffix: 14
1248: requires: triangle
1249: args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1250: test:
1251: suffix: 15
1252: requires: triangle
1253: args: -run_type test -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1254: test:
1255: suffix: 16
1256: requires: triangle
1257: args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1258: test:
1259: suffix: 17
1260: requires: hdf5 ctetgen
1261: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1263: test:
1264: suffix: 18
1265: requires: ctetgen
1266: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1268: test:
1269: suffix: 19
1270: requires: ctetgen
1271: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1273: test:
1274: suffix: 20
1275: requires: ctetgen
1276: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1278: # P1 variable coefficient 21-28
1279: test:
1280: suffix: 21
1281: requires: triangle
1282: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1284: test:
1285: suffix: 22
1286: requires: triangle
1287: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1289: test:
1290: suffix: 23
1291: requires: triangle
1292: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1294: test:
1295: suffix: 24
1296: requires: triangle
1297: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1299: test:
1300: suffix: 25
1301: requires: ctetgen
1302: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1304: test:
1305: suffix: 26
1306: requires: ctetgen
1307: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1309: test:
1310: suffix: 27
1311: requires: ctetgen
1312: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1314: test:
1315: suffix: 28
1316: requires: ctetgen
1317: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1319: # P0 variable coefficient 29-36
1320: test:
1321: suffix: 29
1322: requires: triangle
1323: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1325: test:
1326: suffix: 30
1327: requires: triangle
1328: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1330: test:
1331: suffix: 31
1332: requires: triangle
1333: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1335: test:
1336: requires: triangle
1337: suffix: 32
1338: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1340: test:
1341: requires: ctetgen
1342: suffix: 33
1343: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1345: test:
1346: suffix: 34
1347: requires: ctetgen
1348: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1350: test:
1351: suffix: 35
1352: requires: ctetgen
1353: args: -run_type test -dim 3 -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1355: test:
1356: suffix: 36
1357: requires: ctetgen
1358: args: -run_type test -dim 3 -refinement_limit 0.0125 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1360: # Full solve 39-44
1361: test:
1362: suffix: 39
1363: requires: triangle !single
1364: args: -run_type full -refinement_limit 0.015625 -interpolate 1 -petscspace_degree 2 -pc_type gamg -ksp_rtol 1.0e-10 -ksp_monitor_short -ksp_converged_reason -snes_monitor_short -snes_converged_reason ::ascii_info_detail
1365: test:
1366: suffix: 40
1367: requires: triangle !single
1368: args: -run_type full -refinement_limit 0.015625 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 2 -pc_type svd -ksp_rtol 1.0e-10 -snes_monitor_short -snes_converged_reason ::ascii_info_detail
1369: test:
1370: suffix: 41
1371: requires: triangle !single
1372: args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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: 42
1375: requires: triangle !single
1376: args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1377: test:
1378: suffix: 43
1379: requires: triangle !single
1380: nsize: 2
1381: args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1383: test:
1384: suffix: 44
1385: requires: triangle !single
1386: nsize: 2
1387: args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1389: # These tests use a loose tolerance just to exercise the PtAP operations for MATIS and multiple BDDC setup calls inside PCMG
1390: testset:
1391: requires: triangle !single
1392: nsize: 3
1393: args: -interpolate -run_type full -petscspace_degree 1 -dm_mat_type is -pc_type mg -pc_mg_levels 2 -mg_coarse_pc_type bddc -pc_mg_galerkin pmat -ksp_rtol 1.0e-2 -snes_converged_reason -dm_refine_hierarchy 2 -snes_max_it 4
1394: test:
1395: suffix: gmg_bddc
1396: args: -mg_levels_pc_type jacobi
1397: test:
1398: filter: sed -e "s/iterations 1/iterations 4/g"
1399: suffix: gmg_bddc_lev
1400: args: -mg_levels_pc_type bddc
1402: # Restarting
1403: testset:
1404: suffix: restart
1405: requires: hdf5 triangle !complex
1406: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1
1407: test:
1408: args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1409: test:
1410: args: -f sol.h5 -restart
1412: # Periodicity
1413: test:
1414: suffix: periodic_0
1415: requires: triangle
1416: args: -run_type full -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -snes_converged_reason ::ascii_info_detail
1418: test:
1419: requires: !complex
1420: suffix: periodic_1
1421: args: -quiet -run_type test -simplex 0 -x_periodicity periodic -y_periodicity periodic -vec_view vtk:test.vtu:vtk_vtu -interpolate 1 -petscspace_degree 1 -dm_refine 1
1423: # 2D serial P1 test with field bc
1424: test:
1425: suffix: field_bc_p1_0
1426: requires: triangle
1427: args: -run_type test -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1429: test:
1430: suffix: field_bc_p1_1
1431: requires: triangle
1432: args: -run_type test -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1434: test:
1435: suffix: field_bc_p1_2
1436: requires: triangle
1437: args: -run_type test -interpolate 1 -bc_type neumann -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1439: test:
1440: suffix: field_bc_p1_3
1441: requires: triangle
1442: args: -run_type test -dm_refine 1 -interpolate 1 -bc_type neumann -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1444: # 3D serial P1 test with field bc
1445: test:
1446: suffix: field_bc_p1_4
1447: requires: ctetgen
1448: args: -run_type test -dim 3 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1450: test:
1451: suffix: field_bc_p1_5
1452: requires: ctetgen
1453: args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1455: test:
1456: suffix: field_bc_p1_6
1457: requires: ctetgen
1458: args: -run_type test -dim 3 -interpolate 1 -bc_type neumann -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1460: test:
1461: suffix: field_bc_p1_7
1462: requires: ctetgen
1463: args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type neumann -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1465: # 2D serial P2 test with field bc
1466: test:
1467: suffix: field_bc_p2_0
1468: requires: triangle
1469: args: -run_type test -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1471: test:
1472: suffix: field_bc_p2_1
1473: requires: triangle
1474: args: -run_type test -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1476: test:
1477: suffix: field_bc_p2_2
1478: requires: triangle
1479: args: -run_type test -interpolate 1 -bc_type neumann -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1481: test:
1482: suffix: field_bc_p2_3
1483: requires: triangle
1484: args: -run_type test -dm_refine 1 -interpolate 1 -bc_type neumann -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1486: # 3D serial P2 test with field bc
1487: test:
1488: suffix: field_bc_p2_4
1489: requires: ctetgen
1490: args: -run_type test -dim 3 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1492: test:
1493: suffix: field_bc_p2_5
1494: requires: ctetgen
1495: args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1497: test:
1498: suffix: field_bc_p2_6
1499: requires: ctetgen
1500: args: -run_type test -dim 3 -interpolate 1 -bc_type neumann -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1502: test:
1503: suffix: field_bc_p2_7
1504: requires: ctetgen
1505: args: -run_type test -dim 3 -dm_refine 1 -interpolate 1 -bc_type neumann -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -cells 1,1,1
1507: # Full solve simplex: Convergence
1508: test:
1509: suffix: tet_conv_p1_r0
1510: requires: ctetgen
1511: args: -run_type full -dim 3 -dm_refine 0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1
1512: test:
1513: suffix: tet_conv_p1_r2
1514: requires: ctetgen
1515: args: -run_type full -dim 3 -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1
1516: test:
1517: suffix: tet_conv_p1_r3
1518: requires: ctetgen
1519: args: -run_type full -dim 3 -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1
1520: test:
1521: suffix: tet_conv_p2_r0
1522: requires: ctetgen
1523: args: -run_type full -dim 3 -dm_refine 0 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1
1524: test:
1525: suffix: tet_conv_p2_r2
1526: requires: ctetgen
1527: args: -run_type full -dim 3 -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -dm_view -snes_converged_reason ::ascii_info_detail -pc_type lu -cells 1,1,1
1529: # Full solve simplex: BDDC
1530: test:
1531: suffix: tri_bddc
1532: requires: triangle !single
1533: nsize: 5
1534: args: -run_type full -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0
1536: # Full solve simplex: BDDC
1537: test:
1538: suffix: tri_bddc_parmetis
1539: requires: hdf5 triangle !single parmetis
1540: nsize: 4
1541: args: -run_type full -petscpartitioner_type parmetis -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0
1543: testset:
1544: args: -run_type full -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -dm_mat_type is -pc_type bddc -ksp_type gmres -snes_monitor_short -ksp_monitor_short -snes_view -simplex 0 -petscspace_poly_tensor -pc_bddc_corner_selection -cells 3,3 -ksp_rtol 1.e-9 -pc_bddc_use_edges 0
1545: nsize: 5
1546: output_file: output/ex12_quad_bddc.out
1547: filter: sed -e "s/aijcusparse/aij/g" -e "s/aijviennacl/aij/g" -e "s/factorization: cusparse/factorization: petsc/g"
1548: test:
1549: requires: !single
1550: suffix: quad_bddc
1551: test:
1552: requires: !single veccuda
1553: suffix: quad_bddc_cuda
1554: args: -matis_localmat_type aijcusparse -pc_bddc_dirichlet_pc_factor_mat_solver_type cusparse -pc_bddc_neumann_pc_factor_mat_solver_type cusparse
1555: test:
1556: requires: !single viennacl
1557: suffix: quad_bddc_viennacl
1558: args: -matis_localmat_type aijviennacl
1560: # Full solve simplex: ASM
1561: test:
1562: suffix: tri_q2q1_asm_lu
1563: requires: triangle !single
1564: args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0
1566: test:
1567: suffix: tri_q2q1_msm_lu
1568: requires: triangle !single
1569: args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0
1571: test:
1572: suffix: tri_q2q1_asm_sor
1573: requires: triangle !single
1574: args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0
1576: test:
1577: suffix: tri_q2q1_msm_sor
1578: requires: triangle !single
1579: args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0
1581: # Full solve simplex: FAS
1582: test:
1583: suffix: fas_newton_0
1584: requires: triangle !single
1585: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1587: test:
1588: suffix: fas_newton_1
1589: requires: triangle !single
1590: args: -run_type full -dm_refine_hierarchy 3 -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1592: test:
1593: suffix: fas_ngs_0
1594: requires: triangle !single
1595: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -dm_refine_hierarchy 1 -snes_view -fas_levels_1_snes_type ngs -fas_levels_1_snes_monitor_short
1597: test:
1598: suffix: fas_newton_coarse_0
1599: requires: pragmatic triangle
1600: TODO: broken
1601: args: -run_type full -dm_refine 2 -dm_plex_hash_location -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1603: test:
1604: suffix: mg_newton_coarse_0
1605: requires: triangle pragmatic
1606: args: -run_type full -dm_refine 3 -interpolate 1 -petscspace_degree 1 -snes_monitor_short -ksp_monitor_true_residual -snes_linesearch_type basic -snes_converged_reason ::ascii_info_detail -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
1608: test:
1609: suffix: mg_newton_coarse_1
1610: requires: triangle pragmatic
1611: args: -run_type full -dm_refine 5 -interpolate 1 -petscspace_degree 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_converged_reason ::ascii_info_detail -snes_monitor -ksp_monitor_true_residual -mg_levels_ksp_monitor_true_residual -dm_view -ksp_view
1613: test:
1614: suffix: mg_newton_coarse_2
1615: requires: triangle pragmatic
1616: args: -run_type full -dm_refine 5 -interpolate 1 -petscspace_degree 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_converged_reason ::ascii_info_detail -snes_monitor -ksp_monitor_true_residual -mg_levels_ksp_monitor_true_residual -dm_view -ksp_view
1618: # Full solve tensor
1619: test:
1620: suffix: tensor_plex_2d
1621: args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_refine_hierarchy 2 -cells 2,2
1623: test:
1624: suffix: tensor_p4est_2d
1625: requires: p4est
1626: args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_forest_initial_refinement 2 -dm_forest_minimum_refinement 0 -dm_plex_convert_type p4est -cells 2,2
1628: test:
1629: suffix: tensor_plex_3d
1630: args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dim 3 -dm_refine_hierarchy 1 -cells 2,2,2
1632: test:
1633: suffix: tensor_p4est_3d
1634: requires: p4est
1635: args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_forest_initial_refinement 1 -dm_forest_minimum_refinement 0 -dim 3 -dm_plex_convert_type p8est -cells 2,2,2
1637: test:
1638: suffix: p4est_test_q2_conformal_serial
1639: requires: p4est
1640: args: -run_type test -interpolate 1 -petscspace_degree 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -cells 2,2
1642: test:
1643: suffix: p4est_test_q2_conformal_parallel
1644: requires: p4est
1645: nsize: 7
1646: args: -run_type test -interpolate 1 -petscspace_degree 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type simple -cells 2,2
1648: test:
1649: suffix: p4est_test_q2_conformal_parallel_parmetis
1650: requires: hdf5 p4est
1651: nsize: 4
1652: args: -run_type test -interpolate 1 -petscspace_degree 2 -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type parmetis -cells 2,2
1654: test:
1655: suffix: p4est_test_q2_nonconformal_serial
1656: requires: p4est
1657: filter: grep -v "CG or CGNE: variant"
1658: args: -run_type test -interpolate 1 -petscspace_degree 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
1660: test:
1661: suffix: p4est_test_q2_nonconformal_parallel
1662: requires: p4est
1663: filter: grep -v "CG or CGNE: variant"
1664: nsize: 7
1665: args: -run_type test -interpolate 1 -petscspace_degree 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
1667: test:
1668: suffix: p4est_test_q2_nonconformal_parallel_parmetis
1669: requires: hdf5 p4est
1670: nsize: 4
1671: args: -run_type test -interpolate 1 -petscspace_degree 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 parmetis -cells 2,2
1673: test:
1674: suffix: p4est_exact_q2_conformal_serial
1675: requires: p4est !single
1676: args: -run_type exact -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1678: test:
1679: suffix: p4est_exact_q2_conformal_parallel
1680: TODO: broken
1681: requires: p4est !single
1682: nsize: 7
1683: args: -run_type exact -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1685: test:
1686: suffix: p4est_exact_q2_conformal_parallel_parmetis
1687: requires: hdf5 p4est !single
1688: nsize: 4
1689: args: -run_type exact -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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 -petscpartitioner_type parmetis -cells 2,2
1691: test:
1692: suffix: p4est_exact_q2_nonconformal_serial
1693: requires: p4est
1694: args: -run_type exact -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1696: test:
1697: suffix: p4est_exact_q2_nonconformal_parallel
1698: requires: p4est
1699: nsize: 7
1700: args: -run_type exact -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1702: test:
1703: suffix: p4est_exact_q2_nonconformal_parallel_parmetis
1704: requires: hdf5 p4est
1705: nsize: 4
1706: args: -run_type exact -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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 parmetis -cells 2,2
1708: test:
1709: suffix: p4est_full_q2_nonconformal_serial
1710: requires: p4est !single
1711: filter: grep -v "variant HERMITIAN"
1712: args: -run_type full -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1714: test:
1715: suffix: p4est_full_q2_nonconformal_parallel
1716: requires: p4est !single
1717: filter: grep -v "variant HERMITIAN"
1718: nsize: 7
1719: args: -run_type full -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1721: test:
1722: suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1723: requires: p4est
1724: filter: grep -v "variant HERMITIAN"
1725: nsize: 7
1726: args: -run_type full -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1728: test:
1729: suffix: p4est_full_q2_nonconformal_parallel_bddc
1730: requires: p4est
1731: filter: grep -v "variant HERMITIAN"
1732: nsize: 7
1733: args: -run_type full -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1735: test:
1736: suffix: p4est_fas_q2_conformal_serial
1737: requires: p4est
1738: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1739: TODO: broken
1741: test:
1742: suffix: p4est_fas_q2_nonconformal_serial
1743: requires: p4est broken
1744: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1746: test:
1747: suffix: fas_newton_0_p4est
1748: requires: p4est !single
1749: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 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 ::ascii_info_detail -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
1751: # Full solve simplicial AMR
1752: test:
1753: suffix: tri_p1_adapt_0
1754: requires: pragmatic
1755: args: -run_type exact -dim 2 -dm_refine 5 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -variable_coefficient circle -snes_converged_reason ::ascii_info_detail -pc_type lu -adaptor_refinement_factor 1.0 -dm_view -dm_adapt_view -snes_adapt_initial
1757: test:
1758: suffix: tri_p1_adapt_1
1759: requires: pragmatic
1760: args: -run_type exact -dim 2 -dm_refine 5 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -variable_coefficient circle -snes_converged_reason ::ascii_info_detail -pc_type lu -adaptor_refinement_factor 1.0 -dm_view -dm_adapt_iter_view -dm_adapt_view -snes_adapt_sequence 2
1762: test:
1763: suffix: tri_p1_adapt_analytic_0
1764: requires: pragmatic
1765: args: -run_type exact -dim 2 -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -adaptor_monitor -dm_view -dm_adapt_iter_view
1767: # Full solve tensor AMR
1768: test:
1769: suffix: quad_q1_adapt_0
1770: requires: p4est
1771: args: -run_type exact -dim 2 -simplex 0 -dm_plex_convert_type p4est -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -variable_coefficient circle -snes_converged_reason ::ascii_info_detail -pc_type lu -dm_forest_initial_refinement 4 -snes_adapt_initial -dm_view
1772: filter: grep -v DM_
1774: test:
1775: suffix: amr_0
1776: nsize: 5
1777: args: -run_type test -petscpartitioner_type simple -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_refine 1 -cells 2,2
1779: test:
1780: suffix: amr_1
1781: requires: p4est !complex
1782: args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 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
1784: test:
1785: suffix: p4est_solve_bddc
1786: requires: p4est
1787: args: -run_type full -variable_coefficient nonlinear -nonzero_initial_guess 1 -interpolate 1 -petscspace_degree 2 -snes_max_it 20 -snes_type newtonls -dm_mat_type is -pc_type bddc -ksp_type cg -snes_monitor_short -ksp_monitor -snes_linesearch_type bt -snes_converged_reason -snes_view -simplex 0 -petscspace_poly_tensor -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 -pc_bddc_detect_disconnected
1788: nsize: 4
1790: test:
1791: suffix: p4est_solve_fas
1792: requires: p4est
1793: args: -run_type full -variable_coefficient nonlinear -nonzero_initial_guess 1 -interpolate 1 -petscspace_degree 2 -snes_max_it 10 -snes_type fas -snes_linesearch_type bt -snes_fas_levels 3 -fas_coarse_snes_type newtonls -fas_coarse_snes_linesearch_type basic -fas_coarse_ksp_type cg -fas_coarse_pc_type jacobi -fas_coarse_snes_monitor_short -fas_levels_snes_max_it 4 -fas_levels_snes_type newtonls -fas_levels_snes_linesearch_type bt -fas_levels_ksp_type cg -fas_levels_pc_type jacobi -fas_levels_snes_monitor_short -fas_levels_cycle_snes_linesearch_type bt -snes_monitor_short -snes_converged_reason -snes_view -simplex 0 -petscspace_poly_tensor -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
1794: nsize: 4
1795: TODO: identical machine two runs produce slightly different solver trackers
1797: test:
1798: suffix: p4est_convergence_test_1
1799: requires: p4est
1800: args: -quiet -run_type test -interpolate 1 -petscspace_degree 1 -simplex 0 -petscspace_poly_tensor -dm_plex_convert_type p4est -dm_forest_minimum_refinement 2 -dm_forest_initial_refinement 2 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash
1801: nsize: 4
1803: test:
1804: suffix: p4est_convergence_test_2
1805: requires: p4est
1806: args: -quiet -run_type test -interpolate 1 -petscspace_degree 1 -simplex 0 -petscspace_poly_tensor -dm_plex_convert_type p4est -dm_forest_minimum_refinement 3 -dm_forest_initial_refinement 3 -dm_forest_maximum_refinement 5 -dm_p4est_refine_pattern hash
1808: test:
1809: suffix: p4est_convergence_test_3
1810: requires: p4est
1811: args: -quiet -run_type test -interpolate 1 -petscspace_degree 1 -simplex 0 -petscspace_poly_tensor -dm_plex_convert_type p4est -dm_forest_minimum_refinement 4 -dm_forest_initial_refinement 4 -dm_forest_maximum_refinement 6 -dm_p4est_refine_pattern hash
1813: test:
1814: suffix: p4est_convergence_test_4
1815: requires: p4est
1816: args: -quiet -run_type test -interpolate 1 -petscspace_degree 1 -simplex 0 -petscspace_poly_tensor -dm_plex_convert_type p4est -dm_forest_minimum_refinement 5 -dm_forest_initial_refinement 5 -dm_forest_maximum_refinement 7 -dm_p4est_refine_pattern hash
1817: timeoutfactor: 5
1819: # Serial tests with GLVis visualization
1820: test:
1821: suffix: glvis_2d_tet_p1
1822: args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -f ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh
1823: test:
1824: suffix: glvis_2d_tet_p2
1825: args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -f ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh
1826: test:
1827: suffix: glvis_2d_hex_p1
1828: args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -simplex 0 -dm_refine 1
1829: test:
1830: suffix: glvis_2d_hex_p2
1831: args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -simplex 0 -dm_refine 1
1832: test:
1833: suffix: glvis_2d_hex_p2_p4est
1834: requires: p4est
1835: args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 1 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2 -viewer_glvis_dm_plex_enable_ncmesh
1837: TEST*/