Actual source code: ex12.c
petsc-3.11.0 2019-03-29
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, 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 = PetscSinReal(alpha*xy) * (alpha*PetscAbsReal(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] = PetscSinReal(alpha*xy) * (alpha*PetscAbsReal(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->viewHierarchy = PETSC_FALSE;
440: options->simplex = PETSC_TRUE;
441: options->quiet = PETSC_FALSE;
442: options->nonzInit = PETSC_FALSE;
443: options->bdIntegral = PETSC_FALSE;
445: PetscOptionsBegin(comm, "", "Poisson Problem Options", "DMPLEX");
446: PetscOptionsInt("-debug", "The debugging level", "ex12.c", options->debug, &options->debug, NULL);
447: run = options->runType;
448: PetscOptionsEList("-run_type", "The run type", "ex12.c", runTypes, 4, runTypes[options->runType], &run, NULL);
450: options->runType = (RunType) run;
452: PetscOptionsInt("-dim", "The topological mesh dimension", "ex12.c", options->dim, &options->dim, NULL);
453: bd = options->periodicity[0];
454: PetscOptionsEList("-x_periodicity", "The x-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[0]], &bd, NULL);
455: options->periodicity[0] = (DMBoundaryType) bd;
456: bd = options->periodicity[1];
457: PetscOptionsEList("-y_periodicity", "The y-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[1]], &bd, NULL);
458: options->periodicity[1] = (DMBoundaryType) bd;
459: bd = options->periodicity[2];
460: PetscOptionsEList("-z_periodicity", "The z-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[2]], &bd, NULL);
461: options->periodicity[2] = (DMBoundaryType) bd;
462: n = 3;
463: PetscOptionsIntArray("-cells", "The initial mesh division", "ex12.c", options->cells, &n, NULL);
464: PetscOptionsString("-f", "Mesh filename to read", "ex12.c", options->filename, options->filename, sizeof(options->filename), &flg);
465: PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex12.c", options->interpolate, &options->interpolate, NULL);
466: PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex12.c", options->refinementLimit, &options->refinementLimit, NULL);
467: bc = options->bcType;
468: PetscOptionsEList("-bc_type","Type of boundary condition","ex12.c",bcTypes,3,bcTypes[options->bcType],&bc,NULL);
469: options->bcType = (BCType) bc;
470: coeff = options->variableCoefficient;
471: PetscOptionsEList("-variable_coefficient","Type of variable coefficent","ex12.c",coeffTypes,6,coeffTypes[options->variableCoefficient],&coeff,NULL);
472: options->variableCoefficient = (CoeffType) coeff;
474: PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL);
475: PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL);
476: PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL);
477: PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL);
478: PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL);
479: PetscOptionsBool("-dm_view_hierarchy", "View the coarsened hierarchy", "ex12.c", options->viewHierarchy, &options->viewHierarchy, NULL);
480: PetscOptionsBool("-simplex", "Simplicial (true) or tensor (false) mesh", "ex12.c", options->simplex, &options->simplex, NULL);
481: PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL);
482: PetscOptionsBool("-nonzero_initial_guess", "nonzero intial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL);
483: PetscOptionsBool("-bd_integral", "Compute the integral of the solution on the boundary", "ex12.c", options->bdIntegral, &options->bdIntegral, NULL);
484: PetscOptionsEnd();
485: PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
486: return(0);
487: }
489: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
490: {
491: DMLabel label;
495: DMCreateLabel(dm, name);
496: DMGetLabel(dm, name, &label);
497: DMPlexMarkBoundaryFaces(dm, 1, label);
498: DMPlexLabelComplete(dm, label);
499: return(0);
500: }
502: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
503: {
504: PetscInt dim = user->dim;
505: const char *filename = user->filename;
506: PetscBool interpolate = user->interpolate;
507: PetscReal refinementLimit = user->refinementLimit;
508: size_t len;
512: PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
513: PetscStrlen(filename, &len);
514: if (!len) {
515: PetscInt d;
517: if (user->periodicity[0] || user->periodicity[1] || user->periodicity[2]) for (d = 0; d < dim; ++d) user->cells[d] = PetscMax(user->cells[d], 3);
518: DMPlexCreateBoxMesh(comm, dim, user->simplex, user->cells, NULL, NULL, user->periodicity, interpolate, dm);
519: PetscObjectSetName((PetscObject) *dm, "Mesh");
520: } else {
521: DMPlexCreateFromFile(comm, filename, interpolate, dm);
522: DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
523: }
524: {
525: PetscPartitioner part;
526: DM refinedMesh = NULL;
527: DM distributedMesh = NULL;
529: /* Refine mesh using a volume constraint */
530: if (refinementLimit > 0.0) {
531: DMPlexSetRefinementLimit(*dm, refinementLimit);
532: DMRefine(*dm, comm, &refinedMesh);
533: if (refinedMesh) {
534: const char *name;
536: PetscObjectGetName((PetscObject) *dm, &name);
537: PetscObjectSetName((PetscObject) refinedMesh, name);
538: DMDestroy(dm);
539: *dm = refinedMesh;
540: }
541: }
542: /* Distribute mesh over processes */
543: DMPlexGetPartitioner(*dm,&part);
544: PetscPartitionerSetFromOptions(part);
545: DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
546: if (distributedMesh) {
547: DMDestroy(dm);
548: *dm = distributedMesh;
549: }
550: }
551: if (user->bcType == NEUMANN) {
552: DMLabel label;
554: DMCreateLabel(*dm, "boundary");
555: DMGetLabel(*dm, "boundary", &label);
556: DMPlexMarkBoundaryFaces(*dm, 1, label);
557: } else if (user->bcType == DIRICHLET) {
558: PetscBool hasLabel;
560: DMHasLabel(*dm,"marker",&hasLabel);
561: if (!hasLabel) {CreateBCLabel(*dm, "marker");}
562: }
563: {
564: char convType[256];
565: PetscBool flg;
567: PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
568: PetscOptionsFList("-dm_plex_convert_type","Convert DMPlex to another format","ex12",DMList,DMPLEX,convType,256,&flg);
569: PetscOptionsEnd();
570: if (flg) {
571: DM dmConv;
573: DMConvert(*dm,convType,&dmConv);
574: if (dmConv) {
575: DMDestroy(dm);
576: *dm = dmConv;
577: }
578: }
579: }
580: DMLocalizeCoordinates(*dm); /* needed for periodic */
581: DMSetFromOptions(*dm);
582: DMViewFromOptions(*dm, NULL, "-dm_view");
583: if (user->viewHierarchy) {
584: DM cdm = *dm;
585: PetscInt i = 0;
586: char buf[256];
588: while (cdm) {
589: DMSetUp(cdm);
590: DMGetCoarseDM(cdm, &cdm);
591: ++i;
592: }
593: cdm = *dm;
594: while (cdm) {
595: PetscViewer viewer;
596: PetscBool isHDF5, isVTK;
598: --i;
599: PetscViewerCreate(comm,&viewer);
600: PetscViewerSetType(viewer,PETSCVIEWERHDF5);
601: PetscViewerSetOptionsPrefix(viewer,"hierarchy_");
602: PetscViewerSetFromOptions(viewer);
603: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&isHDF5);
604: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERVTK,&isVTK);
605: if (isHDF5) {
606: PetscSNPrintf(buf, 256, "ex12-%d.h5", i);
607: } else if (isVTK) {
608: PetscSNPrintf(buf, 256, "ex12-%d.vtu", i);
609: PetscViewerPushFormat(viewer,PETSC_VIEWER_VTK_VTU);
610: } else {
611: PetscSNPrintf(buf, 256, "ex12-%d", i);
612: }
613: PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
614: PetscViewerFileSetName(viewer,buf);
615: DMView(cdm, viewer);
616: PetscViewerDestroy(&viewer);
617: DMGetCoarseDM(cdm, &cdm);
618: }
619: }
620: PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
621: return(0);
622: }
624: static PetscErrorCode SetupProblem(DM dm, AppCtx *user)
625: {
626: PetscDS prob;
627: const PetscInt id = 1;
631: DMGetDS(dm, &prob);
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: return(0);
703: }
705: static PetscErrorCode SetupMaterial(DM dm, DM dmAux, AppCtx *user)
706: {
707: PetscErrorCode (*matFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {nu_2d};
708: Vec nu;
712: DMCreateLocalVector(dmAux, &nu);
713: DMProjectFunctionLocal(dmAux, 0.0, matFuncs, NULL, INSERT_ALL_VALUES, nu);
714: PetscObjectCompose((PetscObject) dm, "A", (PetscObject) nu);
715: VecDestroy(&nu);
716: return(0);
717: }
719: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
720: {
721: PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
722: Vec uexact;
723: PetscInt dim;
727: DMGetDimension(dm, &dim);
728: if (dim == 2) bcFuncs[0] = quadratic_u_2d;
729: else bcFuncs[0] = quadratic_u_3d;
730: DMCreateLocalVector(dmAux, &uexact);
731: DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
732: PetscObjectCompose((PetscObject) dm, "A", (PetscObject) uexact);
733: VecDestroy(&uexact);
734: return(0);
735: }
737: static PetscErrorCode SetupAuxDM(DM dm, PetscFE feAux, AppCtx *user)
738: {
739: DM dmAux, coordDM;
743: /* MUST call DMGetCoordinateDM() in order to get p4est setup if present */
744: DMGetCoordinateDM(dm, &coordDM);
745: if (!feAux) return(0);
746: DMClone(dm, &dmAux);
747: PetscObjectCompose((PetscObject) dm, "dmAux", (PetscObject) dmAux);
748: DMSetCoordinateDM(dmAux, coordDM);
749: DMSetField(dmAux, 0, NULL, (PetscObject) feAux);
750: DMCreateDS(dmAux);
751: if (user->fieldBC) {SetupBC(dm, dmAux, user);}
752: else {SetupMaterial(dm, dmAux, user);}
753: DMDestroy(&dmAux);
754: return(0);
755: }
757: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
758: {
759: DM cdm = dm;
760: const PetscInt dim = user->dim;
761: PetscFE fe, feAux = NULL;
762: PetscBool simplex = user->simplex;
763: MPI_Comm comm;
767: /* Create finite element for each field and auxiliary field */
768: PetscObjectGetComm((PetscObject) dm, &comm);
769: PetscFECreateDefault(comm, dim, 1, simplex, NULL, -1, &fe);
770: PetscObjectSetName((PetscObject) fe, "potential");
771: if (user->variableCoefficient == COEFF_FIELD) {
772: PetscQuadrature q;
774: PetscFECreateDefault(comm, dim, 1, simplex, "mat_", -1, &feAux);
775: PetscFEGetQuadrature(fe, &q);
776: PetscFESetQuadrature(feAux, q);
777: } else if (user->fieldBC) {
778: PetscQuadrature q;
780: PetscFECreateDefault(comm, dim, 1, simplex, "bc_", -1, &feAux);
781: PetscFEGetQuadrature(fe, &q);
782: PetscFESetQuadrature(feAux, q);
783: }
784: /* Set discretization and boundary conditions for each mesh */
785: DMSetField(dm, 0, NULL, (PetscObject) fe);
786: DMCreateDS(dm);
787: SetupProblem(dm, user);
788: while (cdm) {
789: DMCopyDisc(dm, cdm);
790: SetupAuxDM(cdm, feAux, user);
791: if (user->bcType == DIRICHLET) {
792: PetscBool hasLabel;
794: DMHasLabel(cdm, "marker", &hasLabel);
795: if (!hasLabel) {CreateBCLabel(cdm, "marker");}
796: }
797: DMGetCoarseDM(cdm, &cdm);
798: }
799: PetscFEDestroy(&fe);
800: PetscFEDestroy(&feAux);
801: return(0);
802: }
804: #include "petsc/private/petscimpl.h"
806: /*@C
807: KSPMonitorError - Outputs the error at each iteration of an iterative solver.
809: Collective on KSP
811: Input Parameters:
812: + ksp - the KSP
813: . its - iteration number
814: . rnorm - 2-norm, preconditioned residual value (may be estimated).
815: - ctx - monitor context
817: Level: intermediate
819: .keywords: KSP, default, monitor, residual
820: .seealso: KSPMonitorSet(), KSPMonitorTrueResidualNorm(), KSPMonitorDefault()
821: @*/
822: static PetscErrorCode KSPMonitorError(KSP ksp, PetscInt its, PetscReal rnorm, void *ctx)
823: {
824: AppCtx *user = (AppCtx *) ctx;
825: DM dm;
826: Vec du = NULL, r;
827: PetscInt level = 0;
828: PetscBool hasLevel;
829: #if defined(PETSC_HAVE_HDF5)
830: PetscViewer viewer;
831: char buf[256];
832: #endif
836: KSPGetDM(ksp, &dm);
837: /* Calculate solution */
838: {
839: PC pc = user->pcmg; /* The MG PC */
840: DM fdm = NULL, cdm = NULL;
841: KSP fksp, cksp;
842: Vec fu, cu = NULL;
843: PetscInt levels, l;
845: KSPBuildSolution(ksp, NULL, &du);
846: PetscObjectComposedDataGetInt((PetscObject) ksp, PetscMGLevelId, level, hasLevel);
847: PCMGGetLevels(pc, &levels);
848: PCMGGetSmoother(pc, levels-1, &fksp);
849: KSPBuildSolution(fksp, NULL, &fu);
850: for (l = levels-1; l > level; --l) {
851: Mat R;
852: Vec s;
854: PCMGGetSmoother(pc, l-1, &cksp);
855: KSPGetDM(cksp, &cdm);
856: DMGetGlobalVector(cdm, &cu);
857: PCMGGetRestriction(pc, l, &R);
858: PCMGGetRScale(pc, l, &s);
859: MatRestrict(R, fu, cu);
860: VecPointwiseMult(cu, cu, s);
861: if (l < levels-1) {DMRestoreGlobalVector(fdm, &fu);}
862: fdm = cdm;
863: fu = cu;
864: }
865: if (levels-1 > level) {
866: VecAXPY(du, 1.0, cu);
867: DMRestoreGlobalVector(cdm, &cu);
868: }
869: }
870: /* Calculate error */
871: DMGetGlobalVector(dm, &r);
872: DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
873: VecAXPY(r,-1.0,du);
874: PetscObjectSetName((PetscObject) r, "solution error");
875: /* View error */
876: #if defined(PETSC_HAVE_HDF5)
877: PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
878: PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
879: VecView(r, viewer);
880: PetscViewerDestroy(&viewer);
881: #endif
882: DMRestoreGlobalVector(dm, &r);
883: return(0);
884: }
886: /*@C
887: SNESMonitorError - Outputs the error at each iteration of an iterative solver.
889: Collective on SNES
891: Input Parameters:
892: + snes - the SNES
893: . its - iteration number
894: . rnorm - 2-norm of residual
895: - ctx - user context
897: Level: intermediate
899: .keywords: SNES, nonlinear, default, monitor, norm
900: .seealso: SNESMonitorDefault(), SNESMonitorSet(), SNESMonitorSolution()
901: @*/
902: static PetscErrorCode SNESMonitorError(SNES snes, PetscInt its, PetscReal rnorm, void *ctx)
903: {
904: AppCtx *user = (AppCtx *) ctx;
905: DM dm;
906: Vec u, r;
907: PetscInt level = -1;
908: PetscBool hasLevel;
909: #if defined(PETSC_HAVE_HDF5)
910: PetscViewer viewer;
911: #endif
912: char buf[256];
916: SNESGetDM(snes, &dm);
917: /* Calculate error */
918: SNESGetSolution(snes, &u);
919: DMGetGlobalVector(dm, &r);
920: PetscObjectSetName((PetscObject) r, "solution error");
921: DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
922: VecAXPY(r, -1.0, u);
923: /* View error */
924: PetscObjectComposedDataGetInt((PetscObject) snes, PetscMGLevelId, level, hasLevel);
925: PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
926: #if defined(PETSC_HAVE_HDF5)
927: PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
928: VecView(r, viewer);
929: PetscViewerDestroy(&viewer);
930: /* Cleanup */
931: DMRestoreGlobalVector(dm, &r);
932: return(0);
933: #else
934: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"You need to configure with --download-hdf5");
935: #endif
936: }
938: int main(int argc, char **argv)
939: {
940: DM dm; /* Problem specification */
941: SNES snes; /* nonlinear solver */
942: Vec u; /* solution vector */
943: Mat A,J; /* Jacobian matrix */
944: MatNullSpace nullSpace; /* May be necessary for Neumann conditions */
945: AppCtx user; /* user-defined work context */
946: JacActionCtx userJ; /* context for Jacobian MF action */
947: PetscReal error = 0.0; /* L_2 error in the solution */
948: PetscBool isFAS;
951: PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
952: ProcessOptions(PETSC_COMM_WORLD, &user);
953: SNESCreate(PETSC_COMM_WORLD, &snes);
954: CreateMesh(PETSC_COMM_WORLD, &user, &dm);
955: SNESSetDM(snes, dm);
956: DMSetApplicationContext(dm, &user);
958: PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
959: SetupDiscretization(dm, &user);
961: DMCreateGlobalVector(dm, &u);
962: PetscObjectSetName((PetscObject) u, "potential");
964: DMCreateMatrix(dm, &J);
965: if (user.jacobianMF) {
966: PetscInt M, m, N, n;
968: MatGetSize(J, &M, &N);
969: MatGetLocalSize(J, &m, &n);
970: MatCreate(PETSC_COMM_WORLD, &A);
971: MatSetSizes(A, m, n, M, N);
972: MatSetType(A, MATSHELL);
973: MatSetUp(A);
974: #if 0
975: MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction);
976: #endif
978: userJ.dm = dm;
979: userJ.J = J;
980: userJ.user = &user;
982: DMCreateLocalVector(dm, &userJ.u);
983: if (user.fieldBC) {DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);}
984: else {DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);}
985: MatShellSetContext(A, &userJ);
986: } else {
987: A = J;
988: }
989: if (user.bcType == NEUMANN) {
990: MatNullSpaceCreate(PetscObjectComm((PetscObject) dm), PETSC_TRUE, 0, NULL, &nullSpace);
991: MatSetNullSpace(A, nullSpace);
992: }
994: DMPlexSetSNESLocalFEM(dm,&user,&user,&user);
995: SNESSetJacobian(snes, A, J, NULL, NULL);
997: SNESSetFromOptions(snes);
999: if (user.fieldBC) {DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u);}
1000: else {DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);}
1001: if (user.restart) {
1002: #if defined(PETSC_HAVE_HDF5)
1003: PetscViewer viewer;
1005: PetscViewerCreate(PETSC_COMM_WORLD, &viewer);
1006: PetscViewerSetType(viewer, PETSCVIEWERHDF5);
1007: PetscViewerFileSetMode(viewer, FILE_MODE_READ);
1008: PetscViewerFileSetName(viewer, user.filename);
1009: PetscViewerHDF5PushGroup(viewer, "/fields");
1010: VecLoad(u, viewer);
1011: PetscViewerHDF5PopGroup(viewer);
1012: PetscViewerDestroy(&viewer);
1013: #endif
1014: }
1015: if (user.showInitial) {
1016: Vec lv;
1017: DMGetLocalVector(dm, &lv);
1018: DMGlobalToLocalBegin(dm, u, INSERT_VALUES, lv);
1019: DMGlobalToLocalEnd(dm, u, INSERT_VALUES, lv);
1020: DMPrintLocalVec(dm, "Local function", 1.0e-10, lv);
1021: DMRestoreLocalVector(dm, &lv);
1022: }
1023: if (user.viewHierarchy) {
1024: SNES lsnes;
1025: KSP ksp;
1026: PC pc;
1027: PetscInt numLevels, l;
1028: PetscBool isMG;
1030: PetscObjectTypeCompare((PetscObject) snes, SNESFAS, &isFAS);
1031: if (isFAS) {
1032: SNESFASGetLevels(snes, &numLevels);
1033: for (l = 0; l < numLevels; ++l) {
1034: SNESFASGetCycleSNES(snes, l, &lsnes);
1035: SNESMonitorSet(lsnes, SNESMonitorError, &user, NULL);
1036: }
1037: } else {
1038: SNESGetKSP(snes, &ksp);
1039: KSPGetPC(ksp, &pc);
1040: PetscObjectTypeCompare((PetscObject) pc, PCMG, &isMG);
1041: if (isMG) {
1042: user.pcmg = pc;
1043: PCMGGetLevels(pc, &numLevels);
1044: for (l = 0; l < numLevels; ++l) {
1045: PCMGGetSmootherDown(pc, l, &ksp);
1046: KSPMonitorSet(ksp, KSPMonitorError, &user, NULL);
1047: }
1048: }
1049: }
1050: }
1051: if (user.runType == RUN_FULL || user.runType == RUN_EXACT) {
1052: PetscErrorCode (*initialGuess[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {zero};
1054: if (user.nonzInit) initialGuess[0] = ecks;
1055: if (user.runType == RUN_FULL) {
1056: DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
1057: }
1058: if (user.debug) {
1059: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1060: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1061: }
1062: SNESSolve(snes, NULL, u);
1063: SNESGetSolution(snes, &u);
1064: SNESGetDM(snes, &dm);
1066: if (user.showSolution) {
1067: PetscPrintf(PETSC_COMM_WORLD, "Solution\n");
1068: VecChop(u, 3.0e-9);
1069: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1070: }
1071: VecViewFromOptions(u, NULL, "-vec_view");
1072: } else if (user.runType == RUN_PERF) {
1073: Vec r;
1074: PetscReal res = 0.0;
1076: SNESGetFunction(snes, &r, NULL, NULL);
1077: SNESComputeFunction(snes, u, r);
1078: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1079: VecChop(r, 1.0e-10);
1080: VecNorm(r, NORM_2, &res);
1081: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1082: } else {
1083: Vec r;
1084: PetscReal res = 0.0, tol = 1.0e-11;
1086: /* Check discretization error */
1087: SNESGetFunction(snes, &r, NULL, NULL);
1088: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1089: if (!user.quiet) {VecView(u, PETSC_VIEWER_STDOUT_WORLD);}
1090: DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
1091: if (error < tol) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < %2.1e\n", tol);}
1092: else {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", error);}
1093: /* Check residual */
1094: SNESComputeFunction(snes, u, r);
1095: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1096: VecChop(r, 1.0e-10);
1097: if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1098: VecNorm(r, NORM_2, &res);
1099: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
1100: /* Check Jacobian */
1101: {
1102: Vec b;
1104: SNESComputeJacobian(snes, u, A, A);
1105: VecDuplicate(u, &b);
1106: VecSet(r, 0.0);
1107: SNESComputeFunction(snes, r, b);
1108: MatMult(A, u, r);
1109: VecAXPY(r, 1.0, b);
1110: VecDestroy(&b);
1111: PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n");
1112: VecChop(r, 1.0e-10);
1113: if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1114: VecNorm(r, NORM_2, &res);
1115: PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", res);
1116: }
1117: }
1118: VecViewFromOptions(u, NULL, "-vec_view");
1120: if (user.bdIntegral) {
1121: DMLabel label;
1122: PetscInt id = 1;
1123: PetscScalar bdInt = 0.0;
1124: PetscReal exact = 3.3333333333;
1126: DMGetLabel(dm, "marker", &label);
1127: DMPlexComputeBdIntegral(dm, u, label, 1, &id, bd_integral_2d, &bdInt, NULL);
1128: PetscPrintf(PETSC_COMM_WORLD, "Solution boundary integral: %.4g\n", (double) PetscAbsScalar(bdInt));
1129: if (PetscAbsReal(PetscAbsScalar(bdInt) - exact) > PETSC_SQRT_MACHINE_EPSILON) SETERRQ2(PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Invalid boundary integral %g != %g", bdInt, exact);
1130: }
1132: if (user.bcType == NEUMANN) {MatNullSpaceDestroy(&nullSpace);}
1133: if (user.jacobianMF) {VecDestroy(&userJ.u);}
1134: if (A != J) {MatDestroy(&A);}
1135: MatDestroy(&J);
1136: VecDestroy(&u);
1137: SNESDestroy(&snes);
1138: DMDestroy(&dm);
1139: PetscFree2(user.exactFuncs, user.exactFields);
1140: PetscFinalize();
1141: return ierr;
1142: }
1144: /*TEST
1145: # 2D serial P1 test 0-4
1146: test:
1147: suffix: 2d_p1_0
1148: requires: triangle
1149: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1151: test:
1152: suffix: 2d_p1_1
1153: requires: triangle
1154: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1156: test:
1157: suffix: 2d_p1_2
1158: requires: triangle
1159: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1161: test:
1162: suffix: 2d_p1_neumann_0
1163: requires: triangle
1164: 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
1166: test:
1167: suffix: 2d_p1_neumann_1
1168: requires: triangle
1169: args: -run_type test -refinement_limit 0.0625 -bc_type neumann -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1171: # 2D serial P2 test 5-8
1172: test:
1173: suffix: 2d_p2_0
1174: requires: triangle
1175: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1177: test:
1178: suffix: 2d_p2_1
1179: requires: triangle
1180: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1182: test:
1183: suffix: 2d_p2_neumann_0
1184: requires: triangle
1185: 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
1187: test:
1188: suffix: 2d_p2_neumann_1
1189: requires: triangle
1190: 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
1192: test:
1193: suffix: bd_int_0
1194: requires: triangle
1195: args: -run_type test -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet
1197: test:
1198: suffix: bd_int_1
1199: requires: triangle
1200: args: -run_type test -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet
1202: # 3D serial P1 test 9-12
1203: test:
1204: suffix: 3d_p1_0
1205: requires: ctetgen
1206: 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
1208: test:
1209: suffix: 3d_p1_1
1210: requires: ctetgen
1211: 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
1213: test:
1214: suffix: 3d_p1_2
1215: requires: ctetgen
1216: 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
1218: test:
1219: suffix: 3d_p1_neumann_0
1220: requires: ctetgen
1221: args: -run_type test -dim 3 -bc_type neumann -interpolate 1 -petscspace_degree 1 -snes_fd -show_initial -dm_plex_print_fem 1 -dm_view -cells 1,1,1
1223: # Analytic variable coefficient 13-20
1224: test:
1225: suffix: 13
1226: requires: triangle
1227: args: -run_type test -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1228: test:
1229: suffix: 14
1230: requires: triangle
1231: args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1232: test:
1233: suffix: 15
1234: requires: triangle
1235: args: -run_type test -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1236: test:
1237: suffix: 16
1238: requires: triangle
1239: args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1240: test:
1241: suffix: 17
1242: requires: ctetgen
1243: 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
1245: test:
1246: suffix: 18
1247: requires: ctetgen
1248: 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
1250: test:
1251: suffix: 19
1252: requires: ctetgen
1253: 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
1255: test:
1256: suffix: 20
1257: requires: ctetgen
1258: 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
1260: # P1 variable coefficient 21-28
1261: test:
1262: suffix: 21
1263: requires: triangle
1264: 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
1266: test:
1267: suffix: 22
1268: requires: triangle
1269: 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
1271: test:
1272: suffix: 23
1273: requires: triangle
1274: 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
1276: test:
1277: suffix: 24
1278: requires: triangle
1279: 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
1281: test:
1282: suffix: 25
1283: requires: ctetgen
1284: 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
1286: test:
1287: suffix: 26
1288: requires: ctetgen
1289: 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
1291: test:
1292: suffix: 27
1293: requires: ctetgen
1294: 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
1296: test:
1297: suffix: 28
1298: requires: ctetgen
1299: 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
1301: # P0 variable coefficient 29-36
1302: test:
1303: suffix: 29
1304: requires: triangle
1305: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1307: test:
1308: suffix: 30
1309: requires: triangle
1310: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1312: test:
1313: suffix: 31
1314: requires: triangle
1315: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1317: test:
1318: requires: triangle
1319: suffix: 32
1320: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1322: test:
1323: requires: ctetgen
1324: suffix: 33
1325: 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
1327: test:
1328: suffix: 34
1329: requires: ctetgen
1330: 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
1332: test:
1333: suffix: 35
1334: requires: ctetgen
1335: 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
1337: test:
1338: suffix: 36
1339: requires: ctetgen
1340: 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
1342: # Full solve 39-44
1343: test:
1344: suffix: 39
1345: requires: triangle !single
1346: 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
1347: test:
1348: suffix: 40
1349: requires: triangle !single
1350: 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
1351: test:
1352: suffix: 41
1353: requires: triangle !single
1354: 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
1355: test:
1356: suffix: 42
1357: requires: triangle !single
1358: 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
1359: test:
1360: suffix: 43
1361: requires: triangle !single
1362: nsize: 2
1363: 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
1365: test:
1366: suffix: 44
1367: requires: triangle !single
1368: nsize: 2
1369: 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
1371: # These tests use a loose tolerance just to exercise the PtAP operations for MATIS and multiple BDDC setup calls inside PCMG
1372: testset:
1373: requires: triangle !single
1374: nsize: 3
1375: 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
1376: test:
1377: suffix: gmg_bddc
1378: args: -mg_levels_pc_type jacobi
1379: test:
1380: filter: sed -e "s/iterations [0-4]/iterations 4/g"
1381: suffix: gmg_bddc_lev
1382: args: -mg_levels_pc_type bddc
1384: # Restarting
1385: testset:
1386: suffix: restart
1387: requires: hdf5 triangle !complex
1388: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1
1389: test:
1390: args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1391: test:
1392: args: -f sol.h5 -restart
1394: # Periodicity
1395: test:
1396: suffix: periodic_0
1397: requires: triangle
1398: args: -run_type full -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -snes_converged_reason ::ascii_info_detail
1400: test:
1401: requires: !complex
1402: suffix: periodic_1
1403: 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
1405: # 2D serial P1 test with field bc
1406: test:
1407: suffix: field_bc_2d_p1_0
1408: requires: triangle
1409: 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
1411: test:
1412: suffix: field_bc_2d_p1_1
1413: requires: triangle
1414: 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
1416: test:
1417: suffix: field_bc_2d_p1_neumann_0
1418: requires: triangle
1419: 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
1421: test:
1422: suffix: field_bc_2d_p1_neumann_1
1423: requires: triangle
1424: 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
1426: # 3D serial P1 test with field bc
1427: test:
1428: suffix: field_bc_3d_p1_0
1429: requires: ctetgen
1430: 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
1432: test:
1433: suffix: field_bc_3d_p1_1
1434: requires: ctetgen
1435: 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
1437: test:
1438: suffix: field_bc_3d_p1_neumann_0
1439: requires: ctetgen
1440: 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
1442: test:
1443: suffix: field_bc_3d_p1_neumann_1
1444: requires: ctetgen
1445: 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
1447: # 2D serial P2 test with field bc
1448: test:
1449: suffix: field_bc_2d_p2_0
1450: requires: triangle
1451: 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
1453: test:
1454: suffix: field_bc_2d_p2_1
1455: requires: triangle
1456: 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
1458: test:
1459: suffix: field_bc_2d_p2_neumann_0
1460: requires: triangle
1461: 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
1463: test:
1464: suffix: field_bc_2d_p2_neumann_1
1465: requires: triangle
1466: 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
1468: # 3D serial P2 test with field bc
1469: test:
1470: suffix: field_bc_3d_p2_0
1471: requires: ctetgen
1472: 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
1474: test:
1475: suffix: field_bc_3d_p2_1
1476: requires: ctetgen
1477: 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
1479: test:
1480: suffix: field_bc_3d_p2_neumann_0
1481: requires: ctetgen
1482: 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
1484: test:
1485: suffix: field_bc_3d_p2_neumann_1
1486: requires: ctetgen
1487: 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
1489: # Full solve simplex: Convergence
1490: test:
1491: suffix: tet_conv_p1_r0
1492: requires: ctetgen
1493: 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
1494: test:
1495: suffix: tet_conv_p1_r2
1496: requires: ctetgen
1497: 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
1498: test:
1499: suffix: tet_conv_p1_r3
1500: requires: ctetgen
1501: 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
1502: test:
1503: suffix: tet_conv_p2_r0
1504: requires: ctetgen
1505: 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
1506: test:
1507: suffix: tet_conv_p2_r2
1508: requires: ctetgen
1509: 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
1511: # Full solve simplex: BDDC
1512: test:
1513: suffix: tri_bddc
1514: requires: triangle !single
1515: nsize: 5
1516: 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
1518: # Full solve simplex: BDDC
1519: test:
1520: suffix: tri_bddc_parmetis
1521: requires: triangle !single parmetis
1522: nsize: 4
1523: 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
1525: testset:
1526: 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
1527: nsize: 5
1528: output_file: output/ex12_quad_bddc.out
1529: filter: sed -e "s/aijcusparse/aij/g" -e "s/aijviennacl/aij/g" -e "s/factorization: cusparse/factorization: petsc/g"
1530: test:
1531: requires: !single
1532: suffix: quad_bddc
1533: test:
1534: requires: !single cuda
1535: suffix: quad_bddc_cuda
1536: args: -matis_localmat_type aijcusparse -pc_bddc_dirichlet_pc_factor_mat_solver_type cusparse -pc_bddc_neumann_pc_factor_mat_solver_type cusparse
1537: test:
1538: requires: !single viennacl
1539: suffix: quad_bddc_viennacl
1540: args: -matis_localmat_type aijviennacl
1542: # Full solve simplex: ASM
1543: test:
1544: suffix: tri_q2q1_asm_lu
1545: requires: triangle !single
1546: 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
1548: test:
1549: suffix: tri_q2q1_msm_lu
1550: requires: triangle !single
1551: 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
1553: test:
1554: suffix: tri_q2q1_asm_sor
1555: requires: triangle !single
1556: 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
1558: test:
1559: suffix: tri_q2q1_msm_sor
1560: requires: triangle !single
1561: 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
1563: # Full solve simplex: FAS
1564: test:
1565: suffix: fas_newton_0
1566: requires: triangle !single
1567: 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
1569: test:
1570: suffix: fas_newton_1
1571: requires: triangle !single
1572: 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
1574: test:
1575: suffix: fas_ngs_0
1576: requires: triangle !single
1577: 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
1579: test:
1580: suffix: fas_newton_coarse_0
1581: requires: pragmatic triangle
1582: TODO: broken
1583: 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
1585: test:
1586: suffix: mg_newton_coarse_0
1587: requires: triangle pragmatic
1588: 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
1590: test:
1591: suffix: mg_newton_coarse_1
1592: requires: triangle pragmatic
1593: 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
1595: test:
1596: suffix: mg_newton_coarse_2
1597: requires: triangle pragmatic
1598: 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
1600: # Full solve tensor
1601: test:
1602: suffix: tensor_plex_2d
1603: args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_refine_hierarchy 2 -cells 2,2
1605: test:
1606: suffix: tensor_p4est_2d
1607: requires: p4est
1608: 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
1610: test:
1611: suffix: tensor_plex_3d
1612: 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
1614: test:
1615: suffix: tensor_p4est_3d
1616: requires: p4est
1617: 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
1619: test:
1620: suffix: p4est_test_q2_conformal_serial
1621: requires: p4est
1622: 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
1624: test:
1625: suffix: p4est_test_q2_conformal_parallel
1626: requires: p4est
1627: nsize: 7
1628: 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
1630: test:
1631: suffix: p4est_test_q2_conformal_parallel_parmetis
1632: requires: parmetis p4est
1633: nsize: 4
1634: 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
1636: test:
1637: suffix: p4est_test_q2_nonconformal_serial
1638: requires: p4est
1639: filter: grep -v "CG or CGNE: variant"
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 -dm_forest_maximum_refinement 4 -dm_p4est_refine_pattern hash -cells 2,2
1642: test:
1643: suffix: p4est_test_q2_nonconformal_parallel
1644: requires: p4est
1645: filter: grep -v "CG or CGNE: variant"
1646: nsize: 7
1647: 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
1649: test:
1650: suffix: p4est_test_q2_nonconformal_parallel_parmetis
1651: requires: parmetis p4est
1652: nsize: 4
1653: 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
1655: test:
1656: suffix: p4est_exact_q2_conformal_serial
1657: requires: p4est !single !complex !__float128
1658: 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
1660: test:
1661: suffix: p4est_exact_q2_conformal_parallel
1662: requires: p4est !single !complex !__float128
1663: nsize: 4
1664: 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
1666: test:
1667: suffix: p4est_exact_q2_conformal_parallel_parmetis
1668: requires: parmetis p4est !single
1669: nsize: 4
1670: 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
1672: test:
1673: suffix: p4est_exact_q2_nonconformal_serial
1674: requires: p4est
1675: 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
1677: test:
1678: suffix: p4est_exact_q2_nonconformal_parallel
1679: requires: p4est
1680: nsize: 7
1681: 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
1683: test:
1684: suffix: p4est_exact_q2_nonconformal_parallel_parmetis
1685: requires: parmetis p4est
1686: nsize: 4
1687: 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
1689: test:
1690: suffix: p4est_full_q2_nonconformal_serial
1691: requires: p4est !single
1692: filter: grep -v "variant HERMITIAN"
1693: 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
1695: test:
1696: suffix: p4est_full_q2_nonconformal_parallel
1697: requires: p4est !single
1698: filter: grep -v "variant HERMITIAN"
1699: nsize: 7
1700: 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
1702: test:
1703: suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1704: requires: p4est
1705: filter: grep -v "variant HERMITIAN"
1706: nsize: 7
1707: 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
1709: test:
1710: suffix: p4est_full_q2_nonconformal_parallel_bddc
1711: requires: p4est
1712: filter: grep -v "variant HERMITIAN"
1713: nsize: 7
1714: 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
1716: test:
1717: TODO: broken
1718: suffix: p4est_fas_q2_conformal_serial
1719: requires: p4est !complex !__float128
1720: 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
1722: test:
1723: TODO: broken
1724: suffix: p4est_fas_q2_nonconformal_serial
1725: requires: p4est
1726: 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
1728: test:
1729: suffix: fas_newton_0_p4est
1730: requires: p4est !single !__float128
1731: 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
1733: # Full solve simplicial AMR
1734: test:
1735: suffix: tri_p1_adapt_0
1736: requires: pragmatic
1737: 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
1739: test:
1740: suffix: tri_p1_adapt_1
1741: requires: pragmatic
1742: 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
1744: test:
1745: suffix: tri_p1_adapt_analytic_0
1746: requires: pragmatic
1747: 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
1749: # Full solve tensor AMR
1750: test:
1751: suffix: quad_q1_adapt_0
1752: requires: p4est
1753: 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
1754: filter: grep -v DM_
1756: test:
1757: suffix: amr_0
1758: nsize: 5
1759: 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
1761: test:
1762: suffix: amr_1
1763: requires: p4est !complex
1764: 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
1766: test:
1767: suffix: p4est_solve_bddc
1768: requires: p4est !complex
1769: 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
1770: nsize: 4
1772: test:
1773: suffix: p4est_solve_fas
1774: requires: p4est
1775: 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
1776: nsize: 4
1777: TODO: identical machine two runs produce slightly different solver trackers
1779: test:
1780: suffix: p4est_convergence_test_1
1781: requires: p4est
1782: 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
1783: nsize: 4
1785: test:
1786: suffix: p4est_convergence_test_2
1787: requires: p4est
1788: 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
1790: test:
1791: suffix: p4est_convergence_test_3
1792: requires: p4est
1793: 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
1795: test:
1796: suffix: p4est_convergence_test_4
1797: requires: p4est
1798: 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
1799: timeoutfactor: 5
1801: # Serial tests with GLVis visualization
1802: test:
1803: suffix: glvis_2d_tet_p1
1804: 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
1805: test:
1806: suffix: glvis_2d_tet_p2
1807: 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
1808: test:
1809: suffix: glvis_2d_hex_p1
1810: args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -simplex 0 -dm_refine 1
1811: test:
1812: suffix: glvis_2d_hex_p2
1813: args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -simplex 0 -dm_refine 1
1814: test:
1815: suffix: glvis_2d_hex_p2_p4est
1816: requires: p4est
1817: 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
1819: TEST*/