Actual source code: ex12.c
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: static void f1_bd_zero(PetscInt dim, PetscInt Nf, PetscInt NfAux,
182: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
183: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
184: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
185: {
186: PetscInt comp;
187: for (comp = 0; comp < dim; ++comp) f1[comp] = 0.0;
188: }
190: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
191: static void f1_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
192: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
193: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
194: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
195: {
196: PetscInt d;
197: for (d = 0; d < dim; ++d) f1[d] = u_x[d];
198: }
200: /* < \nabla v, \nabla u + {\nabla u}^T >
201: This just gives \nabla u, give the perdiagonal for the transpose */
202: static void g3_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
203: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
204: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
205: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
206: {
207: PetscInt d;
208: for (d = 0; d < dim; ++d) g3[d*dim+d] = 1.0;
209: }
211: /*
212: In 2D for x periodicity and y Dirichlet conditions, we use exact solution:
214: u = sin(2 pi x)
215: f = -4 pi^2 sin(2 pi x)
217: so that
219: -\Delta u + f = 4 pi^2 sin(2 pi x) - 4 pi^2 sin(2 pi x) = 0
220: */
221: static PetscErrorCode xtrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
222: {
223: *u = PetscSinReal(2.0*PETSC_PI*x[0]);
224: return 0;
225: }
227: static void f0_xtrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
228: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
229: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
230: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
231: {
232: f0[0] = -4.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
233: }
235: /*
236: In 2D for x-y periodicity, we use exact solution:
238: u = sin(2 pi x) sin(2 pi y)
239: f = -8 pi^2 sin(2 pi x)
241: so that
243: -\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
244: */
245: static PetscErrorCode xytrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
246: {
247: *u = PetscSinReal(2.0*PETSC_PI*x[0])*PetscSinReal(2.0*PETSC_PI*x[1]);
248: return 0;
249: }
251: static void f0_xytrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
252: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
253: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
254: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
255: {
256: f0[0] = -8.0*PetscSqr(PETSC_PI)*PetscSinReal(2.0*PETSC_PI*x[0]);
257: }
259: /*
260: In 2D for Dirichlet conditions with a variable coefficient, we use exact solution:
262: u = x^2 + y^2
263: f = 6 (x + y)
264: nu = (x + y)
266: so that
268: -\div \nu \grad u + f = -6 (x + y) + 6 (x + y) = 0
269: */
270: static PetscErrorCode nu_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
271: {
272: *u = x[0] + x[1];
273: return 0;
274: }
276: static PetscErrorCode checkerboardCoeff(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
277: {
278: AppCtx *user = (AppCtx *) ctx;
279: PetscInt div = user->div;
280: PetscInt k = user->k;
281: PetscInt mask = 0, ind = 0, d;
284: for (d = 0; d < dim; ++d) mask = (mask + (PetscInt) (x[d]*div)) % 2;
285: if (user->kgrid) {
286: for (d = 0; d < dim; ++d) {
287: if (d > 0) ind *= dim;
288: ind += (PetscInt) (x[d]*div);
289: }
290: k = user->kgrid[ind];
291: }
292: u[0] = mask ? 1.0 : PetscPowRealInt(10.0, -k);
293: return(0);
294: }
296: void f0_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 f0[])
300: {
301: f0[0] = 6.0*(x[0] + x[1]);
302: }
304: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
305: void f1_analytic_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] = (x[0] + x[1])*u_x[d];
312: }
314: void f1_field_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
315: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
316: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
317: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
318: {
319: PetscInt d;
320: for (d = 0; d < dim; ++d) f1[d] = a[0]*u_x[d];
321: }
323: /* < \nabla v, \nabla u + {\nabla u}^T >
324: This just gives \nabla u, give the perdiagonal for the transpose */
325: void g3_analytic_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] = x[0] + x[1];
332: }
334: void g3_field_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
335: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
336: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
337: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
338: {
339: PetscInt d;
340: for (d = 0; d < dim; ++d) g3[d*dim+d] = a[0];
341: }
343: /*
344: In 2D for Dirichlet conditions with a nonlinear coefficient (p-Laplacian with p = 4), we use exact solution:
346: u = x^2 + y^2
347: f = 16 (x^2 + y^2)
348: nu = 1/2 |grad u|^2
350: so that
352: -\div \nu \grad u + f = -16 (x^2 + y^2) + 16 (x^2 + y^2) = 0
353: */
354: void f0_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 f0[])
358: {
359: f0[0] = 16.0*(x[0]*x[0] + x[1]*x[1]);
360: }
362: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
363: void f1_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
364: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
365: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
366: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
367: {
368: PetscScalar nu = 0.0;
369: PetscInt d;
370: for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
371: for (d = 0; d < dim; ++d) f1[d] = 0.5*nu*u_x[d];
372: }
374: /*
375: grad (u + eps w) - grad u = eps grad w
377: 1/2 |grad (u + eps w)|^2 grad (u + eps w) - 1/2 |grad u|^2 grad u
378: = 1/2 (|grad u|^2 + 2 eps <grad u,grad w>) (grad u + eps grad w) - 1/2 |grad u|^2 grad u
379: = 1/2 (eps |grad u|^2 grad w + 2 eps <grad u,grad w> grad u)
380: = eps (1/2 |grad u|^2 grad w + grad u <grad u,grad w>)
381: */
382: void g3_analytic_nonlinear_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
383: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
384: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
385: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
386: {
387: PetscScalar nu = 0.0;
388: PetscInt d, e;
389: for (d = 0; d < dim; ++d) nu += u_x[d]*u_x[d];
390: for (d = 0; d < dim; ++d) {
391: g3[d*dim+d] = 0.5*nu;
392: for (e = 0; e < dim; ++e) {
393: g3[d*dim+e] += u_x[d]*u_x[e];
394: }
395: }
396: }
398: /*
399: In 3D for Dirichlet conditions we use exact solution:
401: u = 2/3 (x^2 + y^2 + z^2)
402: f = 4
404: so that
406: -\Delta u + f = -2/3 * 6 + 4 = 0
408: For Neumann conditions, we have
410: -\nabla u \cdot -\hat z |_{z=0} = (2z)|_{z=0} = 0 (bottom)
411: -\nabla u \cdot \hat z |_{z=1} = -(2z)|_{z=1} = -2 (top)
412: -\nabla u \cdot -\hat y |_{y=0} = (2y)|_{y=0} = 0 (front)
413: -\nabla u \cdot \hat y |_{y=1} = -(2y)|_{y=1} = -2 (back)
414: -\nabla u \cdot -\hat x |_{x=0} = (2x)|_{x=0} = 0 (left)
415: -\nabla u \cdot \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)
417: Which we can express as
419: \nabla u \cdot \hat n|_\Gamma = {2 x, 2 y, 2z} \cdot \hat n = 2 (x + y + z)
420: */
421: static PetscErrorCode quadratic_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
422: {
423: *u = 2.0*(x[0]*x[0] + x[1]*x[1] + x[2]*x[2])/3.0;
424: return 0;
425: }
427: static void quadratic_u_field_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
428: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
429: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
430: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
431: {
432: uexact[0] = a[0];
433: }
435: static void bd_integral_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux,
436: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
437: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
438: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar *uint)
439: {
440: uint[0] = u[0];
441: }
443: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
444: {
445: const char *bcTypes[3] = {"neumann", "dirichlet", "none"};
446: const char *runTypes[4] = {"full", "exact", "test", "perf"};
447: const char *coeffTypes[8] = {"none", "analytic", "field", "nonlinear", "circle", "cross", "checkerboard_0", "checkerboard_1"};
448: PetscInt bd, bc, run, coeff, n;
449: PetscBool rand = PETSC_FALSE, flg;
453: options->debug = 0;
454: options->runType = RUN_FULL;
455: options->dim = 2;
456: options->periodicity[0] = DM_BOUNDARY_NONE;
457: options->periodicity[1] = DM_BOUNDARY_NONE;
458: options->periodicity[2] = DM_BOUNDARY_NONE;
459: options->cells[0] = 2;
460: options->cells[1] = 2;
461: options->cells[2] = 2;
462: options->filename[0] = '\0';
463: options->interpolate = PETSC_TRUE;
464: options->refinementLimit = 0.0;
465: options->bcType = DIRICHLET;
466: options->variableCoefficient = COEFF_NONE;
467: options->fieldBC = PETSC_FALSE;
468: options->jacobianMF = PETSC_FALSE;
469: options->showInitial = PETSC_FALSE;
470: options->showSolution = PETSC_FALSE;
471: options->restart = PETSC_FALSE;
472: options->viewHierarchy = PETSC_FALSE;
473: options->simplex = PETSC_TRUE;
474: options->quiet = PETSC_FALSE;
475: options->nonzInit = PETSC_FALSE;
476: options->bdIntegral = PETSC_FALSE;
477: options->checkksp = PETSC_FALSE;
478: options->div = 4;
479: options->k = 1;
480: options->kgrid = NULL;
482: PetscOptionsBegin(comm, "", "Poisson Problem Options", "DMPLEX");
483: PetscOptionsInt("-debug", "The debugging level", "ex12.c", options->debug, &options->debug, NULL);
484: run = options->runType;
485: PetscOptionsEList("-run_type", "The run type", "ex12.c", runTypes, 4, runTypes[options->runType], &run, NULL);
487: options->runType = (RunType) run;
489: PetscOptionsInt("-dim", "The topological mesh dimension", "ex12.c", options->dim, &options->dim, NULL);
490: bd = options->periodicity[0];
491: PetscOptionsEList("-x_periodicity", "The x-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[0]], &bd, NULL);
492: options->periodicity[0] = (DMBoundaryType) bd;
493: bd = options->periodicity[1];
494: PetscOptionsEList("-y_periodicity", "The y-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[1]], &bd, NULL);
495: options->periodicity[1] = (DMBoundaryType) bd;
496: bd = options->periodicity[2];
497: PetscOptionsEList("-z_periodicity", "The z-boundary periodicity", "ex12.c", DMBoundaryTypes, 5, DMBoundaryTypes[options->periodicity[2]], &bd, NULL);
498: options->periodicity[2] = (DMBoundaryType) bd;
499: n = 3;
500: PetscOptionsIntArray("-cells", "The initial mesh division", "ex12.c", options->cells, &n, NULL);
501: PetscOptionsString("-f", "Mesh filename to read", "ex12.c", options->filename, options->filename, sizeof(options->filename), &flg);
502: PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex12.c", options->interpolate, &options->interpolate, NULL);
503: PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex12.c", options->refinementLimit, &options->refinementLimit, NULL);
504: bc = options->bcType;
505: PetscOptionsEList("-bc_type","Type of boundary condition","ex12.c",bcTypes,3,bcTypes[options->bcType],&bc,NULL);
506: options->bcType = (BCType) bc;
507: coeff = options->variableCoefficient;
508: PetscOptionsEList("-variable_coefficient","Type of variable coefficent","ex12.c",coeffTypes,8,coeffTypes[options->variableCoefficient],&coeff,NULL);
509: options->variableCoefficient = (CoeffType) coeff;
511: PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL);
512: PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL);
513: PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL);
514: PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL);
515: PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL);
516: PetscOptionsBool("-dm_view_hierarchy", "View the coarsened hierarchy", "ex12.c", options->viewHierarchy, &options->viewHierarchy, NULL);
517: PetscOptionsBool("-simplex", "Simplicial (true) or tensor (false) mesh", "ex12.c", options->simplex, &options->simplex, NULL);
518: PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL);
519: PetscOptionsBool("-nonzero_initial_guess", "nonzero initial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL);
520: PetscOptionsBool("-bd_integral", "Compute the integral of the solution on the boundary", "ex12.c", options->bdIntegral, &options->bdIntegral, NULL);
521: if (options->runType == RUN_TEST) {
522: PetscOptionsBool("-run_test_check_ksp", "Check solution of KSP", "ex12.c", options->checkksp, &options->checkksp, NULL);
523: }
524: PetscOptionsInt("-div", "The number of division for the checkerboard coefficient", "ex12.c", options->div, &options->div, NULL);
525: PetscOptionsInt("-k", "The exponent for the checkerboard coefficient", "ex12.c", options->k, &options->k, NULL);
526: PetscOptionsBool("-k_random", "Assign random k values to checkerboard", "ex12.c", rand, &rand, NULL);
527: PetscOptionsEnd();
528: PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
530: if (rand) {
531: PetscRandom r;
532: PetscReal val;
533: PetscInt N = PetscPowInt(options->div, options->dim), i;
535: PetscMalloc1(N, &options->kgrid);
536: PetscRandomCreate(PETSC_COMM_SELF, &r);
537: PetscRandomSetFromOptions(r);
538: PetscRandomSetInterval(r, 0.0, options->k);
539: PetscRandomSetSeed(r, 1973);
540: PetscRandomSeed(r);
541: for (i = 0; i < N; ++i) {
542: PetscRandomGetValueReal(r, &val);
543: options->kgrid[i] = 1 + (PetscInt) val;
544: }
545: PetscRandomDestroy(&r);
546: }
547: return(0);
548: }
550: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
551: {
552: DM plex;
553: DMLabel label;
557: DMCreateLabel(dm, name);
558: DMGetLabel(dm, name, &label);
559: DMConvert(dm, DMPLEX, &plex);
560: DMPlexMarkBoundaryFaces(plex, 1, label);
561: DMDestroy(&plex);
562: return(0);
563: }
565: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
566: {
567: PetscInt dim = user->dim;
568: const char *filename = user->filename;
569: PetscBool interpolate = user->interpolate;
570: PetscReal refinementLimit = user->refinementLimit;
571: size_t len;
575: PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
576: PetscStrlen(filename, &len);
577: if (!len) {
578: PetscInt d;
580: if (user->periodicity[0] || user->periodicity[1] || user->periodicity[2]) for (d = 0; d < dim; ++d) user->cells[d] = PetscMax(user->cells[d], 3);
581: DMPlexCreateBoxMesh(comm, dim, user->simplex, user->cells, NULL, NULL, user->periodicity, interpolate, dm);
582: PetscObjectSetName((PetscObject) *dm, "Mesh");
583: } else {
584: DMPlexCreateFromFile(comm, filename, interpolate, dm);
585: DMPlexSetRefinementUniform(*dm, PETSC_FALSE);
586: }
587: {
588: PetscPartitioner part;
589: DM refinedMesh = NULL;
590: DM distributedMesh = NULL;
592: /* Refine mesh using a volume constraint */
593: if (refinementLimit > 0.0) {
594: DMPlexSetRefinementLimit(*dm, refinementLimit);
595: DMRefine(*dm, comm, &refinedMesh);
596: if (refinedMesh) {
597: const char *name;
599: PetscObjectGetName((PetscObject) *dm, &name);
600: PetscObjectSetName((PetscObject) refinedMesh, name);
601: DMDestroy(dm);
602: *dm = refinedMesh;
603: }
604: }
605: /* Distribute mesh over processes */
606: DMPlexGetPartitioner(*dm,&part);
607: PetscPartitionerSetFromOptions(part);
608: DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
609: if (distributedMesh) {
610: DMDestroy(dm);
611: *dm = distributedMesh;
612: }
613: }
614: if (interpolate) {
615: if (user->bcType == NEUMANN) {
616: DMLabel label;
618: DMCreateLabel(*dm, "boundary");
619: DMGetLabel(*dm, "boundary", &label);
620: DMPlexMarkBoundaryFaces(*dm, 1, label);
621: } else if (user->bcType == DIRICHLET) {
622: PetscBool hasLabel;
624: DMHasLabel(*dm,"marker",&hasLabel);
625: if (!hasLabel) {CreateBCLabel(*dm, "marker");}
626: }
627: }
628: {
629: char convType[256];
630: PetscBool flg;
632: PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
633: PetscOptionsFList("-dm_plex_convert_type","Convert DMPlex to another format","ex12",DMList,DMPLEX,convType,256,&flg);
634: PetscOptionsEnd();
635: if (flg) {
636: DM dmConv;
638: DMConvert(*dm,convType,&dmConv);
639: if (dmConv) {
640: DMDestroy(dm);
641: *dm = dmConv;
642: }
643: }
644: }
645: DMLocalizeCoordinates(*dm); /* needed for periodic */
646: DMSetFromOptions(*dm);
647: DMViewFromOptions(*dm, NULL, "-dm_view");
648: if (user->viewHierarchy) {
649: DM cdm = *dm;
650: PetscInt i = 0;
651: char buf[256];
653: while (cdm) {
654: DMSetUp(cdm);
655: DMGetCoarseDM(cdm, &cdm);
656: ++i;
657: }
658: cdm = *dm;
659: while (cdm) {
660: PetscViewer viewer;
661: PetscBool isHDF5, isVTK;
663: --i;
664: PetscViewerCreate(comm,&viewer);
665: PetscViewerSetType(viewer,PETSCVIEWERHDF5);
666: PetscViewerSetOptionsPrefix(viewer,"hierarchy_");
667: PetscViewerSetFromOptions(viewer);
668: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&isHDF5);
669: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERVTK,&isVTK);
670: if (isHDF5) {
671: PetscSNPrintf(buf, 256, "ex12-%d.h5", i);
672: } else if (isVTK) {
673: PetscSNPrintf(buf, 256, "ex12-%d.vtu", i);
674: PetscViewerPushFormat(viewer,PETSC_VIEWER_VTK_VTU);
675: } else {
676: PetscSNPrintf(buf, 256, "ex12-%d", i);
677: }
678: PetscViewerFileSetMode(viewer,FILE_MODE_WRITE);
679: PetscViewerFileSetName(viewer,buf);
680: DMView(cdm, viewer);
681: PetscViewerDestroy(&viewer);
682: DMGetCoarseDM(cdm, &cdm);
683: }
684: }
685: PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
686: return(0);
687: }
689: static PetscErrorCode SetupProblem(DM dm, AppCtx *user)
690: {
691: PetscDS prob;
692: const PetscInt id = 1;
696: DMGetDS(dm, &prob);
697: switch (user->variableCoefficient) {
698: case COEFF_NONE:
699: if (user->periodicity[0]) {
700: if (user->periodicity[1]) {
701: PetscDSSetResidual(prob, 0, f0_xytrig_u, f1_u);
702: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
703: } else {
704: PetscDSSetResidual(prob, 0, f0_xtrig_u, f1_u);
705: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
706: }
707: } else {
708: PetscDSSetResidual(prob, 0, f0_u, f1_u);
709: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
710: }
711: break;
712: case COEFF_ANALYTIC:
713: PetscDSSetResidual(prob, 0, f0_analytic_u, f1_analytic_u);
714: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_analytic_uu);
715: break;
716: case COEFF_FIELD:
717: PetscDSSetResidual(prob, 0, f0_analytic_u, f1_field_u);
718: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_field_uu);
719: break;
720: case COEFF_NONLINEAR:
721: PetscDSSetResidual(prob, 0, f0_analytic_nonlinear_u, f1_analytic_nonlinear_u);
722: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_analytic_nonlinear_uu);
723: break;
724: case COEFF_CIRCLE:
725: PetscDSSetResidual(prob, 0, f0_circle_u, f1_u);
726: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
727: break;
728: case COEFF_CROSS:
729: PetscDSSetResidual(prob, 0, f0_cross_u, f1_u);
730: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
731: break;
732: case COEFF_CHECKERBOARD_0:
733: PetscDSSetResidual(prob, 0, f0_checkerboard_0_u, f1_field_u);
734: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_field_uu);
735: break;
736: default: SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid variable coefficient type %d", user->variableCoefficient);
737: }
738: switch (user->dim) {
739: case 2:
740: switch (user->variableCoefficient) {
741: case COEFF_CIRCLE:
742: user->exactFuncs[0] = circle_u_2d;break;
743: case COEFF_CROSS:
744: user->exactFuncs[0] = cross_u_2d;break;
745: case COEFF_CHECKERBOARD_0:
746: user->exactFuncs[0] = zero;break;
747: default:
748: if (user->periodicity[0]) {
749: if (user->periodicity[1]) {
750: user->exactFuncs[0] = xytrig_u_2d;
751: } else {
752: user->exactFuncs[0] = xtrig_u_2d;
753: }
754: } else {
755: user->exactFuncs[0] = quadratic_u_2d;
756: user->exactFields[0] = quadratic_u_field_2d;
757: }
758: }
759: if (user->bcType == NEUMANN) {PetscDSSetBdResidual(prob, 0, f0_bd_u, f1_bd_zero);}
760: break;
761: case 3:
762: user->exactFuncs[0] = quadratic_u_3d;
763: user->exactFields[0] = quadratic_u_field_3d;
764: if (user->bcType == NEUMANN) {PetscDSSetBdResidual(prob, 0, f0_bd_u, f1_bd_zero);}
765: break;
766: default:
767: SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %d", user->dim);
768: }
769: /* Setup constants */
770: switch (user->variableCoefficient) {
771: case COEFF_CHECKERBOARD_0:
772: {
773: PetscScalar constants[2];
775: constants[0] = user->div;
776: constants[1] = user->k;
777: PetscDSSetConstants(prob, 2, constants);
778: }
779: break;
780: default: break;
781: }
782: PetscDSSetExactSolution(prob, 0, user->exactFuncs[0], user);
783: /* Setup Boundary Conditions */
784: if (user->bcType != NONE) {
785: DMAddBoundary(dm, user->bcType == DIRICHLET ? (user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL) : DM_BC_NATURAL,
786: "wall", user->bcType == DIRICHLET ? "marker" : "boundary", 0, 0, NULL,
787: user->fieldBC ? (void (*)(void)) user->exactFields[0] : (void (*)(void)) user->exactFuncs[0], NULL, 1, &id, user);
788: }
789: return(0);
790: }
792: static PetscErrorCode SetupMaterial(DM dm, DM dmAux, AppCtx *user)
793: {
794: PetscErrorCode (*matFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {nu_2d};
795: void *ctx[1];
796: Vec nu;
797: PetscErrorCode ierr;
800: ctx[0] = user;
801: if (user->variableCoefficient == COEFF_CHECKERBOARD_0) {matFuncs[0] = checkerboardCoeff;}
802: DMCreateLocalVector(dmAux, &nu);
803: PetscObjectSetName((PetscObject) nu, "Coefficient");
804: DMProjectFunctionLocal(dmAux, 0.0, matFuncs, ctx, INSERT_ALL_VALUES, nu);
805: PetscObjectCompose((PetscObject) dm, "A", (PetscObject) nu);
806: VecDestroy(&nu);
807: return(0);
808: }
810: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
811: {
812: PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
813: Vec uexact;
814: PetscInt dim;
818: DMGetDimension(dm, &dim);
819: if (dim == 2) bcFuncs[0] = quadratic_u_2d;
820: else bcFuncs[0] = quadratic_u_3d;
821: DMCreateLocalVector(dmAux, &uexact);
822: DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
823: PetscObjectCompose((PetscObject) dm, "A", (PetscObject) uexact);
824: VecDestroy(&uexact);
825: return(0);
826: }
828: static PetscErrorCode SetupAuxDM(DM dm, PetscFE feAux, AppCtx *user)
829: {
830: DM dmAux, coordDM;
834: /* MUST call DMGetCoordinateDM() in order to get p4est setup if present */
835: DMGetCoordinateDM(dm, &coordDM);
836: if (!feAux) return(0);
837: DMClone(dm, &dmAux);
838: PetscObjectCompose((PetscObject) dm, "dmAux", (PetscObject) dmAux);
839: DMSetCoordinateDM(dmAux, coordDM);
840: DMSetField(dmAux, 0, NULL, (PetscObject) feAux);
841: DMCreateDS(dmAux);
842: if (user->fieldBC) {SetupBC(dm, dmAux, user);}
843: else {SetupMaterial(dm, dmAux, user);}
844: DMDestroy(&dmAux);
845: return(0);
846: }
848: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
849: {
850: DM cdm = dm;
851: const PetscInt dim = user->dim;
852: PetscFE fe, feAux = NULL;
853: PetscBool simplex = user->simplex;
854: MPI_Comm comm;
858: /* Create finite element for each field and auxiliary field */
859: PetscObjectGetComm((PetscObject) dm, &comm);
860: PetscFECreateDefault(comm, dim, 1, simplex, NULL, -1, &fe);
861: PetscObjectSetName((PetscObject) fe, "potential");
862: if (user->variableCoefficient == COEFF_FIELD || user->variableCoefficient == COEFF_CHECKERBOARD_0) {
863: PetscFECreateDefault(comm, dim, 1, simplex, "mat_", -1, &feAux);
864: PetscObjectSetName((PetscObject) feAux, "coefficient");
865: PetscFECopyQuadrature(fe, feAux);
866: } else if (user->fieldBC) {
867: PetscFECreateDefault(comm, dim, 1, simplex, "bc_", -1, &feAux);
868: PetscFECopyQuadrature(fe, feAux);
869: }
870: /* Set discretization and boundary conditions for each mesh */
871: DMSetField(dm, 0, NULL, (PetscObject) fe);
872: DMCreateDS(dm);
873: SetupProblem(dm, user);
874: while (cdm) {
875: SetupAuxDM(cdm, feAux, user);
876: if (user->bcType == DIRICHLET && user->interpolate) {
877: PetscBool hasLabel;
879: DMHasLabel(cdm, "marker", &hasLabel);
880: if (!hasLabel) {CreateBCLabel(cdm, "marker");}
881: }
882: DMCopyDisc(dm, cdm);
883: DMGetCoarseDM(cdm, &cdm);
884: }
885: PetscFEDestroy(&fe);
886: PetscFEDestroy(&feAux);
887: return(0);
888: }
890: #include "petsc/private/petscimpl.h"
892: /*
893: MonitorError - Outputs the error at each iteration of an iterative solver.
895: Collective on KSP
897: Input Parameters:
898: + ksp - the KSP
899: . its - iteration number
900: . rnorm - 2-norm, preconditioned residual value (may be estimated).
901: - ctx - monitor context
903: Level: intermediate
905: .seealso: KSPMonitorSet(), KSPMonitorTrueResidual(), KSPMonitorResidual()
906: */
907: static PetscErrorCode MonitorError(KSP ksp, PetscInt its, PetscReal rnorm, void *ctx)
908: {
909: AppCtx *user = (AppCtx *) ctx;
910: DM dm;
911: Vec du = NULL, r;
912: PetscInt level = 0;
913: PetscBool hasLevel;
914: #if defined(PETSC_HAVE_HDF5)
915: PetscViewer viewer;
916: char buf[256];
917: #endif
921: KSPGetDM(ksp, &dm);
922: /* Calculate solution */
923: {
924: PC pc = user->pcmg; /* The MG PC */
925: DM fdm = NULL, cdm = NULL;
926: KSP fksp, cksp;
927: Vec fu, cu = NULL;
928: PetscInt levels, l;
930: KSPBuildSolution(ksp, NULL, &du);
931: PetscObjectComposedDataGetInt((PetscObject) ksp, PetscMGLevelId, level, hasLevel);
932: PCMGGetLevels(pc, &levels);
933: PCMGGetSmoother(pc, levels-1, &fksp);
934: KSPBuildSolution(fksp, NULL, &fu);
935: for (l = levels-1; l > level; --l) {
936: Mat R;
937: Vec s;
939: PCMGGetSmoother(pc, l-1, &cksp);
940: KSPGetDM(cksp, &cdm);
941: DMGetGlobalVector(cdm, &cu);
942: PCMGGetRestriction(pc, l, &R);
943: PCMGGetRScale(pc, l, &s);
944: MatRestrict(R, fu, cu);
945: VecPointwiseMult(cu, cu, s);
946: if (l < levels-1) {DMRestoreGlobalVector(fdm, &fu);}
947: fdm = cdm;
948: fu = cu;
949: }
950: if (levels-1 > level) {
951: VecAXPY(du, 1.0, cu);
952: DMRestoreGlobalVector(cdm, &cu);
953: }
954: }
955: /* Calculate error */
956: DMGetGlobalVector(dm, &r);
957: DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
958: VecAXPY(r,-1.0,du);
959: PetscObjectSetName((PetscObject) r, "solution error");
960: /* View error */
961: #if defined(PETSC_HAVE_HDF5)
962: PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
963: PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
964: VecView(r, viewer);
965: PetscViewerDestroy(&viewer);
966: #endif
967: DMRestoreGlobalVector(dm, &r);
968: return(0);
969: }
971: /*@C
972: SNESMonitorError - Outputs the error at each iteration of an iterative solver.
974: Collective on SNES
976: Input Parameters:
977: + snes - the SNES
978: . its - iteration number
979: . rnorm - 2-norm of residual
980: - ctx - user context
982: Level: intermediate
984: .seealso: SNESMonitorDefault(), SNESMonitorSet(), SNESMonitorSolution()
985: @*/
986: static PetscErrorCode SNESMonitorError(SNES snes, PetscInt its, PetscReal rnorm, void *ctx)
987: {
988: AppCtx *user = (AppCtx *) ctx;
989: DM dm;
990: Vec u, r;
991: PetscInt level = -1;
992: PetscBool hasLevel;
993: #if defined(PETSC_HAVE_HDF5)
994: PetscViewer viewer;
995: #endif
996: char buf[256];
1000: SNESGetDM(snes, &dm);
1001: /* Calculate error */
1002: SNESGetSolution(snes, &u);
1003: DMGetGlobalVector(dm, &r);
1004: PetscObjectSetName((PetscObject) r, "solution error");
1005: DMProjectFunction(dm, 0.0, user->exactFuncs, NULL, INSERT_ALL_VALUES, r);
1006: VecAXPY(r, -1.0, u);
1007: /* View error */
1008: PetscObjectComposedDataGetInt((PetscObject) snes, PetscMGLevelId, level, hasLevel);
1009: PetscSNPrintf(buf, 256, "ex12-%D.h5", level);
1010: #if defined(PETSC_HAVE_HDF5)
1011: PetscViewerHDF5Open(PETSC_COMM_WORLD, buf, FILE_MODE_APPEND, &viewer);
1012: VecView(r, viewer);
1013: PetscViewerDestroy(&viewer);
1014: /* Cleanup */
1015: DMRestoreGlobalVector(dm, &r);
1016: return(0);
1017: #else
1018: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"You need to configure with --download-hdf5");
1019: #endif
1020: }
1022: int main(int argc, char **argv)
1023: {
1024: DM dm; /* Problem specification */
1025: SNES snes; /* nonlinear solver */
1026: Vec u; /* solution vector */
1027: Mat A,J; /* Jacobian matrix */
1028: MatNullSpace nullSpace; /* May be necessary for Neumann conditions */
1029: AppCtx user; /* user-defined work context */
1030: JacActionCtx userJ; /* context for Jacobian MF action */
1031: PetscReal error = 0.0; /* L_2 error in the solution */
1032: PetscBool isFAS;
1035: PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
1036: ProcessOptions(PETSC_COMM_WORLD, &user);
1037: SNESCreate(PETSC_COMM_WORLD, &snes);
1038: CreateMesh(PETSC_COMM_WORLD, &user, &dm);
1039: SNESSetDM(snes, dm);
1040: DMSetApplicationContext(dm, &user);
1042: PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
1043: SetupDiscretization(dm, &user);
1045: DMCreateGlobalVector(dm, &u);
1046: PetscObjectSetName((PetscObject) u, "potential");
1048: DMCreateMatrix(dm, &J);
1049: if (user.jacobianMF) {
1050: PetscInt M, m, N, n;
1052: MatGetSize(J, &M, &N);
1053: MatGetLocalSize(J, &m, &n);
1054: MatCreate(PETSC_COMM_WORLD, &A);
1055: MatSetSizes(A, m, n, M, N);
1056: MatSetType(A, MATSHELL);
1057: MatSetUp(A);
1058: #if 0
1059: MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction);
1060: #endif
1062: userJ.dm = dm;
1063: userJ.J = J;
1064: userJ.user = &user;
1066: DMCreateLocalVector(dm, &userJ.u);
1067: if (user.fieldBC) {DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);}
1068: else {DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);}
1069: MatShellSetContext(A, &userJ);
1070: } else {
1071: A = J;
1072: }
1074: nullSpace = NULL;
1075: if (user.bcType != DIRICHLET) {
1076: MatNullSpaceCreate(PetscObjectComm((PetscObject) dm), PETSC_TRUE, 0, NULL, &nullSpace);
1077: MatSetNullSpace(A, nullSpace);
1078: }
1080: DMPlexSetSNESLocalFEM(dm,&user,&user,&user);
1081: SNESSetJacobian(snes, A, J, NULL, NULL);
1083: SNESSetFromOptions(snes);
1085: if (user.fieldBC) {DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u);}
1086: else {DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);}
1087: if (user.restart) {
1088: #if defined(PETSC_HAVE_HDF5)
1089: PetscViewer viewer;
1091: PetscViewerCreate(PETSC_COMM_WORLD, &viewer);
1092: PetscViewerSetType(viewer, PETSCVIEWERHDF5);
1093: PetscViewerFileSetMode(viewer, FILE_MODE_READ);
1094: PetscViewerFileSetName(viewer, user.filename);
1095: PetscViewerHDF5PushGroup(viewer, "/fields");
1096: VecLoad(u, viewer);
1097: PetscViewerHDF5PopGroup(viewer);
1098: PetscViewerDestroy(&viewer);
1099: #endif
1100: }
1101: if (user.showInitial) {
1102: Vec lv;
1103: DMGetLocalVector(dm, &lv);
1104: DMGlobalToLocalBegin(dm, u, INSERT_VALUES, lv);
1105: DMGlobalToLocalEnd(dm, u, INSERT_VALUES, lv);
1106: DMPrintLocalVec(dm, "Local function", 1.0e-10, lv);
1107: DMRestoreLocalVector(dm, &lv);
1108: }
1109: if (user.viewHierarchy) {
1110: SNES lsnes;
1111: KSP ksp;
1112: PC pc;
1113: PetscInt numLevels, l;
1114: PetscBool isMG;
1116: PetscObjectTypeCompare((PetscObject) snes, SNESFAS, &isFAS);
1117: if (isFAS) {
1118: SNESFASGetLevels(snes, &numLevels);
1119: for (l = 0; l < numLevels; ++l) {
1120: SNESFASGetCycleSNES(snes, l, &lsnes);
1121: SNESMonitorSet(lsnes, SNESMonitorError, &user, NULL);
1122: }
1123: } else {
1124: SNESGetKSP(snes, &ksp);
1125: KSPGetPC(ksp, &pc);
1126: PetscObjectTypeCompare((PetscObject) pc, PCMG, &isMG);
1127: if (isMG) {
1128: user.pcmg = pc;
1129: PCMGGetLevels(pc, &numLevels);
1130: for (l = 0; l < numLevels; ++l) {
1131: PCMGGetSmootherDown(pc, l, &ksp);
1132: KSPMonitorSet(ksp, MonitorError, &user, NULL);
1133: }
1134: }
1135: }
1136: }
1137: if (user.runType == RUN_FULL || user.runType == RUN_EXACT) {
1138: PetscErrorCode (*initialGuess[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {zero};
1140: if (user.nonzInit) initialGuess[0] = ecks;
1141: if (user.runType == RUN_FULL) {
1142: DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
1143: }
1144: if (user.debug) {
1145: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1146: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1147: }
1148: VecViewFromOptions(u, NULL, "-guess_vec_view");
1149: SNESSolve(snes, NULL, u);
1150: SNESGetSolution(snes, &u);
1151: SNESGetDM(snes, &dm);
1153: if (user.showSolution) {
1154: PetscPrintf(PETSC_COMM_WORLD, "Solution\n");
1155: VecChop(u, 3.0e-9);
1156: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
1157: }
1158: } else if (user.runType == RUN_PERF) {
1159: Vec r;
1160: PetscReal res = 0.0;
1162: SNESGetFunction(snes, &r, NULL, NULL);
1163: SNESComputeFunction(snes, u, r);
1164: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1165: VecChop(r, 1.0e-10);
1166: VecNorm(r, NORM_2, &res);
1167: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res);
1168: } else {
1169: Vec r;
1170: PetscReal res = 0.0, tol = 1.0e-11;
1172: /* Check discretization error */
1173: SNESGetFunction(snes, &r, NULL, NULL);
1174: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
1175: if (!user.quiet) {VecView(u, PETSC_VIEWER_STDOUT_WORLD);}
1176: DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
1177: if (error < tol) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < %2.1e\n", (double)tol);}
1178: else {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", (double)error);}
1179: /* Check residual */
1180: SNESComputeFunction(snes, u, r);
1181: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
1182: VecChop(r, 1.0e-10);
1183: if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1184: VecNorm(r, NORM_2, &res);
1185: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res);
1186: /* Check Jacobian */
1187: {
1188: Vec b;
1190: SNESComputeJacobian(snes, u, A, A);
1191: VecDuplicate(u, &b);
1192: VecSet(r, 0.0);
1193: SNESComputeFunction(snes, r, b);
1194: MatMult(A, u, r);
1195: VecAXPY(r, 1.0, b);
1196: PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n");
1197: VecChop(r, 1.0e-10);
1198: if (!user.quiet) {VecView(r, PETSC_VIEWER_STDOUT_WORLD);}
1199: VecNorm(r, NORM_2, &res);
1200: PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", (double)res);
1201: /* check solver */
1202: if (user.checkksp) {
1203: KSP ksp;
1205: if (nullSpace) {
1206: MatNullSpaceRemove(nullSpace, u);
1207: }
1208: SNESComputeJacobian(snes, u, A, J);
1209: MatMult(A, u, b);
1210: SNESGetKSP(snes, &ksp);
1211: KSPSetOperators(ksp, A, J);
1212: KSPSolve(ksp, b, r);
1213: VecAXPY(r, -1.0, u);
1214: VecNorm(r, NORM_2, &res);
1215: PetscPrintf(PETSC_COMM_WORLD, "KSP Error: %g\n", (double)res);
1216: }
1217: VecDestroy(&b);
1218: }
1219: }
1220: VecViewFromOptions(u, NULL, "-vec_view");
1221: {
1222: Vec nu;
1224: PetscObjectQuery((PetscObject) dm, "A", (PetscObject *) &nu);
1225: if (nu) {VecViewFromOptions(nu, NULL, "-coeff_view");}
1226: }
1228: if (user.bdIntegral) {
1229: DMLabel label;
1230: PetscInt id = 1;
1231: PetscScalar bdInt = 0.0;
1232: PetscReal exact = 3.3333333333;
1234: DMGetLabel(dm, "marker", &label);
1235: DMPlexComputeBdIntegral(dm, u, label, 1, &id, bd_integral_2d, &bdInt, NULL);
1236: PetscPrintf(PETSC_COMM_WORLD, "Solution boundary integral: %.4g\n", (double) PetscAbsScalar(bdInt));
1237: 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);
1238: }
1240: MatNullSpaceDestroy(&nullSpace);
1241: if (user.jacobianMF) {VecDestroy(&userJ.u);}
1242: if (A != J) {MatDestroy(&A);}
1243: MatDestroy(&J);
1244: VecDestroy(&u);
1245: SNESDestroy(&snes);
1246: DMDestroy(&dm);
1247: PetscFree2(user.exactFuncs, user.exactFields);
1248: PetscFree(user.kgrid);
1249: PetscFinalize();
1250: return ierr;
1251: }
1253: /*TEST
1254: # 2D serial P1 test 0-4
1255: test:
1256: suffix: 2d_p1_0
1257: requires: triangle
1258: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1260: test:
1261: suffix: 2d_p1_1
1262: requires: triangle
1263: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1265: test:
1266: suffix: 2d_p1_2
1267: requires: triangle
1268: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1270: test:
1271: suffix: 2d_p1_neumann_0
1272: requires: triangle
1273: 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
1275: test:
1276: suffix: 2d_p1_neumann_1
1277: requires: triangle
1278: args: -run_type test -refinement_limit 0.0625 -bc_type neumann -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1280: # 2D serial P2 test 5-8
1281: test:
1282: suffix: 2d_p2_0
1283: requires: triangle
1284: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1286: test:
1287: suffix: 2d_p2_1
1288: requires: triangle
1289: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1291: test:
1292: suffix: 2d_p2_neumann_0
1293: requires: triangle
1294: 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
1296: test:
1297: suffix: 2d_p2_neumann_1
1298: requires: triangle
1299: 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
1301: test:
1302: suffix: bd_int_0
1303: requires: triangle
1304: args: -run_type test -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet
1306: test:
1307: suffix: bd_int_1
1308: requires: triangle
1309: args: -run_type test -dm_refine 2 -bc_type dirichlet -interpolate 1 -petscspace_degree 2 -bd_integral -dm_view -quiet
1311: # 3D serial P1 test 9-12
1312: test:
1313: suffix: 3d_p1_0
1314: requires: ctetgen
1315: 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
1317: test:
1318: suffix: 3d_p1_1
1319: requires: ctetgen
1320: 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
1322: test:
1323: suffix: 3d_p1_2
1324: requires: ctetgen
1325: 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
1327: test:
1328: suffix: 3d_p1_neumann_0
1329: requires: ctetgen
1330: 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
1332: # Analytic variable coefficient 13-20
1333: test:
1334: suffix: 13
1335: requires: triangle
1336: args: -run_type test -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1337: test:
1338: suffix: 14
1339: requires: triangle
1340: args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1341: test:
1342: suffix: 15
1343: requires: triangle
1344: args: -run_type test -refinement_limit 0.0 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1345: test:
1346: suffix: 16
1347: requires: triangle
1348: args: -run_type test -refinement_limit 0.0625 -variable_coefficient analytic -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1349: test:
1350: suffix: 17
1351: requires: ctetgen
1352: 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
1354: test:
1355: suffix: 18
1356: requires: ctetgen
1357: 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
1359: test:
1360: suffix: 19
1361: requires: ctetgen
1362: 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
1364: test:
1365: suffix: 20
1366: requires: ctetgen
1367: 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
1369: # P1 variable coefficient 21-28
1370: test:
1371: suffix: 21
1372: requires: triangle
1373: 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
1375: test:
1376: suffix: 22
1377: requires: triangle
1378: 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
1380: test:
1381: suffix: 23
1382: requires: triangle
1383: 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
1385: test:
1386: suffix: 24
1387: requires: triangle
1388: 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
1390: test:
1391: suffix: 25
1392: requires: ctetgen
1393: 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
1395: test:
1396: suffix: 26
1397: requires: ctetgen
1398: 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
1400: test:
1401: suffix: 27
1402: requires: ctetgen
1403: 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
1405: test:
1406: suffix: 28
1407: requires: ctetgen
1408: 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
1410: # P0 variable coefficient 29-36
1411: test:
1412: suffix: 29
1413: requires: triangle
1414: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1416: test:
1417: suffix: 30
1418: requires: triangle
1419: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1421: test:
1422: suffix: 31
1423: requires: triangle
1424: args: -run_type test -refinement_limit 0.0 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1426: test:
1427: requires: triangle
1428: suffix: 32
1429: args: -run_type test -refinement_limit 0.0625 -variable_coefficient field -interpolate 1 -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1431: test:
1432: requires: ctetgen
1433: suffix: 33
1434: 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
1436: test:
1437: suffix: 34
1438: requires: ctetgen
1439: 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
1441: test:
1442: suffix: 35
1443: requires: ctetgen
1444: 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
1446: test:
1447: suffix: 36
1448: requires: ctetgen
1449: 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
1451: # Full solve 39-44
1452: test:
1453: suffix: 39
1454: requires: triangle !single
1455: 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
1456: test:
1457: suffix: 40
1458: requires: triangle !single
1459: 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
1460: test:
1461: suffix: 41
1462: requires: triangle !single
1463: 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
1464: test:
1465: suffix: 42
1466: requires: triangle !single
1467: 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
1468: test:
1469: suffix: 43
1470: requires: triangle !single
1471: nsize: 2
1472: 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: 44
1476: requires: triangle !single
1477: nsize: 2
1478: 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
1480: # These tests use a loose tolerance just to exercise the PtAP operations for MATIS and multiple PCBDDC setup calls inside PCMG
1481: testset:
1482: requires: triangle !single
1483: nsize: 3
1484: 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
1485: test:
1486: suffix: gmg_bddc
1487: filter: sed -e "s/CONVERGED_FNORM_RELATIVE iterations 3/CONVERGED_FNORM_RELATIVE iterations 4/g"
1488: args: -mg_levels_pc_type jacobi
1489: test:
1490: filter: sed -e "s/iterations [0-4]/iterations 4/g"
1491: suffix: gmg_bddc_lev
1492: args: -mg_levels_pc_type bddc
1494: # Restarting
1495: testset:
1496: suffix: restart
1497: requires: hdf5 triangle !complex
1498: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1
1499: test:
1500: args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1501: test:
1502: args: -f sol.h5 -restart
1504: # Periodicity
1505: test:
1506: suffix: periodic_0
1507: requires: triangle
1508: args: -run_type full -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 -snes_converged_reason ::ascii_info_detail
1510: test:
1511: requires: !complex
1512: suffix: periodic_1
1513: 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
1515: # 2D serial P1 test with field bc
1516: test:
1517: suffix: field_bc_2d_p1_0
1518: requires: triangle
1519: 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
1521: test:
1522: suffix: field_bc_2d_p1_1
1523: requires: triangle
1524: 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
1526: test:
1527: suffix: field_bc_2d_p1_neumann_0
1528: requires: triangle
1529: 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
1531: test:
1532: suffix: field_bc_2d_p1_neumann_1
1533: requires: triangle
1534: 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
1536: # 3D serial P1 test with field bc
1537: test:
1538: suffix: field_bc_3d_p1_0
1539: requires: ctetgen
1540: 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
1542: test:
1543: suffix: field_bc_3d_p1_1
1544: requires: ctetgen
1545: 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
1547: test:
1548: suffix: field_bc_3d_p1_neumann_0
1549: requires: ctetgen
1550: 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
1552: test:
1553: suffix: field_bc_3d_p1_neumann_1
1554: requires: ctetgen
1555: 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
1557: # 2D serial P2 test with field bc
1558: test:
1559: suffix: field_bc_2d_p2_0
1560: requires: triangle
1561: 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
1563: test:
1564: suffix: field_bc_2d_p2_1
1565: requires: triangle
1566: 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
1568: test:
1569: suffix: field_bc_2d_p2_neumann_0
1570: requires: triangle
1571: 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
1573: test:
1574: suffix: field_bc_2d_p2_neumann_1
1575: requires: triangle
1576: 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
1578: # 3D serial P2 test with field bc
1579: test:
1580: suffix: field_bc_3d_p2_0
1581: requires: ctetgen
1582: 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
1584: test:
1585: suffix: field_bc_3d_p2_1
1586: requires: ctetgen
1587: 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
1589: test:
1590: suffix: field_bc_3d_p2_neumann_0
1591: requires: ctetgen
1592: 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
1594: test:
1595: suffix: field_bc_3d_p2_neumann_1
1596: requires: ctetgen
1597: 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
1599: # Full solve simplex: Convergence
1600: test:
1601: suffix: 3d_p1_conv
1602: requires: ctetgen
1603: args: -run_type full -dim 3 -cells 1,1,1 -dm_refine 1 -bc_type dirichlet -interpolate 1 -petscspace_degree 1 \
1604: -snes_convergence_estimate -convest_num_refine 1 -pc_type lu
1606: # Full solve simplex: PCBDDC
1607: test:
1608: suffix: tri_bddc
1609: requires: triangle !single
1610: nsize: 5
1611: 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
1613: # Full solve simplex: PCBDDC
1614: test:
1615: suffix: tri_parmetis_bddc
1616: requires: triangle !single parmetis
1617: nsize: 4
1618: 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
1620: testset:
1621: 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
1622: nsize: 5
1623: output_file: output/ex12_quad_bddc.out
1624: filter: sed -e "s/aijcusparse/aij/g" -e "s/aijviennacl/aij/g" -e "s/factorization: cusparse/factorization: petsc/g"
1625: test:
1626: requires: !single
1627: suffix: quad_bddc
1628: test:
1629: requires: !single cuda
1630: suffix: quad_bddc_cuda
1631: args: -matis_localmat_type aijcusparse -pc_bddc_dirichlet_pc_factor_mat_solver_type cusparse -pc_bddc_neumann_pc_factor_mat_solver_type cusparse
1632: test:
1633: requires: !single viennacl
1634: suffix: quad_bddc_viennacl
1635: args: -matis_localmat_type aijviennacl
1637: # Full solve simplex: ASM
1638: test:
1639: suffix: tri_q2q1_asm_lu
1640: requires: triangle !single
1641: 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
1643: test:
1644: suffix: tri_q2q1_msm_lu
1645: requires: triangle !single
1646: 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
1648: test:
1649: suffix: tri_q2q1_asm_sor
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 sor -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_sor
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 sor -snes_monitor_short -ksp_monitor_short -snes_converged_reason ::ascii_info_detail -ksp_converged_reason -snes_view -show_solution 0
1658: # Full solve simplex: FAS
1659: test:
1660: suffix: fas_newton_0
1661: requires: triangle !single
1662: 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
1664: test:
1665: suffix: fas_newton_1
1666: requires: triangle !single
1667: 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
1668: filter: sed -e "s/total number of linear solver iterations=14/total number of linear solver iterations=15/g"
1670: test:
1671: suffix: fas_ngs_0
1672: requires: triangle !single
1673: 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
1675: test:
1676: suffix: fas_newton_coarse_0
1677: requires: pragmatic triangle
1678: TODO: broken
1679: 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
1681: test:
1682: suffix: mg_newton_coarse_0
1683: requires: triangle pragmatic
1684: TODO: broken
1685: 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
1687: test:
1688: suffix: mg_newton_coarse_1
1689: requires: triangle pragmatic
1690: TODO: broken
1691: 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
1693: test:
1694: suffix: mg_newton_coarse_2
1695: requires: triangle pragmatic
1696: TODO: broken
1697: 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
1699: # Full solve tensor
1700: test:
1701: suffix: tensor_plex_2d
1702: args: -run_type test -refinement_limit 0.0 -simplex 0 -interpolate -bc_type dirichlet -petscspace_degree 1 -dm_refine_hierarchy 2 -cells 2,2
1704: test:
1705: suffix: tensor_p4est_2d
1706: requires: p4est
1707: 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
1709: test:
1710: suffix: tensor_plex_3d
1711: 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
1713: test:
1714: suffix: tensor_p4est_3d
1715: requires: p4est
1716: 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
1718: test:
1719: suffix: p4est_test_q2_conformal_serial
1720: requires: p4est
1721: 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
1723: test:
1724: suffix: p4est_test_q2_conformal_parallel
1725: requires: p4est
1726: nsize: 7
1727: 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
1729: test:
1730: suffix: p4est_test_q2_conformal_parallel_parmetis
1731: requires: parmetis p4est
1732: nsize: 4
1733: 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
1735: test:
1736: suffix: p4est_test_q2_nonconformal_serial
1737: requires: p4est
1738: filter: grep -v "CG or CGNE: variant"
1739: 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
1741: test:
1742: suffix: p4est_test_q2_nonconformal_parallel
1743: requires: p4est
1744: filter: grep -v "CG or CGNE: variant"
1745: nsize: 7
1746: 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
1748: test:
1749: suffix: p4est_test_q2_nonconformal_parallel_parmetis
1750: requires: parmetis p4est
1751: nsize: 4
1752: 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
1754: test:
1755: suffix: p4est_exact_q2_conformal_serial
1756: requires: p4est !single !complex !__float128
1757: 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
1759: test:
1760: suffix: p4est_exact_q2_conformal_parallel
1761: requires: p4est !single !complex !__float128
1762: nsize: 4
1763: 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
1765: test:
1766: suffix: p4est_exact_q2_conformal_parallel_parmetis
1767: requires: parmetis p4est !single
1768: nsize: 4
1769: 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
1771: test:
1772: suffix: p4est_exact_q2_nonconformal_serial
1773: requires: p4est
1774: 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
1776: test:
1777: suffix: p4est_exact_q2_nonconformal_parallel
1778: requires: p4est
1779: nsize: 7
1780: 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
1782: test:
1783: suffix: p4est_exact_q2_nonconformal_parallel_parmetis
1784: requires: parmetis p4est
1785: nsize: 4
1786: 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
1788: test:
1789: suffix: p4est_full_q2_nonconformal_serial
1790: requires: p4est !single
1791: filter: grep -v "variant HERMITIAN"
1792: 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
1794: test:
1795: suffix: p4est_full_q2_nonconformal_parallel
1796: requires: p4est !single
1797: filter: grep -v "variant HERMITIAN"
1798: nsize: 7
1799: 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
1801: test:
1802: suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1803: requires: p4est !single
1804: filter: grep -v "variant HERMITIAN"
1805: nsize: 7
1806: 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
1808: test:
1809: suffix: p4est_full_q2_nonconformal_parallel_bddc
1810: requires: p4est !single
1811: filter: grep -v "variant HERMITIAN"
1812: nsize: 7
1813: 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
1815: test:
1816: TODO: broken
1817: suffix: p4est_fas_q2_conformal_serial
1818: requires: p4est !complex !__float128
1819: 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
1821: test:
1822: TODO: broken
1823: suffix: p4est_fas_q2_nonconformal_serial
1824: requires: p4est
1825: 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
1827: test:
1828: suffix: fas_newton_0_p4est
1829: requires: p4est !single !__float128
1830: 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
1832: # Full solve simplicial AMR
1833: test:
1834: suffix: tri_p1_adapt_0
1835: requires: pragmatic
1836: TODO: broken
1837: 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
1839: test:
1840: suffix: tri_p1_adapt_1
1841: requires: pragmatic
1842: TODO: broken
1843: 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
1845: test:
1846: suffix: tri_p1_adapt_analytic_0
1847: requires: pragmatic
1848: TODO: broken
1849: 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
1851: # Full solve tensor AMR
1852: test:
1853: suffix: quad_q1_adapt_0
1854: requires: p4est
1855: 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
1856: filter: grep -v DM_
1858: test:
1859: suffix: amr_0
1860: nsize: 5
1861: 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
1863: test:
1864: suffix: amr_1
1865: requires: p4est !complex
1866: 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
1868: test:
1869: suffix: p4est_solve_bddc
1870: requires: p4est !complex
1871: 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
1872: nsize: 4
1874: test:
1875: suffix: p4est_solve_fas
1876: requires: p4est
1877: 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
1878: nsize: 4
1879: TODO: identical machine two runs produce slightly different solver trackers
1881: test:
1882: suffix: p4est_convergence_test_1
1883: requires: p4est
1884: 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
1885: nsize: 4
1887: test:
1888: suffix: p4est_convergence_test_2
1889: requires: p4est
1890: 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
1892: test:
1893: suffix: p4est_convergence_test_3
1894: requires: p4est
1895: 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
1897: test:
1898: suffix: p4est_convergence_test_4
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 5 -dm_forest_initial_refinement 5 -dm_forest_maximum_refinement 7 -dm_p4est_refine_pattern hash
1901: timeoutfactor: 5
1903: # Serial tests with GLVis visualization
1904: test:
1905: suffix: glvis_2d_tet_p1
1906: 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
1907: test:
1908: suffix: glvis_2d_tet_p2
1909: 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
1910: test:
1911: suffix: glvis_2d_hex_p1
1912: args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -simplex 0 -dm_refine 1
1913: test:
1914: suffix: glvis_2d_hex_p2
1915: args: -quiet -run_type test -interpolate 1 -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -simplex 0 -dm_refine 1
1916: test:
1917: suffix: glvis_2d_hex_p2_p4est
1918: requires: p4est
1919: 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
1920: test:
1921: suffix: glvis_2d_tet_p0
1922: 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
1923: test:
1924: suffix: glvis_2d_hex_p0
1925: args: -run_type exact -interpolate 1 -guess_vec_view glvis: -nonzero_initial_guess 1 -cells 5,7 -simplex 0 -petscspace_degree 0
1927: # PCHPDDM tests
1928: testset:
1929: nsize: 4
1930: requires: hpddm slepc !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1931: 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
1932: test:
1933: suffix: quad_singular_hpddm
1934: args: -cells 6,7
1935: test:
1936: requires: p4est
1937: suffix: p4est_singular_2d_hpddm
1938: args: -dm_plex_convert_type p4est -dm_forest_minimum_refinement 1 -dm_forest_initial_refinement 3 -dm_forest_maximum_refinement 3
1939: test:
1940: requires: p4est
1941: suffix: p4est_nc_singular_2d_hpddm
1942: 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
1943: testset:
1944: nsize: 4
1945: requires: hpddm slepc triangle !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1946: 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
1947: test:
1948: args: -pc_hpddm_coarse_mat_type baij -options_left no
1949: suffix: tri_hpddm_reuse_baij
1950: test:
1951: requires: !complex
1952: suffix: tri_hpddm_reuse
1953: testset:
1954: nsize: 4
1955: requires: hpddm slepc !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1956: 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
1957: test:
1958: args: -pc_hpddm_coarse_mat_type baij -options_left no
1959: suffix: quad_hpddm_reuse_baij
1960: test:
1961: requires: !complex
1962: suffix: quad_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 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
1967: test:
1968: args: -pc_hpddm_coarse_mat_type baij -options_left no
1969: suffix: quad_hpddm_reuse_threshold_baij
1970: test:
1971: requires: !complex
1972: suffix: quad_hpddm_reuse_threshold
1973: testset:
1974: nsize: 4
1975: requires: hpddm slepc parmetis !single define(PETSC_HAVE_DYNAMIC_LIBRARIES) define(PETSC_USE_SHARED_LIBRARIES)
1976: filter: sed -e "s/linear solver iterations=17/linear solver iterations=16/g"
1977: 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
1978: test:
1979: args: -pc_hpddm_coarse_mat_type baij -options_left no
1980: 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"
1981: suffix: tri_parmetis_hpddm_baij
1982: test:
1983: 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"
1984: requires: !complex
1985: suffix: tri_parmetis_hpddm
1987: # 2D serial P1 tests for adaptive MG
1988: test:
1989: suffix: 2d_p1_adaptmg_0
1990: requires: triangle bamg
1991: args: -dm_refine_hierarchy 3 -cells 4,4 -interpolate -bc_type dirichlet -petscspace_degree 1 \
1992: -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
1993: -snes_max_it 1 -ksp_converged_reason \
1994: -ksp_rtol 1e-8 -pc_type mg
1995: # -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
1996: test:
1997: suffix: 2d_p1_adaptmg_1
1998: requires: triangle bamg
1999: args: -dm_refine_hierarchy 3 -cells 4,4 -interpolate -bc_type dirichlet -petscspace_degree 1 \
2000: -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
2001: -snes_max_it 1 -ksp_converged_reason \
2002: -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 \
2003: -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
2005: TEST*/