Actual source code: ex12.c
petsc-main 2021-04-15
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 :~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, COEFF_CHECKERBOARD_0, COEFF_CHECKERBOARD_1} 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: /* Reproducing tests from SISC 40(3), pp. A1473-A1493, 2018 */
53: PetscInt div; /* Number of divisions */
54: PetscInt k; /* Parameter for checkerboard coefficient */
55: PetscInt *kgrid; /* Random parameter grid */
56: /* Solver */
57: PC pcmg; /* This is needed for error monitoring */
58: PetscBool checkksp; /* Whether to check the KSPSolve for runType == RUN_TEST */
59: } AppCtx;
61: static PetscErrorCode zero(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
62: {
63: u[0] = 0.0;
64: return 0;
65: }
67: static PetscErrorCode ecks(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
68: {
69: u[0] = x[0];
70: return 0;
71: }
73: /*
74: In 2D for Dirichlet conditions, we use exact solution:
76: u = x^2 + y^2
77: f = 4
79: so that
81: -\Delta u + f = -4 + 4 = 0
83: For Neumann conditions, we have
85: -\nabla u \cdot -\hat y |_{y=0} = (2y)|_{y=0} = 0 (bottom)
86: -\nabla u \cdot \hat y |_{y=1} = -(2y)|_{y=1} = -2 (top)
87: -\nabla u \cdot -\hat x |_{x=0} = (2x)|_{x=0} = 0 (left)
88: -\nabla u \cdot \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)
90: Which we can express as
92: \nabla u \cdot \hat n|_\Gamma = {2 x, 2 y} \cdot \hat n = 2 (x + y)
94: The boundary integral of this solution is (assuming we are not orienting the edges)
96: \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
97: */
98: static PetscErrorCode quadratic_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
99: {
100: *u = x[0]*x[0] + x[1]*x[1];
101: return 0;
102: }
104: static void quadratic_u_field_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
105: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
106: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
107: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
108: {
109: uexact[0] = a[0];
110: }
112: static PetscErrorCode circle_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
113: {
114: const PetscReal alpha = 500.;
115: const PetscReal radius2 = PetscSqr(0.15);
116: const PetscReal r2 = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5);
117: const PetscReal xi = alpha*(radius2 - r2);
119: *u = PetscTanhScalar(xi) + 1.0;
120: return 0;
121: }
123: static PetscErrorCode cross_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
124: {
125: const PetscReal alpha = 50*4;
126: const PetscReal xy = (x[0]-0.5)*(x[1]-0.5);
128: *u = PetscSinReal(alpha*xy) * (alpha*PetscAbsReal(xy) < 2*PETSC_PI ? 1 : 0.01);
129: return 0;
130: }
132: static void f0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
133: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
134: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
135: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
136: {
137: f0[0] = 4.0;
138: }
140: static void f0_circle_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
141: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
142: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
143: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
144: {
145: const PetscReal alpha = 500.;
146: const PetscReal radius2 = PetscSqr(0.15);
147: const PetscReal r2 = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5);
148: const PetscReal xi = alpha*(radius2 - r2);
150: f0[0] = (-4.0*alpha - 8.0*PetscSqr(alpha)*r2*PetscTanhReal(xi)) * PetscSqr(1.0/PetscCoshReal(xi));
151: }
153: static void f0_cross_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
154: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
155: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
156: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
157: {
158: const PetscReal alpha = 50*4;
159: const PetscReal xy = (x[0]-0.5)*(x[1]-0.5);
161: f0[0] = PetscSinReal(alpha*xy) * (alpha*PetscAbsReal(xy) < 2*PETSC_PI ? 1 : 0.01);
162: }
164: static void f0_checkerboard_0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
165: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
166: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
167: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
168: {
169: f0[0] = -20.0*PetscExpReal(-(PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5)));
170: }
172: static void f0_bd_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
173: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
174: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
175: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
176: {
177: PetscInt d;
178: for (d = 0, f0[0] = 0.0; d < dim; ++d) f0[0] += -n[d]*2.0*x[d];
179: }
181: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
182: static void f1_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
183: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
184: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
185: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
186: {
187: PetscInt d;
188: for (d = 0; d < dim; ++d) f1[d] = u_x[d];
189: }
191: /* < \nabla v, \nabla u + {\nabla u}^T >
192: This just gives \nabla u, give the perdiagonal for the transpose */
193: static void g3_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
194: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
195: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
196: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
197: {
198: PetscInt d;
199: for (d = 0; d < dim; ++d) g3[d*dim+d] = 1.0;
200: }
202: /*
203: In 2D for x periodicity and y Dirichlet conditions, we use exact solution:
205: u = sin(2 pi x)
206: f = -4 pi^2 sin(2 pi x)
208: so that
210: -\Delta u + f = 4 pi^2 sin(2 pi x) - 4 pi^2 sin(2 pi x) = 0
211: */
212: static PetscErrorCode xtrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
213: {
214: *u = PetscSinReal(2.0*PETSC_PI*x[0]);
215: return 0;
216: }
218: static void f0_xtrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
219: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
220: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
221: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
222: {
223: f0[0] = -4.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
224: }
226: /*
227: In 2D for x-y periodicity, we use exact solution:
229: u = sin(2 pi x) sin(2 pi y)
230: f = -8 pi^2 sin(2 pi x)
232: so that
234: -\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
235: */
236: static PetscErrorCode xytrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
237: {
238: *u = PetscSinReal(2.0*PETSC_PI*x[0])*PetscSinReal(2.0*PETSC_PI*x[1]);
239: return 0;
240: }
242: static void f0_xytrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
243: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
244: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
245: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
246: {
247: f0[0] = -8.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
248: }
250: /*
251: In 2D for Dirichlet conditions with a variable coefficient, we use exact solution:
253: u = x^2 + y^2
254: f = 6 (x + y)
255: nu = (x + y)
257: so that
259: -\div \nu \grad u + f = -6 (x + y) + 6 (x + y) = 0
260: */
261: static PetscErrorCode nu_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
262: {
263: *u = x[0] + x[1];
264: return 0;
265: }
267: static PetscErrorCode checkerboardCoeff(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
268: {
269: AppCtx *user = (AppCtx *) ctx;
270: PetscInt div = user->div;
271: PetscInt k = user->k;
272: PetscInt mask = 0, ind = 0, d;
275: for (d = 0; d < dim; ++d) mask = (mask + (PetscInt) (x[d]*div)) % 2;
276: if (user->kgrid) {
277: for (d = 0; d < dim; ++d) {
278: if (d > 0) ind *= dim;
279: ind += (PetscInt) (x[d]*div);
280: }
281: k = user->kgrid[ind];
282: }
283: u[0] = mask ? 1.0 : PetscPowRealInt(10.0, -k);
284: return(0);
285: }
287: void f0_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
288: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
289: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
290: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
291: {
292: f0[0] = 6.0*(x[0] + x[1]);
293: }
295: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
296: void f1_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
297: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
298: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
299: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
300: {
301: PetscInt d;
302: for (d = 0; d < dim; ++d) f1[d] = (x[0] + x[1])*u_x[d];
303: }
305: void f1_field_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
306: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
307: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
308: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
309: {
310: PetscInt d;
311: for (d = 0; d < dim; ++d) f1[d] = a[0]*u_x[d];
312: }
314: /* < \nabla v, \nabla u + {\nabla u}^T >
315: This just gives \nabla u, give the perdiagonal for the transpose */
316: void g3_analytic_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
317: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
318: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
319: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
320: {
321: PetscInt d;
322: for (d = 0; d < dim; ++d) g3[d*dim+d] = x[0] + x[1];
323: }
325: void g3_field_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
326: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
327: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
328: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
329: {
330: PetscInt d;
331: for (d = 0; d < dim; ++d) g3[d*dim+d] = a[0];
332: }
334: /*
335: In 2D for Dirichlet conditions with a nonlinear coefficient (p-Laplacian with p = 4), we use exact solution:
337: u = x^2 + y^2
338: f = 16 (x^2 + y^2)
339: nu = 1/2 |grad u|^2
341: so that
343: -\div \nu \grad u + f = -16 (x^2 + y^2) + 16 (x^2 + y^2) = 0
344: */
345: void f0_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
346: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
347: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
348: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
349: {
350: f0[0] = 16.0*(x[0]*x[0] + x[1]*x[1]);
351: }
353: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
354: void f1_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
355: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
356: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
357: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
358: {
359: PetscScalar nu = 0.0;
360: PetscInt d;
361: for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
362: for (d = 0; d < dim; ++d) f1[d] = 0.5*nu*u_x[d];
363: }
365: /*
366: grad (u + eps w) - grad u = eps grad w
368: 1/2 |grad (u + eps w)|^2 grad (u + eps w) - 1/2 |grad u|^2 grad u
369: = 1/2 (|grad u|^2 + 2 eps <grad u,grad w>) (grad u + eps grad w) - 1/2 |grad u|^2 grad u
370: = 1/2 (eps |grad u|^2 grad w + 2 eps <grad u,grad w> grad u)
371: = eps (1/2 |grad u|^2 grad w + grad u <grad u,grad w>)
372: */
373: void g3_analytic_nonlinear_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
374: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
375: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
376: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
377: {
378: PetscScalar nu = 0.0;
379: PetscInt d, e;
380: for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
381: for (d = 0; d < dim; ++d) {
382: g3[d*dim+d] = 0.5*nu;
383: for (e = 0; e < dim; ++e) {
384: g3[d*dim+e] += u_x[d]*u_x[e];
385: }
386: }
387: }
389: /*
390: In 3D for Dirichlet conditions we use exact solution:
392: u = 2/3 (x^2 + y^2 + z^2)
393: f = 4
395: so that
397: -\Delta u + f = -2/3 * 6 + 4 = 0
399: For Neumann conditions, we have
401: -\nabla u \cdot -\hat z |_{z=0} = (2z)|_{z=0} = 0 (bottom)
402: -\nabla u \cdot \hat z |_{z=1} = -(2z)|_{z=1} = -2 (top)
403: -\nabla u \cdot -\hat y |_{y=0} = (2y)|_{y=0} = 0 (front)
404: -\nabla u \cdot \hat y |_{y=1} = -(2y)|_{y=1} = -2 (back)
405: -\nabla u \cdot -\hat x |_{x=0} = (2x)|_{x=0} = 0 (left)
406: -\nabla u \cdot \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)
408: Which we can express as
410: \nabla u \cdot \hat n|_\Gamma = {2 x, 2 y, 2z} \cdot \hat n = 2 (x + y + z)
411: */
412: static PetscErrorCode quadratic_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
413: {
414: *u = 2.0*(x[0]*x[0] + x[1]*x[1] + x[2]*x[2])/3.0;
415: return 0;
416: }
418: static void quadratic_u_field_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
419: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
420: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
421: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
422: {
423: uexact[0] = a[0];
424: }
426: static void bd_integral_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
427: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
428: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
429: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar *uint)
430: {
431: uint[0] = u[0];
432: }
434: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
435: {
436: const char *bcTypes[3] = {"neumann", "dirichlet", "none"};
437: const char *runTypes[4] = {"full", "exact", "test", "perf"};
438: const char *coeffTypes[8] = {"none", "analytic", "field", "nonlinear", "circle", "cross", "checkerboard_0", "checkerboard_1"};
439: PetscInt bd, bc, run, coeff, n;
440: PetscBool rand = PETSC_FALSE, flg;
444: options->debug = 0;
445: options->runType = RUN_FULL;
446: options->dim = 2;
447: options->periodicity[0] = DM_BOUNDARY_NONE;
448: options->periodicity[1] = DM_BOUNDARY_NONE;
449: options->periodicity[2] = DM_BOUNDARY_NONE;
450: options->cells[0] = 2;
451: options->cells[1] = 2;
452: options->cells[2] = 2;
453: options->filename[0] = '\0';
454: options->interpolate = PETSC_TRUE;
455: options->refinementLimit = 0.0;
456: options->bcType = DIRICHLET;
457: options->variableCoefficient = COEFF_NONE;
458: options->fieldBC = PETSC_FALSE;
459: options->jacobianMF = PETSC_FALSE;
460: options->showInitial = PETSC_FALSE;
461: options->showSolution = PETSC_FALSE;
462: options->restart = PETSC_FALSE;
463: options->viewHierarchy = PETSC_FALSE;
464: options->simplex = PETSC_TRUE;
465: options->quiet = PETSC_FALSE;
466: options->nonzInit = PETSC_FALSE;
467: options->bdIntegral = PETSC_FALSE;
468: options->checkksp = PETSC_FALSE;
469: options->div = 4;
470: options->k = 1;
471: options->kgrid = NULL;
473: PetscOptionsBegin(comm, "", "Poisson Problem Options", "DMPLEX");
474: PetscOptionsInt("-debug", "The debugging level", "ex12.c", options->debug, &options->debug, NULL);
475: run = options->runType;
476: PetscOptionsEList("-run_type", "The run type", "ex12.c", runTypes, 4, runTypes[options->runType], &run, NULL);
478: options->runType = (RunType) run;
480: PetscOptionsInt("-dim", "The topological mesh dimension", "ex12.c", options->dim, &options->dim, NULL);
481: bd = options->periodicity[0];
482: PetscOptionsEList("-x_periodicity", "The x-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[0]], &bd, NULL);
483: options->periodicity[0] = (DMBoundaryType) bd;
484: bd = options->periodicity[1];
485: PetscOptionsEList("-y_periodicity", "The y-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[1]], &bd, NULL);
486: options->periodicity[1] = (DMBoundaryType) bd;
487: bd = options->periodicity[2];
488: PetscOptionsEList("-z_periodicity", "The z-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[2]], &bd, NULL);
489: options->periodicity[2] = (DMBoundaryType) bd;
490: n = 3;
491: PetscOptionsIntArray("-cells", "The initial mesh division", "ex12.c", options->cells, &n, NULL);
492: PetscOptionsString("-f", "Mesh filename to read", "ex12.c", options->filename, options->filename, sizeof(options->filename), &flg);
493: PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex12.c", options->interpolate, &options->interpolate, NULL);
494: PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex12.c", options->refinementLimit, &options->refinementLimit, NULL);
495: bc = options->bcType;
496: PetscOptionsEList("-bc_type","Type of boundary condition","ex12.c",bcTypes,3,bcTypes[options->bcType],&bc,NULL);
497: options->bcType = (BCType) bc;
498: coeff = options->variableCoefficient;
499: PetscOptionsEList("-variable_coefficient","Type of variable coefficent","ex12.c",coeffTypes,8,coeffTypes[options->variableCoefficient],&coeff,NULL);
500: options->variableCoefficient = (CoeffType) coeff;
502: PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL);
503: PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL);
504: PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL);
505: PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL);
506: PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL);
507: PetscOptionsBool("-dm_view_hierarchy", "View the coarsened hierarchy", "ex12.c", options->viewHierarchy, &options->viewHierarchy, NULL);
508: PetscOptionsBool("-simplex", "Simplicial (true) or tensor (false) mesh", "ex12.c", options->simplex, &options->simplex, NULL);
509: PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL);
510: PetscOptionsBool("-nonzero_initial_guess", "nonzero initial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL);
511: PetscOptionsBool("-bd_integral", "Compute the integral of the solution on the boundary", "ex12.c", options->bdIntegral, &options->bdIntegral, NULL);
512: if (options->runType == RUN_TEST) {
513: PetscOptionsBool("-run_test_check_ksp", "Check solution of KSP", "ex12.c", options->checkksp, &options->checkksp, NULL);
514: }
515: PetscOptionsInt("-div", "The number of division for the checkerboard coefficient", "ex12.c", options->div, &options->div, NULL);
516: PetscOptionsInt("-k", "The exponent for the checkerboard coefficient", "ex12.c", options->k, &options->k, NULL);
517: PetscOptionsBool("-k_random", "Assign random k values to checkerboard", "ex12.c", rand, &rand, NULL);
518: PetscOptionsEnd();
519: PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
521: if (rand) {
522: PetscRandom r;
523: PetscReal val;
524: PetscInt N = PetscPowInt(options->div, options->dim), i;
526: PetscMalloc1(N, &options->kgrid);
527: PetscRandomCreate(PETSC_COMM_SELF, &r);
528: PetscRandomSetFromOptions(r);
529: PetscRandomSetInterval(r, 0.0, options->k);
530: PetscRandomSetSeed(r, 1973);
531: PetscRandomSeed(r);
532: for (i = 0; i < N; ++i) {
533: PetscRandomGetValueReal(r, &val);
534: options->kgrid[i] = 1 + (PetscInt) val;
535: }
536: PetscRandomDestroy(&r);
537: }
538: return(0);
539: }
541: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
542: {
543: DM plex;
544: DMLabel label;
548: DMCreateLabel(dm, name);
549: DMGetLabel(dm, name, &label);
550: DMConvert(dm, DMPLEX, &plex);
551: DMPlexMarkBoundaryFaces(plex, 1, label);
552: DMDestroy(&plex);
553: return(0);
554: }
556: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
557: {
558: PetscInt dim = user->dim;
559: const char *filename = user->filename;
560: PetscBool interpolate = user->interpolate;
561: PetscReal refinementLimit = user->refinementLimit;
562: size_t len;
566: PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
567: PetscStrlen(filename, &len);
568: if (!len) {
569: PetscInt d;
571: if (user->periodicity[0] || user->periodicity[1] || user->periodicity[2]) for (d = 0; d < dim; ++d) user->cells[d] = PetscMax(user->cells[d], 3);
572: DMPlexCreateBoxMesh(comm, dim, user->simplex, user->cells, NULL, NULL, user->periodicity, interpolate, dm);
573: PetscObjectSetName((PetscObject) *dm, "Mesh");
574: } else {
575: DMPlexCreateFromFile(comm, filename, interpolate, dm);
576: DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
577: }
578: {
579: PetscPartitioner part;
580: DM refinedMesh = NULL;
581: DM distributedMesh = NULL;
583: /* Refine mesh using a volume constraint */
584: if (refinementLimit > 0.0) {
585: DMPlexSetRefinementLimit(*dm, refinementLimit);
586: DMRefine(*dm, comm, &refinedMesh);
587: if (refinedMesh) {
588: const char *name;
590: PetscObjectGetName((PetscObject) *dm, &name);
591: PetscObjectSetName((PetscObject) refinedMesh, name);
592: DMDestroy(dm);
593: *dm = refinedMesh;
594: }
595: }
596: /* Distribute mesh over processes */
597: DMPlexGetPartitioner(*dm,&part);
598: PetscPartitionerSetFromOptions(part);
599: DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
600: if (distributedMesh) {
601: DMDestroy(dm);
602: *dm = distributedMesh;
603: }
604: }
605: if (interpolate) {
606: if (user->bcType == NEUMANN) {
607: DMLabel label;
609: DMCreateLabel(*dm, "boundary");
610: DMGetLabel(*dm, "boundary", &label);
611: DMPlexMarkBoundaryFaces(*dm, 1, label);
612: } else if (user->bcType == DIRICHLET) {
613: PetscBool hasLabel;
615: DMHasLabel(*dm,"marker",&hasLabel);
616: if (!hasLabel) {
617: //DMSetUp(*dm);
618: CreateBCLabel(*dm, "marker");
619: }
620: }
621: }
622: {
623: char convType[256];
624: PetscBool flg;
626: PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
627: PetscOptionsFList("-dm_plex_convert_type","Convert DMPlex to another format","ex12",DMList,DMPLEX,convType,256,&flg);
628: PetscOptionsEnd();
629: if (flg) {
630: DM dmConv;
632: DMConvert(*dm,convType,&dmConv);
633: if (dmConv) {
634: DMDestroy(dm);
635: *dm = dmConv;
636: }
637: }
638: }
639: DMLocalizeCoordinates(*dm); /* needed for periodic */
640: DMSetFromOptions(*dm);
641: DMViewFromOptions(*dm, NULL, "-dm_view");
642: if (user->viewHierarchy) {
643: DM cdm = *dm;
644: PetscInt i = 0;
645: char buf[256];
647: while (cdm) {
648: DMSetUp(cdm);
649: DMGetCoarseDM(cdm, &cdm);
650: ++i;
651: }
652: cdm = *dm;
653: while (cdm) {
654: PetscViewer viewer;
655: PetscBool isHDF5, isVTK;
657: --i;
658: PetscViewerCreate(comm,&viewer);
659: PetscViewerSetType(viewer,PETSCVIEWERHDF5);
660: PetscViewerSetOptionsPrefix(viewer,"hierarchy_");
661: PetscViewerSetFromOptions(viewer);
662: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&isHDF5);
663: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERVTK,&isVTK);
664: if (isHDF5) {
665: PetscSNPrintf(buf, 256, "ex12-%d.h5", i);
666: } else if (isVTK) {
667: PetscSNPrintf(buf, 256, "ex12-%d.vtu", i);
668: PetscViewerPushFormat(viewer,PETSC_VIEWER_VTK_VTU);
669: } else {
670: PetscSNPrintf(buf, 256, "ex12-%d", i);
671: }
672: PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
673: PetscViewerFileSetName(viewer,buf);
674: DMView(cdm, viewer);
675: PetscViewerDestroy(&viewer);
676: DMGetCoarseDM(cdm, &cdm);
677: }
678: }
679: PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
680: return(0);
681: }
683: static PetscErrorCode SetupProblem(DM dm, AppCtx *user)
684: {
685: PetscDS ds;
686: DMLabel label;
687: PetscWeakForm wf;
688: const PetscInt id = 1;
689: PetscInt bd;
693: DMGetDS(dm, &ds);
694: switch (user->variableCoefficient) {
695: case COEFF_NONE:
696: if (user->periodicity[0]) {
697: if (user->periodicity[1]) {
698: PetscDSSetResidual(ds, 0, f0_xytrig_u, f1_u);
699: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
700: } else {
701: PetscDSSetResidual(ds, 0, f0_xtrig_u, f1_u);
702: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
703: }
704: } else {
705: PetscDSSetResidual(ds, 0, f0_u, f1_u);
706: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
707: }
708: break;
709: case COEFF_ANALYTIC:
710: PetscDSSetResidual(ds, 0, f0_analytic_u, f1_analytic_u);
711: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_analytic_uu);
712: break;
713: case COEFF_FIELD:
714: PetscDSSetResidual(ds, 0, f0_analytic_u, f1_field_u);
715: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_field_uu);
716: break;
717: case COEFF_NONLINEAR:
718: PetscDSSetResidual(ds, 0, f0_analytic_nonlinear_u, f1_analytic_nonlinear_u);
719: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_analytic_nonlinear_uu);
720: break;
721: case COEFF_CIRCLE:
722: PetscDSSetResidual(ds, 0, f0_circle_u, f1_u);
723: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
724: break;
725: case COEFF_CROSS:
726: PetscDSSetResidual(ds, 0, f0_cross_u, f1_u);
727: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
728: break;
729: case COEFF_CHECKERBOARD_0:
730: PetscDSSetResidual(ds, 0, f0_checkerboard_0_u, f1_field_u);
731: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_field_uu);
732: break;
733: default: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid variable coefficient type %d", user->variableCoefficient);
734: }
735: switch (user->dim) {
736: case 2:
737: switch (user->variableCoefficient) {
738: case COEFF_CIRCLE:
739: user->exactFuncs[0] = circle_u_2d;break;
740: case COEFF_CROSS:
741: user->exactFuncs[0] = cross_u_2d;break;
742: case COEFF_CHECKERBOARD_0:
743: user->exactFuncs[0] = zero;break;
744: default:
745: if (user->periodicity[0]) {
746: if (user->periodicity[1]) {
747: user->exactFuncs[0] = xytrig_u_2d;
748: } else {
749: user->exactFuncs[0] = xtrig_u_2d;
750: }
751: } else {
752: user->exactFuncs[0] = quadratic_u_2d;
753: user->exactFields[0] = quadratic_u_field_2d;
754: }
755: }
756: if (user->bcType == NEUMANN) {
757: DMGetLabel(dm, "boundary", &label);
758: DMAddBoundary(dm, DM_BC_NATURAL, "wall", label, 1, &id, 0, 0, NULL, NULL, NULL, user, &bd);
759: PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL);
760: PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, f0_bd_u, 0, NULL);
761: }
762: break;
763: case 3:
764: user->exactFuncs[0] = quadratic_u_3d;
765: user->exactFields[0] = quadratic_u_field_3d;
766: if (user->bcType == NEUMANN) {
767: DMGetLabel(dm, "boundary", &label);
768: DMAddBoundary(dm, DM_BC_NATURAL, "wall", label, 1, &id, 0, 0, NULL, NULL, NULL, user, &bd);
769: PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL);
770: PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, f0_bd_u, 0, NULL);
771: }
772: break;
773: default:
774: SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %d", user->dim);
775: }
776: /* Setup constants */
777: switch (user->variableCoefficient) {
778: case COEFF_CHECKERBOARD_0:
779: {
780: PetscScalar constants[2];
782: constants[0] = user->div;
783: constants[1] = user->k;
784: PetscDSSetConstants(ds, 2, constants);
785: }
786: break;
787: default: break;
788: }
789: PetscDSSetExactSolution(ds, 0, user->exactFuncs[0], user);
790: /* Setup Boundary Conditions */
791: if (user->bcType == DIRICHLET) {
792: DMGetLabel(dm, "marker", &label);
793: if (!label) {
794: /* Right now, p4est cannot create labels immediately */
795: PetscDSAddBoundaryByName(ds, user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL, "wall", "marker", 1, &id, 0, 0, NULL, user->fieldBC ? (void (*)(void)) user->exactFields[0] : (void (*)(void)) user->exactFuncs[0], NULL, user, NULL);
796: } else {
797: DMAddBoundary(dm, user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL, "wall", label, 1, &id, 0, 0, NULL, user->fieldBC ? (void (*)(void)) user->exactFields[0] : (void (*)(void)) user->exactFuncs[0], NULL, user, NULL);
798: }
799: }
800: return(0);
801: }
803: static PetscErrorCode SetupMaterial(DM dm, DM dmAux, AppCtx *user)
804: {
805: PetscErrorCode (*matFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {nu_2d};
806: void *ctx[1];
807: Vec nu;
808: PetscErrorCode ierr;
811: ctx[0] = user;
812: if (user->variableCoefficient == COEFF_CHECKERBOARD_0) {matFuncs[0] = checkerboardCoeff;}
813: DMCreateLocalVector(dmAux, &nu);
814: PetscObjectSetName((PetscObject) nu, "Coefficient");
815: DMProjectFunctionLocal(dmAux, 0.0, matFuncs, ctx, INSERT_ALL_VALUES, nu);
816: DMSetAuxiliaryVec(dm, NULL, 0, nu);
817: VecDestroy(&nu);
818: return(0);
819: }
821: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
822: {
823: PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
824: Vec uexact;
825: PetscInt dim;
829: DMGetDimension(dm, &dim);
830: if (dim == 2) bcFuncs[0] = quadratic_u_2d;
831: else bcFuncs[0] = quadratic_u_3d;
832: DMCreateLocalVector(dmAux, &uexact);
833: DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
834: DMSetAuxiliaryVec(dm, NULL, 0, uexact);
835: VecDestroy(&uexact);
836: return(0);
837: }
839: static PetscErrorCode SetupAuxDM(DM dm, PetscFE feAux, AppCtx *user)
840: {
841: DM dmAux, coordDM;
845: /* MUST call DMGetCoordinateDM() in order to get p4est setup if present */
846: DMGetCoordinateDM(dm, &coordDM);
847: if (!feAux) return(0);
848: DMClone(dm, &dmAux);
849: DMSetCoordinateDM(dmAux, coordDM);
850: DMSetField(dmAux, 0, NULL, (PetscObject) feAux);
851: DMCreateDS(dmAux);
852: if (user->fieldBC) {SetupBC(dm, dmAux, user);}
853: else {SetupMaterial(dm, dmAux, user);}
854: DMDestroy(&dmAux);
855: return(0);
856: }
858: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
859: {
860: DM cdm = dm;
861: const PetscInt dim = user->dim;
862: PetscFE fe, feAux = NULL;
863: PetscBool simplex = user->simplex;
864: MPI_Comm comm;
868: /* Create finite element for each field and auxiliary field */
869: PetscObjectGetComm((PetscObject) dm, &comm);
870: PetscFECreateDefault(comm, dim, 1, simplex, NULL, -1, &fe);
871: PetscObjectSetName((PetscObject) fe, "potential");
872: if (user->variableCoefficient == COEFF_FIELD || user->variableCoefficient == COEFF_CHECKERBOARD_0) {
873: PetscFECreateDefault(comm, dim, 1, simplex, "mat_", -1, &feAux);
874: PetscObjectSetName((PetscObject) feAux, "coefficient");
875: PetscFECopyQuadrature(fe, feAux);
876: } else if (user->fieldBC) {
877: PetscFECreateDefault(comm, dim, 1, simplex, "bc_", -1, &feAux);
878: PetscFECopyQuadrature(fe, feAux);
879: }
880: /* Set discretization and boundary conditions for each mesh */
881: DMSetField(dm, 0, NULL, (PetscObject) fe);
882: DMCreateDS(dm);
883: SetupProblem(dm, user);
884: while (cdm) {
885: SetupAuxDM(cdm, feAux, user);
886: if (user->bcType == DIRICHLET && user->interpolate) {
887: PetscBool hasLabel;
889: DMHasLabel(cdm, "marker", &hasLabel);
890: if (!hasLabel) {CreateBCLabel(cdm, "marker");}
891: }
892: DMCopyDisc(dm, cdm);
893: DMGetCoarseDM(cdm, &cdm);
894: }
895: PetscFEDestroy(&fe);
896: PetscFEDestroy(&feAux);
897: return(0);
898: }
900: #include "petsc/private/petscimpl.h"
902: /*
903: MonitorError - Outputs the error at each iteration of an iterative solver.
905: Collective on KSP
907: Input Parameters:
908: + ksp - the KSP
909: . its - iteration number
910: . rnorm - 2-norm, preconditioned residual value (may be estimated).
911: - ctx - monitor context
913: Level: intermediate
915: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual(), KSPMonitorResidual()
916: */
917: static PetscErrorCode MonitorError(KSP ksp, PetscInt its, PetscReal rnorm, void *ctx)
918: {
919: AppCtx *user = (AppCtx *) ctx;
920: DM dm;
921: Vec du = NULL, r;
922: PetscInt level = 0;
923: PetscBool hasLevel;
924: #if defined(PETSC_HAVE_HDF5)
925: PetscViewer viewer;
926: char buf[256];
927: #endif
931: KSPGetDM(ksp, &dm);
932: /* Calculate solution */
933: {
934: PC pc = user->pcmg; /* The MG PC */
935: DM fdm = NULL, cdm = NULL;
936: KSP fksp, cksp;
937: Vec fu, cu = NULL;
938: PetscInt levels, l;
940: KSPBuildSolution(ksp, NULL, &du);
941: PetscObjectComposedDataGetInt((PetscObject) ksp, PetscMGLevelId, level, hasLevel);
942: PCMGGetLevels(pc, &levels);
943: PCMGGetSmoother(pc, levels-1, &fksp);
944: KSPBuildSolution(fksp, NULL, &fu);
945: for (l = levels-1; l > level; --l) {
946: Mat R;
947: Vec s;
949: PCMGGetSmoother(pc, l-1, &cksp);
950: KSPGetDM(cksp, &cdm);
951: DMGetGlobalVector(cdm, &cu);
952: PCMGGetRestriction(pc, l, &R);
953: PCMGGetRScale(pc, l, &s);
954: MatRestrict(R, fu, cu);
955: VecPointwiseMult(cu, cu, s);
956: if (l < levels-1) {DMRestoreGlobalVector(fdm, &fu);}
957: fdm = cdm;
958: fu = cu;
959: }
960: if (levels-1 > level) {
961: VecAXPY(du, 1.0, cu);
962: DMRestoreGlobalVector(cdm, &cu);
963: }
964: }
965: /* Calculate error */
966: DMGetGlobalVector(dm, &r);
967: DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
968: VecAXPY(r,-1.0,du);
969: PetscObjectSetName((PetscObject) r, "solution error");
970: /* View error */
971: #if defined(PETSC_HAVE_HDF5)
972: PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
973: PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
974: VecView(r, viewer);
975: PetscViewerDestroy(&viewer);
976: #endif
977: DMRestoreGlobalVector(dm, &r);
978: return(0);
979: }
981: /*@C
982: SNESMonitorError - Outputs the error at each iteration of an iterative solver.
984: Collective on SNES
986: Input Parameters:
987: + snes - the SNES
988: . its - iteration number
989: . rnorm - 2-norm of residual
990: - ctx - user context
992: Level: intermediate
994: .seealso: SNESMonitorDefault(), SNESMonitorSet(), SNESMonitorSolution()
995: @*/
996: static PetscErrorCode SNESMonitorError(SNES snes, PetscInt its, PetscReal rnorm, void *ctx)
997: {
998: AppCtx *user = (AppCtx *) ctx;
999: DM dm;
1000: Vec u, r;
1001: PetscInt level = -1;
1002: PetscBool hasLevel;
1003: #if defined(PETSC_HAVE_HDF5)
1004: PetscViewer viewer;
1005: #endif
1006: char buf[256];
1010: SNESGetDM(snes, &dm);
1011: /* Calculate error */
1012: SNESGetSolution(snes, &u);
1013: DMGetGlobalVector(dm, &r);
1014: PetscObjectSetName((PetscObject) r, "solution error");
1015: DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
1016: VecAXPY(r, -1.0, u);
1017: /* View error */
1018: PetscObjectComposedDataGetInt((PetscObject) snes, PetscMGLevelId, level, hasLevel);
1019: PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
1020: #if defined(PETSC_HAVE_HDF5)
1021: PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
1022: VecView(r, viewer);
1023: PetscViewerDestroy(&viewer);
1024: /* Cleanup */
1025: DMRestoreGlobalVector(dm, &r);
1026: return(0);
1027: #else
1028: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"You need to configure with --download-hdf5");
1029: #endif
1030: }
1032: int main(int argc, char **argv)
1033: {
1034: DM dm; /* Problem specification */
1035: SNES snes; /* nonlinear solver */
1036: Vec u; /* solution vector */
1037: Mat A,J; /* Jacobian matrix */
1038: MatNullSpace nullSpace; /* May be necessary for Neumann conditions */
1039: AppCtx user; /* user-defined work context */
1040: JacActionCtx userJ; /* context for Jacobian MF action */
1041: PetscReal error = 0.0; /* L_2 error in the solution */
1042: PetscBool isFAS;
1045: PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
1046: ProcessOptions(PETSC_COMM_WORLD, &user);
1047: SNESCreate(PETSC_COMM_WORLD, &snes);
1048: CreateMesh(PETSC_COMM_WORLD, &user, &dm);
1049: SNESSetDM(snes, dm);
1050: DMSetApplicationContext(dm, &user);
1052: PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
1053: SetupDiscretization(dm, &user);
1055: DMCreateGlobalVector(dm, &u);
1056: PetscObjectSetName((PetscObject) u, "potential");
1058: DMCreateMatrix(dm, &J);
1059: if (user.jacobianMF) {
1060: PetscInt M, m, N, n;
1062: MatGetSize(J, &M, &N);
1063: MatGetLocalSize(J, &m, &n);
1064: MatCreate(PETSC_COMM_WORLD, &A);
1065: MatSetSizes(A, m, n, M, N);
1066: MatSetType(A, MATSHELL);
1067: MatSetUp(A);
1068: #if 0
1069: MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction);
1070: #endif
1072: userJ.dm = dm;
1073: userJ.J = J;
1074: userJ.user = &user;
1076: DMCreateLocalVector(dm, &userJ.u);
1077: if (user.fieldBC) {DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);}
1078: else {DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);}
1079: MatShellSetContext(A, &userJ);
1080: } else {
1081: A = J;
1082: }
1084: nullSpace = NULL;
1085: if (user.bcType != DIRICHLET) {
1086: MatNullSpaceCreate(PetscObjectComm((PetscObject) dm), PETSC_TRUE, 0, NULL, &nullSpace);
1087: MatSetNullSpace(A, nullSpace);
1088: }
1090: DMPlexSetSNESLocalFEM(dm,&user,&user,&user);
1091: SNESSetJacobian(snes, A, J, NULL, NULL);
1093: SNESSetFromOptions(snes);
1095: if (user.fieldBC) {DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u);}
1096: else {DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);}
1097: if (user.restart) {
1098: #if defined(PETSC_HAVE_HDF5)
1099: PetscViewer viewer;
1101: PetscViewerCreate(PETSC_COMM_WORLD, &viewer);
1102: PetscViewerSetType(viewer, PETSCVIEWERHDF5);
1103: PetscViewerFileSetMode(viewer, FILE_MODE_READ);
1104: PetscViewerFileSetName(viewer, user.filename);
1105: PetscViewerHDF5PushGroup(viewer, "/fields");
1106: VecLoad(u, viewer);
1107: PetscViewerHDF5PopGroup(viewer);
1108: PetscViewerDestroy(&viewer);
1109: #endif
1110: }
1111: if (user.showInitial) {
1112: Vec lv;
1113: DMGetLocalVector(dm, &lv);
1114: DMGlobalToLocalBegin(dm, u, INSERT_VALUES, lv);
1115: DMGlobalToLocalEnd(dm, u, INSERT_VALUES, lv);
1116: DMPrintLocalVec(dm, "Local function", 1.0e-10, lv);
1117: DMRestoreLocalVector(dm, &lv);
1118: }
1119: if (user.viewHierarchy) {
1120: SNES lsnes;
1121: KSP ksp;
1122: PC pc;
1123: PetscInt numLevels, l;
1124: PetscBool isMG;
1126: PetscObjectTypeCompare((PetscObject) snes, SNESFAS, &isFAS);
1127: if (isFAS) {
1128: SNESFASGetLevels(snes, &numLevels);
1129: for (l = 0; l < numLevels; ++l) {
1130: SNESFASGetCycleSNES(snes, l, &lsnes);
1131: SNESMonitorSet(lsnes, SNESMonitorError, &user, NULL);
1132: }
1133: } else {
1134: SNESGetKSP(snes, &ksp);
1135: KSPGetPC(ksp, &pc);
1136: PetscObjectTypeCompare((PetscObject) pc, PCMG, &isMG);
1137: if (isMG) {
1138: user.pcmg = pc;
1139: PCMGGetLevels(pc, &numLevels);
1140: for (l = 0; l < numLevels; ++l) {
1141: PCMGGetSmootherDown(pc, l, &ksp);
1142: KSPMonitorSet(ksp, MonitorError, &user, NULL);
1143: }
1144: }
1145: }
1146: }
1147: if (user.runType == RUN_FULL || user.runType == RUN_EXACT) {
1148: PetscErrorCode (*initialGuess[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {zero};
1150: if (user.nonzInit) initialGuess[0] = ecks;
1151: if (user.runType == RUN_FULL) {
1152: DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
1153: }
1154: if (user.debug) {
1155: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1156: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1157: }
1158: VecViewFromOptions(u, NULL, "-guess_vec_view");
1159: SNESSolve(snes, NULL, u);
1160: SNESGetSolution(snes, &u);
1161: SNESGetDM(snes, &dm);
1163: if (user.showSolution) {
1164: PetscPrintf(PETSC_COMM_WORLD, "Solution\n");
1165: VecChop(u, 3.0e-9);
1166: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1167: }
1168: } else if (user.runType == RUN_PERF) {
1169: Vec r;
1170: PetscReal res = 0.0;
1172: SNESGetFunction(snes, &r, NULL, NULL);
1173: SNESComputeFunction(snes, u, r);
1174: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1175: VecChop(r, 1.0e-10);
1176: VecNorm(r, NORM_2, &res);
1177: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res);
1178: } else {
1179: Vec r;
1180: PetscReal res = 0.0, tol = 1.0e-11;
1182: /* Check discretization error */
1183: SNESGetFunction(snes, &r, NULL, NULL);
1184: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1185: if (!user.quiet) {VecView(u, PETSC_VIEWER_STDOUT_WORLD);}
1186: DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
1187: if (error < tol) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < %2.1e\n", (double)tol);}
1188: else {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", (double)error);}
1189: /* Check residual */
1190: SNESComputeFunction(snes, u, r);
1191: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1192: VecChop(r, 1.0e-10);
1193: if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1194: VecNorm(r, NORM_2, &res);
1195: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res);
1196: /* Check Jacobian */
1197: {
1198: Vec b;
1200: SNESComputeJacobian(snes, u, A, A);
1201: VecDuplicate(u, &b);
1202: VecSet(r, 0.0);
1203: SNESComputeFunction(snes, r, b);
1204: MatMult(A, u, r);
1205: VecAXPY(r, 1.0, b);
1206: PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n");
1207: VecChop(r, 1.0e-10);
1208: if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1209: VecNorm(r, NORM_2, &res);
1210: PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", (double)res);
1211: /* check solver */
1212: if (user.checkksp) {
1213: KSP ksp;
1215: if (nullSpace) {
1216: MatNullSpaceRemove(nullSpace, u);
1217: }
1218: SNESComputeJacobian(snes, u, A, J);
1219: MatMult(A, u, b);
1220: SNESGetKSP(snes, &ksp);
1221: KSPSetOperators(ksp, A, J);
1222: KSPSolve(ksp, b, r);
1223: VecAXPY(r, -1.0, u);
1224: VecNorm(r, NORM_2, &res);
1225: PetscPrintf(PETSC_COMM_WORLD, "KSP Error: %g\n", (double)res);
1226: }
1227: VecDestroy(&b);
1228: }
1229: }
1230: VecViewFromOptions(u, NULL, "-vec_view");
1231: {
1232: Vec nu;
1234: DMGetAuxiliaryVec(dm, NULL, 0, &nu);
1235: if (nu) {VecViewFromOptions(nu, NULL, "-coeff_view");}
1236: }
1238: if (user.bdIntegral) {
1239: DMLabel label;
1240: PetscInt id = 1;
1241: PetscScalar bdInt = 0.0;
1242: PetscReal exact = 3.3333333333;
1244: DMGetLabel(dm, "marker", &label);
1245: DMPlexComputeBdIntegral(dm, u, label, 1, &id, bd_integral_2d, &bdInt, NULL);
1246: PetscPrintf(PETSC_COMM_WORLD, "Solution boundary integral: %.4g\n", (double) PetscAbsScalar(bdInt));
1247: if (PetscAbsReal(PetscAbsScalar(bdInt) - exact) > PETSC_SQRT_MACHINE_EPSILON) SETERRQ2(PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Invalid boundary integral %g != %g", (double) PetscAbsScalar(bdInt), (double)exact);
1248: }
1250: MatNullSpaceDestroy(&nullSpace);
1251: if (user.jacobianMF) {VecDestroy(&userJ.u);}
1252: if (A != J) {MatDestroy(&A);}
1253: MatDestroy(&J);
1254: VecDestroy(&u);
1255: SNESDestroy(&snes);
1256: DMDestroy(&dm);
1257: PetscFree2(user.exactFuncs, user.exactFields);
1258: PetscFree(user.kgrid);
1259: PetscFinalize();
1260: return ierr;
1261: }
1263: /*TEST
1264: # 2D serial P1 test 0-4
1265: test:
1266: suffix: 2d_p1_0
1267: requires: triangle
1268: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1270: test:
1271: suffix: 2d_p1_1
1272: requires: triangle
1273: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1275: test:
1276: suffix: 2d_p1_2
1277: requires: triangle
1278: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1280: test:
1281: suffix: 2d_p1_neumann_0
1282: requires: triangle
1283: 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
1285: test:
1286: suffix: 2d_p1_neumann_1
1287: requires: triangle
1288: args: -run_type test -refinement_limit 0.0625 -bc_type neumann -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1290: # 2D serial P2 test 5-8
1291: test:
1292: suffix: 2d_p2_0
1293: requires: triangle
1294: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1296: test:
1297: suffix: 2d_p2_1
1298: requires: triangle
1299: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1301: test:
1302: suffix: 2d_p2_neumann_0
1303: requires: triangle
1304: 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
1306: test:
1307: suffix: 2d_p2_neumann_1
1308: requires: triangle
1309: 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
1311: test:
1312: suffix: bd_int_0
1313: requires: triangle
1314: args: -run_type test -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet
1316: test:
1317: suffix: bd_int_1
1318: requires: triangle
1319: args: -run_type test -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet
1321: # 3D serial P1 test 9-12
1322: test:
1323: suffix: 3d_p1_0
1324: requires: ctetgen
1325: 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
1327: test:
1328: suffix: 3d_p1_1
1329: requires: ctetgen
1330: 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
1332: test:
1333: suffix: 3d_p1_2
1334: requires: ctetgen
1335: 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
1337: test:
1338: suffix: 3d_p1_neumann_0
1339: requires: ctetgen
1340: 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
1342: # Analytic variable coefficient 13-20
1343: test:
1344: suffix: 13
1345: requires: triangle
1346: args: -run_type test -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1347: test:
1348: suffix: 14
1349: requires: triangle
1350: args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1351: test:
1352: suffix: 15
1353: requires: triangle
1354: args: -run_type test -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1355: test:
1356: suffix: 16
1357: requires: triangle
1358: args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1359: test:
1360: suffix: 17
1361: requires: ctetgen
1362: 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
1364: test:
1365: suffix: 18
1366: requires: ctetgen
1367: 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
1369: test:
1370: suffix: 19
1371: requires: ctetgen
1372: 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
1374: test:
1375: suffix: 20
1376: requires: ctetgen
1377: 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
1379: # P1 variable coefficient 21-28
1380: test:
1381: suffix: 21
1382: requires: triangle
1383: 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
1385: test:
1386: suffix: 22
1387: requires: triangle
1388: 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
1390: test:
1391: suffix: 23
1392: requires: triangle
1393: 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
1395: test:
1396: suffix: 24
1397: requires: triangle
1398: 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
1400: test:
1401: suffix: 25
1402: requires: ctetgen
1403: 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
1405: test:
1406: suffix: 26
1407: requires: ctetgen
1408: 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
1410: test:
1411: suffix: 27
1412: requires: ctetgen
1413: 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
1415: test:
1416: suffix: 28
1417: requires: ctetgen
1418: 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
1420: # P0 variable coefficient 29-36
1421: test:
1422: suffix: 29
1423: requires: triangle
1424: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1426: test:
1427: suffix: 30
1428: requires: triangle
1429: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1431: test:
1432: suffix: 31
1433: requires: triangle
1434: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1436: test:
1437: requires: triangle
1438: suffix: 32
1439: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1441: test:
1442: requires: ctetgen
1443: suffix: 33
1444: 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
1446: test:
1447: suffix: 34
1448: requires: ctetgen
1449: 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
1451: test:
1452: suffix: 35
1453: requires: ctetgen
1454: 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
1456: test:
1457: suffix: 36
1458: requires: ctetgen
1459: 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
1461: # Full solve 39-44
1462: test:
1463: suffix: 39
1464: requires: triangle !single
1465: 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
1466: test:
1467: suffix: 40
1468: requires: triangle !single
1469: 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
1470: test:
1471: suffix: 41
1472: requires: triangle !single
1473: args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 2 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -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
1474: test:
1475: suffix: 42
1476: requires: triangle !single
1477: args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -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
1478: test:
1479: suffix: 43
1480: requires: triangle !single
1481: nsize: 2
1482: args: -run_type full -refinement_limit 0.03125 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 2 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -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
1484: test:
1485: suffix: 44
1486: requires: triangle !single
1487: nsize: 2
1488: args: -run_type full -refinement_limit 0.0625 -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -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
1490: # These tests use a loose tolerance just to exercise the PtAP operations for MATIS and multiple PCBDDC setup calls inside PCMG
1491: testset:
1492: requires: triangle !single
1493: nsize: 3
1494: 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
1495: test:
1496: suffix: gmg_bddc
1497: filter: sed -e "s/CONVERGED_FNORM_RELATIVE iterations 3/CONVERGED_FNORM_RELATIVE iterations 4/g"
1498: args: -mg_levels_pc_type jacobi
1499: test:
1500: filter: sed -e "s/iterations [0-4]/iterations 4/g"
1501: suffix: gmg_bddc_lev
1502: args: -mg_levels_pc_type bddc
1504: # Restarting
1505: testset:
1506: suffix: restart
1507: requires: hdf5 triangle !complex
1508: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1
1509: test:
1510: args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1511: test:
1512: args: -f sol.h5 -restart
1514: # Periodicity
1515: test:
1516: suffix: periodic_0
1517: requires: triangle
1518: args: -run_type full -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -snes_converged_reason ::ascii_info_detail
1520: test:
1521: requires: !complex
1522: suffix: periodic_1
1523: 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
1525: # 2D serial P1 test with field bc
1526: test:
1527: suffix: field_bc_2d_p1_0
1528: requires: triangle
1529: 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
1531: test:
1532: suffix: field_bc_2d_p1_1
1533: requires: triangle
1534: 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
1536: test:
1537: suffix: field_bc_2d_p1_neumann_0
1538: requires: triangle
1539: 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
1541: test:
1542: suffix: field_bc_2d_p1_neumann_1
1543: requires: triangle
1544: 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
1546: # 3D serial P1 test with field bc
1547: test:
1548: suffix: field_bc_3d_p1_0
1549: requires: ctetgen
1550: 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
1552: test:
1553: suffix: field_bc_3d_p1_1
1554: requires: ctetgen
1555: 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
1557: test:
1558: suffix: field_bc_3d_p1_neumann_0
1559: requires: ctetgen
1560: 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
1562: test:
1563: suffix: field_bc_3d_p1_neumann_1
1564: requires: ctetgen
1565: 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
1567: # 2D serial P2 test with field bc
1568: test:
1569: suffix: field_bc_2d_p2_0
1570: requires: triangle
1571: 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
1573: test:
1574: suffix: field_bc_2d_p2_1
1575: requires: triangle
1576: 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
1578: test:
1579: suffix: field_bc_2d_p2_neumann_0
1580: requires: triangle
1581: 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
1583: test:
1584: suffix: field_bc_2d_p2_neumann_1
1585: requires: triangle
1586: 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
1588: # 3D serial P2 test with field bc
1589: test:
1590: suffix: field_bc_3d_p2_0
1591: requires: ctetgen
1592: 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
1594: test:
1595: suffix: field_bc_3d_p2_1
1596: requires: ctetgen
1597: 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
1599: test:
1600: suffix: field_bc_3d_p2_neumann_0
1601: requires: ctetgen
1602: 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
1604: test:
1605: suffix: field_bc_3d_p2_neumann_1
1606: requires: ctetgen
1607: 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
1609: # Full solve simplex: Convergence
1610: test:
1611: suffix: 3d_p1_conv
1612: requires: ctetgen
1613: args: -run_type full -dim 3 -cells 1,1,1 -dm_refine 1 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 \
1614: -snes_convergence_estimate -convest_num_refine 1 -pc_type lu
1616: # Full solve simplex: PCBDDC
1617: test:
1618: suffix: tri_bddc
1619: requires: triangle !single
1620: nsize: 5
1621: 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
1623: # Full solve simplex: PCBDDC
1624: test:
1625: suffix: tri_parmetis_bddc
1626: requires: triangle !single parmetis
1627: nsize: 4
1628: 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
1630: testset:
1631: 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
1632: nsize: 5
1633: output_file: output/ex12_quad_bddc.out
1634: filter: sed -e "s/aijcusparse/aij/g" -e "s/aijviennacl/aij/g" -e "s/factorization: cusparse/factorization: petsc/g"
1635: test:
1636: requires: !single
1637: suffix: quad_bddc
1638: test:
1639: requires: !single cuda
1640: suffix: quad_bddc_cuda
1641: args: -matis_localmat_type aijcusparse -pc_bddc_dirichlet_pc_factor_mat_solver_type cusparse -pc_bddc_neumann_pc_factor_mat_solver_type cusparse
1642: test:
1643: requires: !single viennacl
1644: suffix: quad_bddc_viennacl
1645: args: -matis_localmat_type aijviennacl
1647: # Full solve simplex: ASM
1648: test:
1649: suffix: tri_q2q1_asm_lu
1650: requires: triangle !single
1651: 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
1653: test:
1654: suffix: tri_q2q1_msm_lu
1655: requires: triangle !single
1656: 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
1658: test:
1659: suffix: tri_q2q1_asm_sor
1660: requires: triangle !single
1661: 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
1663: test:
1664: suffix: tri_q2q1_msm_sor
1665: requires: triangle !single
1666: 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
1668: # Full solve simplex: FAS
1669: test:
1670: suffix: fas_newton_0
1671: requires: triangle !single
1672: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 2 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -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
1674: test:
1675: suffix: fas_newton_1
1676: requires: triangle !single
1677: args: -run_type full -dm_refine_hierarchy 3 -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type lu -fas_coarse_snes_monitor_short -snes_monitor_short -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
1678: filter: sed -e "s/total number of linear solver iterations=14/total number of linear solver iterations=15/g"
1680: test:
1681: suffix: fas_ngs_0
1682: requires: triangle !single
1683: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 2 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -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
1685: test:
1686: suffix: fas_newton_coarse_0
1687: requires: pragmatic triangle
1688: TODO: broken
1689: 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 -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
1691: test:
1692: suffix: mg_newton_coarse_0
1693: requires: triangle pragmatic
1694: TODO: broken
1695: args: -run_type full -dm_refine 3 -interpolate 1 -petscspace_degree 1 -snes_monitor_short -ksp_monitor_true_residual -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
1697: test:
1698: suffix: mg_newton_coarse_1
1699: requires: triangle pragmatic
1700: TODO: broken
1701: 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
1703: test:
1704: suffix: mg_newton_coarse_2
1705: requires: triangle pragmatic
1706: TODO: broken
1707: 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
1709: # Full solve tensor
1710: test:
1711: suffix: tensor_plex_2d
1712: args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_refine_hierarchy 2 -cells 2,2
1714: test:
1715: suffix: tensor_p4est_2d
1716: requires: p4est
1717: 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
1719: test:
1720: suffix: tensor_plex_3d
1721: 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
1723: test:
1724: suffix: tensor_p4est_3d
1725: requires: p4est
1726: 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
1728: test:
1729: suffix: p4est_test_q2_conformal_serial
1730: requires: p4est
1731: 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
1733: test:
1734: suffix: p4est_test_q2_conformal_parallel
1735: requires: p4est
1736: nsize: 7
1737: 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
1739: test:
1740: suffix: p4est_test_q2_conformal_parallel_parmetis
1741: requires: parmetis p4est
1742: nsize: 4
1743: 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
1745: test:
1746: suffix: p4est_test_q2_nonconformal_serial
1747: requires: p4est
1748: filter: grep -v "CG or CGNE: variant"
1749: 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
1751: test:
1752: suffix: p4est_test_q2_nonconformal_parallel
1753: requires: p4est
1754: filter: grep -v "CG or CGNE: variant"
1755: nsize: 7
1756: 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
1758: test:
1759: suffix: p4est_test_q2_nonconformal_parallel_parmetis
1760: requires: parmetis p4est
1761: nsize: 4
1762: 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
1764: test:
1765: suffix: p4est_exact_q2_conformal_serial
1766: requires: p4est !single !complex !__float128
1767: args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -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
1769: test:
1770: suffix: p4est_exact_q2_conformal_parallel
1771: requires: p4est !single !complex !__float128
1772: nsize: 4
1773: args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -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
1775: test:
1776: suffix: p4est_exact_q2_conformal_parallel_parmetis
1777: requires: parmetis p4est !single
1778: nsize: 4
1779: args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -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
1781: test:
1782: suffix: p4est_exact_q2_nonconformal_serial
1783: requires: p4est
1784: args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -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
1786: test:
1787: suffix: p4est_exact_q2_nonconformal_parallel
1788: requires: p4est
1789: nsize: 7
1790: args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -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
1792: test:
1793: suffix: p4est_exact_q2_nonconformal_parallel_parmetis
1794: requires: parmetis p4est
1795: nsize: 4
1796: args: -run_type exact -interpolate 1 -petscspace_degree 2 -fas_levels_snes_atol 1.e-10 -snes_max_it 1 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type none -fas_coarse_ksp_type preonly -fas_coarse_snes_monitor_short -snes_monitor_short -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
1798: test:
1799: suffix: p4est_full_q2_nonconformal_serial
1800: requires: p4est !single
1801: filter: grep -v "variant HERMITIAN"
1802: args: -run_type full -interpolate 1 -petscspace_degree 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type jacobi -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -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
1804: test:
1805: suffix: p4est_full_q2_nonconformal_parallel
1806: requires: p4est !single
1807: filter: grep -v "variant HERMITIAN"
1808: nsize: 7
1809: args: -run_type full -interpolate 1 -petscspace_degree 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -fas_coarse_pc_type jacobi -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -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
1811: test:
1812: suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1813: requires: p4est !single
1814: filter: grep -v "variant HERMITIAN"
1815: nsize: 7
1816: args: -run_type full -interpolate 1 -petscspace_degree 2 -snes_max_it 20 -snes_type fas -snes_fas_levels 3 -dm_mat_type is -fas_coarse_pc_type bddc -fas_coarse_ksp_type cg -fas_coarse_snes_monitor_short -snes_monitor_short -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
1818: test:
1819: suffix: p4est_full_q2_nonconformal_parallel_bddc
1820: requires: p4est !single
1821: filter: grep -v "variant HERMITIAN"
1822: nsize: 7
1823: 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
1825: test:
1826: TODO: broken
1827: suffix: p4est_fas_q2_conformal_serial
1828: requires: p4est !complex !__float128
1829: 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 -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
1831: test:
1832: TODO: broken
1833: suffix: p4est_fas_q2_nonconformal_serial
1834: requires: p4est
1835: 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 -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
1837: test:
1838: suffix: fas_newton_0_p4est
1839: requires: p4est !single !__float128
1840: args: -run_type full -variable_coefficient nonlinear -interpolate 1 -petscspace_degree 1 -snes_type fas -snes_fas_levels 2 -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -snes_monitor_short -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
1842: # Full solve simplicial AMR
1843: test:
1844: suffix: tri_p1_adapt_0
1845: requires: pragmatic
1846: TODO: broken
1847: 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 1
1849: test:
1850: suffix: tri_p1_adapt_1
1851: requires: pragmatic
1852: TODO: broken
1853: 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
1855: test:
1856: suffix: tri_p1_adapt_analytic_0
1857: requires: pragmatic
1858: TODO: broken
1859: 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
1861: # Full solve tensor AMR
1862: test:
1863: suffix: quad_q1_adapt_0
1864: requires: p4est
1865: 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 1 -dm_view
1866: filter: grep -v DM_
1868: test:
1869: suffix: amr_0
1870: nsize: 5
1871: 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
1873: test:
1874: suffix: amr_1
1875: requires: p4est !complex
1876: 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
1878: test:
1879: suffix: p4est_solve_bddc
1880: requires: p4est !complex
1881: 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
1882: nsize: 4
1884: test:
1885: suffix: p4est_solve_fas
1886: requires: p4est
1887: 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
1888: nsize: 4
1889: TODO: identical machine two runs produce slightly different solver trackers
1891: test:
1892: suffix: p4est_convergence_test_1
1893: requires: p4est
1894: 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
1895: nsize: 4
1897: test:
1898: suffix: p4est_convergence_test_2
1899: requires: p4est
1900: 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
1902: test:
1903: suffix: p4est_convergence_test_3
1904: requires: p4est
1905: 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
1907: test:
1908: suffix: p4est_convergence_test_4
1909: requires: p4est
1910: 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
1911: timeoutfactor: 5
1913: # Serial tests with GLVis visualization
1914: test:
1915: suffix: glvis_2d_tet_p1
1916: 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 -dm_plex_gmsh_periodic 0
1917: test:
1918: suffix: glvis_2d_tet_p2
1919: 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 -dm_plex_gmsh_periodic 0
1920: test:
1921: suffix: glvis_2d_hex_p1
1922: args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -simplex 0 -dm_refine 1
1923: test:
1924: suffix: glvis_2d_hex_p2
1925: args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -simplex 0 -dm_refine 1
1926: test:
1927: suffix: glvis_2d_hex_p2_p4est
1928: requires: p4est
1929: 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
1930: test:
1931: suffix: glvis_2d_tet_p0
1932: args: -run_type exact -interpolate 1 -guess_vec_view glvis: -nonzero_initial_guess 1 -f ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -petscspace_degree 0
1933: test:
1934: suffix: glvis_2d_hex_p0
1935: args: -run_type exact -interpolate 1 -guess_vec_view glvis: -nonzero_initial_guess 1 -cells 5,7 -simplex 0 -petscspace_degree 0
1937: # PCHPDDM tests
1938: testset:
1939: nsize: 4
1940: requires: hpddm slepc !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1941: args: -run_type test -run_test_check_ksp -quiet -petscspace_degree 1 -interpolate 1 -petscpartitioner_type simple -bc_type none -simplex 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 2 -pc_hpddm_coarse_p 1 -pc_hpddm_coarse_pc_type svd -ksp_rtol 1.e-10 -pc_hpddm_levels_1_st_pc_factor_shift_type INBLOCKS -ksp_converged_reason
1942: test:
1943: suffix: quad_singular_hpddm
1944: args: -cells 6,7
1945: test:
1946: requires: p4est
1947: suffix: p4est_singular_2d_hpddm
1948: args: -dm_plex_convert_type p4est -dm_forest_minimum_refinement 1 -dm_forest_initial_refinement 3 -dm_forest_maximum_refinement 3
1949: test:
1950: requires: p4est
1951: suffix: p4est_nc_singular_2d_hpddm
1952: args: -dm_plex_convert_type p4est -dm_forest_minimum_refinement 1 -dm_forest_initial_refinement 1 -dm_forest_maximum_refinement 3 -dm_p4est_refine_pattern hash
1953: testset:
1954: nsize: 4
1955: requires: hpddm slepc triangle !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1956: args: -run_type full -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 4 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-1
1957: test:
1958: args: -pc_hpddm_coarse_mat_type baij -options_left no
1959: suffix: tri_hpddm_reuse_baij
1960: test:
1961: requires: !complex
1962: suffix: tri_hpddm_reuse
1963: testset:
1964: nsize: 4
1965: requires: hpddm slepc !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1966: args: -run_type full -petscpartitioner_type simple -cells 7,5 -dm_refine 2 -simplex 0 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_nev 4 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-1
1967: test:
1968: args: -pc_hpddm_coarse_mat_type baij -options_left no
1969: suffix: quad_hpddm_reuse_baij
1970: test:
1971: requires: !complex
1972: suffix: quad_hpddm_reuse
1973: testset:
1974: nsize: 4
1975: requires: hpddm slepc !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1976: args: -run_type full -petscpartitioner_type simple -cells 7,5 -dm_refine 2 -simplex 0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type lu -pc_hpddm_levels_1_eps_threshold 0.1 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-1
1977: test:
1978: args: -pc_hpddm_coarse_mat_type baij -options_left no
1979: suffix: quad_hpddm_reuse_threshold_baij
1980: test:
1981: requires: !complex
1982: suffix: quad_hpddm_reuse_threshold
1983: testset:
1984: nsize: 4
1985: requires: hpddm slepc parmetis !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1986: filter: sed -e "s/linear solver iterations=17/linear solver iterations=16/g"
1987: args: -run_type full -petscpartitioner_type parmetis -dm_refine 3 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -ksp_type gmres -ksp_gmres_restart 100 -pc_type hpddm -snes_monitor_short -snes_converged_reason ::ascii_info_detail -snes_view -show_solution 0 -pc_type hpddm -pc_hpddm_levels_1_sub_pc_type icc -pc_hpddm_levels_1_eps_nev 20 -pc_hpddm_coarse_p 2 -pc_hpddm_coarse_pc_type redundant -ksp_rtol 1.e-10 -f ${PETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -pc_hpddm_levels_1_sub_pc_factor_levels 3 -variable_coefficient circle -dm_plex_gmsh_periodic 0
1988: test:
1989: args: -pc_hpddm_coarse_mat_type baij -options_left no
1990: filter: grep -v " total: nonzeros=" | grep -v " rows=" | sed -e "s/total number of linear solver iterations=1[5-7]/total number of linear solver iterations=16/g"
1991: suffix: tri_parmetis_hpddm_baij
1992: test:
1993: filter: grep -v " total: nonzeros=" | grep -v " rows=" | sed -e "s/total number of linear solver iterations=1[5-7]/total number of linear solver iterations=16/g"
1994: requires: !complex
1995: suffix: tri_parmetis_hpddm
1997: # 2D serial P1 tests for adaptive MG
1998: test:
1999: suffix: 2d_p1_adaptmg_0
2000: requires: triangle bamg
2001: args: -dm_refine_hierarchy 3 -cells 4,4 -interpolate -bc_type dirichlet -petscspace_degree 1 \
2002: -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
2003: -snes_max_it 1 -ksp_converged_reason \
2004: -ksp_rtol 1e-8 -pc_type mg
2005: # -ksp_monitor_true_residual -ksp_converged_reason -mg_levels_ksp_monitor_true_residual -pc_mg_mesp_monitor -dm_adapt_interp_view_fine draw -dm_adapt_interp_view_coarse draw -draw_pause 1
2006: test:
2007: suffix: 2d_p1_adaptmg_1
2008: requires: triangle bamg
2009: args: -dm_refine_hierarchy 3 -cells 4,4 -interpolate -bc_type dirichlet -petscspace_degree 1 \
2010: -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
2011: -snes_max_it 1 -ksp_converged_reason \
2012: -ksp_rtol 1e-8 -pc_type mg -pc_mg_galerkin -pc_mg_adapt_interp -pc_mg_adapt_interp_coarse_space eigenvector -pc_mg_adapt_interp_n 1 \
2013: -pc_mg_mesp_ksp_type richardson -pc_mg_mesp_ksp_richardson_self_scale -pc_mg_mesp_ksp_max_it 100 -pc_mg_mesp_pc_type none
2015: TEST*/