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 {
24: NEUMANN,
25: DIRICHLET,
26: NONE
27: } BCType;
28: typedef enum {
29: RUN_FULL,
30: RUN_EXACT,
31: RUN_TEST,
32: RUN_PERF
33: } RunType;
34: typedef enum {
35: COEFF_NONE,
36: COEFF_ANALYTIC,
37: COEFF_FIELD,
38: COEFF_NONLINEAR,
39: COEFF_BALL,
40: COEFF_CROSS,
41: COEFF_CHECKERBOARD_0,
42: COEFF_CHECKERBOARD_1
43: } CoeffType;
45: typedef struct {
46: RunType runType; /* Whether to run tests, or solve the full problem */
47: PetscBool jacobianMF; /* Whether to calculate the Jacobian action on the fly */
48: PetscBool showInitial, showSolution, restart, quiet, nonzInit;
49: /* Problem definition */
50: BCType bcType;
51: CoeffType variableCoefficient;
52: PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx);
53: PetscBool fieldBC;
54: void (**exactFields)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]);
55: PetscBool bdIntegral; /* Compute the integral of the solution on the boundary */
56: /* Reproducing tests from SISC 40(3), pp. A1473-A1493, 2018 */
57: PetscInt div; /* Number of divisions */
58: PetscInt k; /* Parameter for checkerboard coefficient */
59: PetscInt *kgrid; /* Random parameter grid */
60: PetscBool rand; /* Make random assignments */
61: /* Solver */
62: PC pcmg; /* This is needed for error monitoring */
63: PetscBool checkksp; /* Whether to check the KSPSolve for runType == RUN_TEST */
64: } AppCtx;
66: static PetscErrorCode zero(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
67: {
68: u[0] = 0.0;
69: return 0;
70: }
72: static PetscErrorCode ecks(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
73: {
74: u[0] = x[0];
75: return 0;
76: }
78: /*
79: In 2D for Dirichlet conditions, we use exact solution:
81: u = x^2 + y^2
82: f = 4
84: so that
86: -\Delta u + f = -4 + 4 = 0
88: For Neumann conditions, we have
90: -\nabla u \cdot -\hat y |_{y=0} = (2y)|_{y=0} = 0 (bottom)
91: -\nabla u \cdot \hat y |_{y=1} = -(2y)|_{y=1} = -2 (top)
92: -\nabla u \cdot -\hat x |_{x=0} = (2x)|_{x=0} = 0 (left)
93: -\nabla u \cdot \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)
95: Which we can express as
97: \nabla u \cdot \hat n|_\Gamma = {2 x, 2 y} \cdot \hat n = 2 (x + y)
99: The boundary integral of this solution is (assuming we are not orienting the edges)
101: \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
102: */
103: static PetscErrorCode quadratic_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
104: {
105: *u = x[0] * x[0] + x[1] * x[1];
106: return 0;
107: }
109: static void quadratic_u_field_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
110: {
111: uexact[0] = a[0];
112: }
114: static PetscErrorCode ball_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
115: {
116: const PetscReal alpha = 500.;
117: const PetscReal radius2 = PetscSqr(0.15);
118: const PetscReal r2 = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5);
119: const PetscReal xi = alpha * (radius2 - r2);
121: *u = PetscTanhScalar(xi) + 1.0;
122: return 0;
123: }
125: static PetscErrorCode cross_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
126: {
127: const PetscReal alpha = 50 * 4;
128: const PetscReal xy = (x[0] - 0.5) * (x[1] - 0.5);
130: *u = PetscSinReal(alpha * xy) * (alpha * PetscAbsReal(xy) < 2 * PETSC_PI ? 1 : 0.01);
131: return 0;
132: }
134: static void f0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
135: {
136: f0[0] = 4.0;
137: }
139: static void f0_ball_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
140: {
141: PetscInt d;
142: const PetscReal alpha = 500., radius2 = PetscSqr(0.15);
143: PetscReal r2, xi;
145: for (d = 0, r2 = 0.0; d < dim; ++d) r2 += PetscSqr(x[d] - 0.5);
146: xi = alpha * (radius2 - r2);
147: f0[0] = (-2.0 * dim * alpha - 8.0 * PetscSqr(alpha) * r2 * PetscTanhReal(xi)) * PetscSqr(1.0 / PetscCoshReal(xi));
148: }
150: static void f0_cross_u_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
151: {
152: const PetscReal alpha = 50 * 4;
153: const PetscReal xy = (x[0] - 0.5) * (x[1] - 0.5);
155: f0[0] = PetscSinReal(alpha * xy) * (alpha * PetscAbsReal(xy) < 2 * PETSC_PI ? 1 : 0.01);
156: }
158: static void f0_checkerboard_0_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
159: {
160: f0[0] = -20.0 * PetscExpReal(-(PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5)));
161: }
163: static void f0_bd_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
164: {
165: PetscInt d;
166: for (d = 0, f0[0] = 0.0; d < dim; ++d) f0[0] += -n[d] * 2.0 * x[d];
167: }
169: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
170: static void f1_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
171: {
172: PetscInt d;
173: for (d = 0; d < dim; ++d) f1[d] = u_x[d];
174: }
176: /* < \nabla v, \nabla u + {\nabla u}^T >
177: This just gives \nabla u, give the perdiagonal for the transpose */
178: static void g3_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
179: {
180: PetscInt d;
181: for (d = 0; d < dim; ++d) g3[d * dim + d] = 1.0;
182: }
184: /*
185: In 2D for x periodicity and y Dirichlet conditions, we use exact solution:
187: u = sin(2 pi x)
188: f = -4 pi^2 sin(2 pi x)
190: so that
192: -\Delta u + f = 4 pi^2 sin(2 pi x) - 4 pi^2 sin(2 pi x) = 0
193: */
194: static PetscErrorCode xtrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
195: {
196: *u = PetscSinReal(2.0 * PETSC_PI * x[0]);
197: return 0;
198: }
200: static void f0_xtrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
201: {
202: f0[0] = -4.0 * PetscSqr(PETSC_PI) * PetscSinReal(2.0 * PETSC_PI * x[0]);
203: }
205: /*
206: In 2D for x-y periodicity, we use exact solution:
208: u = sin(2 pi x) sin(2 pi y)
209: f = -8 pi^2 sin(2 pi x)
211: so that
213: -\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
214: */
215: static PetscErrorCode xytrig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
216: {
217: *u = PetscSinReal(2.0 * PETSC_PI * x[0]) * PetscSinReal(2.0 * PETSC_PI * x[1]);
218: return 0;
219: }
221: static void f0_xytrig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
222: {
223: f0[0] = -8.0 * PetscSqr(PETSC_PI) * PetscSinReal(2.0 * PETSC_PI * x[0]);
224: }
226: /*
227: In 2D for Dirichlet conditions with a variable coefficient, we use exact solution:
229: u = x^2 + y^2
230: f = 6 (x + y)
231: nu = (x + y)
233: so that
235: -\div \nu \grad u + f = -6 (x + y) + 6 (x + y) = 0
236: */
237: static PetscErrorCode nu_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
238: {
239: *u = x[0] + x[1];
240: return 0;
241: }
243: static PetscErrorCode checkerboardCoeff(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
244: {
245: AppCtx *user = (AppCtx *)ctx;
246: PetscInt div = user->div;
247: PetscInt k = user->k;
248: PetscInt mask = 0, ind = 0, d;
251: for (d = 0; d < dim; ++d) mask = (mask + (PetscInt)(x[d] * div)) % 2;
252: if (user->kgrid) {
253: for (d = 0; d < dim; ++d) {
254: if (d > 0) ind *= dim;
255: ind += (PetscInt)(x[d] * div);
256: }
257: k = user->kgrid[ind];
258: }
259: u[0] = mask ? 1.0 : PetscPowRealInt(10.0, -k);
260: return 0;
261: }
263: void f0_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
264: {
265: f0[0] = 6.0 * (x[0] + x[1]);
266: }
268: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
269: void f1_analytic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
270: {
271: PetscInt d;
272: for (d = 0; d < dim; ++d) f1[d] = (x[0] + x[1]) * u_x[d];
273: }
275: void f1_field_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
276: {
277: PetscInt d;
278: for (d = 0; d < dim; ++d) f1[d] = a[0] * u_x[d];
279: }
281: /* < \nabla v, \nabla u + {\nabla u}^T >
282: This just gives \nabla u, give the perdiagonal for the transpose */
283: void g3_analytic_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
284: {
285: PetscInt d;
286: for (d = 0; d < dim; ++d) g3[d * dim + d] = x[0] + x[1];
287: }
289: void g3_field_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
290: {
291: PetscInt d;
292: for (d = 0; d < dim; ++d) g3[d * dim + d] = a[0];
293: }
295: /*
296: In 2D for Dirichlet conditions with a nonlinear coefficient (p-Laplacian with p = 4), we use exact solution:
298: u = x^2 + y^2
299: f = 16 (x^2 + y^2)
300: nu = 1/2 |grad u|^2
302: so that
304: -\div \nu \grad u + f = -16 (x^2 + y^2) + 16 (x^2 + y^2) = 0
305: */
306: void f0_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
307: {
308: f0[0] = 16.0 * (x[0] * x[0] + x[1] * x[1]);
309: }
311: /* gradU[comp*dim+d] = {u_x, u_y} or {u_x, u_y, u_z} */
312: void f1_analytic_nonlinear_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
313: {
314: PetscScalar nu = 0.0;
315: PetscInt d;
316: for (d = 0; d < dim; ++d) nu += u_x[d] * u_x[d];
317: for (d = 0; d < dim; ++d) f1[d] = 0.5 * nu * u_x[d];
318: }
320: /*
321: grad (u + eps w) - grad u = eps grad w
323: 1/2 |grad (u + eps w)|^2 grad (u + eps w) - 1/2 |grad u|^2 grad u
324: = 1/2 (|grad u|^2 + 2 eps <grad u,grad w>) (grad u + eps grad w) - 1/2 |grad u|^2 grad u
325: = 1/2 (eps |grad u|^2 grad w + 2 eps <grad u,grad w> grad u)
326: = eps (1/2 |grad u|^2 grad w + grad u <grad u,grad w>)
327: */
328: void g3_analytic_nonlinear_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
329: {
330: PetscScalar nu = 0.0;
331: PetscInt d, e;
332: for (d = 0; d < dim; ++d) nu += u_x[d] * u_x[d];
333: for (d = 0; d < dim; ++d) {
334: g3[d * dim + d] = 0.5 * nu;
335: for (e = 0; e < dim; ++e) g3[d * dim + e] += u_x[d] * u_x[e];
336: }
337: }
339: /*
340: In 3D for Dirichlet conditions we use exact solution:
342: u = 2/3 (x^2 + y^2 + z^2)
343: f = 4
345: so that
347: -\Delta u + f = -2/3 * 6 + 4 = 0
349: For Neumann conditions, we have
351: -\nabla u \cdot -\hat z |_{z=0} = (2z)|_{z=0} = 0 (bottom)
352: -\nabla u \cdot \hat z |_{z=1} = -(2z)|_{z=1} = -2 (top)
353: -\nabla u \cdot -\hat y |_{y=0} = (2y)|_{y=0} = 0 (front)
354: -\nabla u \cdot \hat y |_{y=1} = -(2y)|_{y=1} = -2 (back)
355: -\nabla u \cdot -\hat x |_{x=0} = (2x)|_{x=0} = 0 (left)
356: -\nabla u \cdot \hat x |_{x=1} = -(2x)|_{x=1} = -2 (right)
358: Which we can express as
360: \nabla u \cdot \hat n|_\Gamma = {2 x, 2 y, 2z} \cdot \hat n = 2 (x + y + z)
361: */
362: static PetscErrorCode quadratic_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
363: {
364: *u = 2.0 * (x[0] * x[0] + x[1] * x[1] + x[2] * x[2]) / 3.0;
365: return 0;
366: }
368: static PetscErrorCode ball_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
369: {
370: const PetscReal alpha = 500.;
371: const PetscReal radius2 = PetscSqr(0.15);
372: const PetscReal r2 = PetscSqr(x[0] - 0.5) + PetscSqr(x[1] - 0.5) + PetscSqr(x[2] - 0.5);
373: const PetscReal xi = alpha * (radius2 - r2);
375: *u = PetscTanhScalar(xi) + 1.0;
376: return 0;
377: }
379: static void quadratic_u_field_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uexact[])
380: {
381: uexact[0] = a[0];
382: }
384: static PetscErrorCode cross_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar *u, void *ctx)
385: {
386: const PetscReal alpha = 50 * 4;
387: const PetscReal xyz = (x[0] - 0.5) * (x[1] - 0.5) * (x[2] - 0.5);
389: *u = PetscSinReal(alpha * xyz) * (alpha * PetscAbsReal(xyz) < 2 * PETSC_PI ? (alpha * PetscAbsReal(xyz) > -2 * PETSC_PI ? 1.0 : 0.01) : 0.01);
390: return 0;
391: }
393: static void f0_cross_u_3d(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
394: {
395: const PetscReal alpha = 50 * 4;
396: const PetscReal xyz = (x[0] - 0.5) * (x[1] - 0.5) * (x[2] - 0.5);
398: f0[0] = PetscSinReal(alpha * xyz) * (alpha * PetscAbsReal(xyz) < 2 * PETSC_PI ? (alpha * PetscAbsReal(xyz) > -2 * PETSC_PI ? 1.0 : 0.01) : 0.01);
399: }
401: static void bd_integral_2d(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar *uint)
402: {
403: uint[0] = u[0];
404: }
406: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
407: {
408: const char *bcTypes[3] = {"neumann", "dirichlet", "none"};
409: const char *runTypes[4] = {"full", "exact", "test", "perf"};
410: const char *coeffTypes[8] = {"none", "analytic", "field", "nonlinear", "ball", "cross", "checkerboard_0", "checkerboard_1"};
411: PetscInt bc, run, coeff;
414: options->runType = RUN_FULL;
415: options->bcType = DIRICHLET;
416: options->variableCoefficient = COEFF_NONE;
417: options->fieldBC = PETSC_FALSE;
418: options->jacobianMF = PETSC_FALSE;
419: options->showInitial = PETSC_FALSE;
420: options->showSolution = PETSC_FALSE;
421: options->restart = PETSC_FALSE;
422: options->quiet = PETSC_FALSE;
423: options->nonzInit = PETSC_FALSE;
424: options->bdIntegral = PETSC_FALSE;
425: options->checkksp = PETSC_FALSE;
426: options->div = 4;
427: options->k = 1;
428: options->kgrid = NULL;
429: options->rand = PETSC_FALSE;
431: PetscOptionsBegin(comm, "", "Poisson Problem Options", "DMPLEX");
432: run = options->runType;
433: PetscOptionsEList("-run_type", "The run type", "ex12.c", runTypes, 4, runTypes[options->runType], &run, NULL);
434: options->runType = (RunType)run;
435: bc = options->bcType;
436: PetscOptionsEList("-bc_type", "Type of boundary condition", "ex12.c", bcTypes, 3, bcTypes[options->bcType], &bc, NULL);
437: options->bcType = (BCType)bc;
438: coeff = options->variableCoefficient;
439: PetscOptionsEList("-variable_coefficient", "Type of variable coefficient", "ex12.c", coeffTypes, 8, coeffTypes[options->variableCoefficient], &coeff, NULL);
440: options->variableCoefficient = (CoeffType)coeff;
442: PetscOptionsBool("-field_bc", "Use a field representation for the BC", "ex12.c", options->fieldBC, &options->fieldBC, NULL);
443: PetscOptionsBool("-jacobian_mf", "Calculate the action of the Jacobian on the fly", "ex12.c", options->jacobianMF, &options->jacobianMF, NULL);
444: PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex12.c", options->showInitial, &options->showInitial, NULL);
445: PetscOptionsBool("-show_solution", "Output the solution for verification", "ex12.c", options->showSolution, &options->showSolution, NULL);
446: PetscOptionsBool("-restart", "Read in the mesh and solution from a file", "ex12.c", options->restart, &options->restart, NULL);
447: PetscOptionsBool("-quiet", "Don't print any vecs", "ex12.c", options->quiet, &options->quiet, NULL);
448: PetscOptionsBool("-nonzero_initial_guess", "nonzero initial guess", "ex12.c", options->nonzInit, &options->nonzInit, NULL);
449: PetscOptionsBool("-bd_integral", "Compute the integral of the solution on the boundary", "ex12.c", options->bdIntegral, &options->bdIntegral, NULL);
450: if (options->runType == RUN_TEST) PetscOptionsBool("-run_test_check_ksp", "Check solution of KSP", "ex12.c", options->checkksp, &options->checkksp, NULL);
451: PetscOptionsInt("-div", "The number of division for the checkerboard coefficient", "ex12.c", options->div, &options->div, NULL);
452: PetscOptionsInt("-k", "The exponent for the checkerboard coefficient", "ex12.c", options->k, &options->k, NULL);
453: PetscOptionsBool("-k_random", "Assign random k values to checkerboard", "ex12.c", options->rand, &options->rand, NULL);
454: PetscOptionsEnd();
455: return 0;
456: }
458: static PetscErrorCode CreateBCLabel(DM dm, const char name[])
459: {
460: DM plex;
461: DMLabel label;
464: DMCreateLabel(dm, name);
465: DMGetLabel(dm, name, &label);
466: DMConvert(dm, DMPLEX, &plex);
467: DMPlexMarkBoundaryFaces(plex, 1, label);
468: DMDestroy(&plex);
469: return 0;
470: }
472: static PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
473: {
475: DMCreate(comm, dm);
476: DMSetType(*dm, DMPLEX);
477: DMSetFromOptions(*dm);
478: {
479: char convType[256];
480: PetscBool flg;
482: PetscOptionsBegin(comm, "", "Mesh conversion options", "DMPLEX");
483: PetscOptionsFList("-dm_plex_convert_type", "Convert DMPlex to another format", "ex12", DMList, DMPLEX, convType, 256, &flg);
484: PetscOptionsEnd();
485: if (flg) {
486: DM dmConv;
488: DMConvert(*dm, convType, &dmConv);
489: if (dmConv) {
490: DMDestroy(dm);
491: *dm = dmConv;
492: }
493: DMSetFromOptions(*dm);
494: DMSetUp(*dm);
495: }
496: }
497: DMViewFromOptions(*dm, NULL, "-dm_view");
498: if (user->rand) {
499: PetscRandom r;
500: PetscReal val;
501: PetscInt dim, N, i;
503: DMGetDimension(*dm, &dim);
504: N = PetscPowInt(user->div, dim);
505: PetscMalloc1(N, &user->kgrid);
506: PetscRandomCreate(PETSC_COMM_SELF, &r);
507: PetscRandomSetFromOptions(r);
508: PetscRandomSetInterval(r, 0.0, user->k);
509: PetscRandomSetSeed(r, 1973);
510: PetscRandomSeed(r);
511: for (i = 0; i < N; ++i) {
512: PetscRandomGetValueReal(r, &val);
513: user->kgrid[i] = 1 + (PetscInt)val;
514: }
515: PetscRandomDestroy(&r);
516: }
517: return 0;
518: }
520: static PetscErrorCode SetupProblem(DM dm, AppCtx *user)
521: {
522: PetscDS ds;
523: DMLabel label;
524: PetscWeakForm wf;
525: const PetscReal *L;
526: const PetscInt id = 1;
527: PetscInt bd, dim;
530: DMGetDS(dm, &ds);
531: DMGetDimension(dm, &dim);
532: DMGetPeriodicity(dm, NULL, NULL, &L);
533: switch (user->variableCoefficient) {
534: case COEFF_NONE:
535: if (L && L[0]) {
536: if (L && L[1]) {
537: PetscDSSetResidual(ds, 0, f0_xytrig_u, f1_u);
538: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
539: } else {
540: PetscDSSetResidual(ds, 0, f0_xtrig_u, f1_u);
541: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
542: }
543: } else {
544: PetscDSSetResidual(ds, 0, f0_u, f1_u);
545: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
546: }
547: break;
548: case COEFF_ANALYTIC:
549: PetscDSSetResidual(ds, 0, f0_analytic_u, f1_analytic_u);
550: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_analytic_uu);
551: break;
552: case COEFF_FIELD:
553: PetscDSSetResidual(ds, 0, f0_analytic_u, f1_field_u);
554: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_field_uu);
555: break;
556: case COEFF_NONLINEAR:
557: PetscDSSetResidual(ds, 0, f0_analytic_nonlinear_u, f1_analytic_nonlinear_u);
558: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_analytic_nonlinear_uu);
559: break;
560: case COEFF_BALL:
561: PetscDSSetResidual(ds, 0, f0_ball_u, f1_u);
562: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
563: break;
564: case COEFF_CROSS:
565: switch (dim) {
566: case 2:
567: PetscDSSetResidual(ds, 0, f0_cross_u_2d, f1_u);
568: break;
569: case 3:
570: PetscDSSetResidual(ds, 0, f0_cross_u_3d, f1_u);
571: break;
572: default:
573: SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %" PetscInt_FMT, dim);
574: }
575: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);
576: break;
577: case COEFF_CHECKERBOARD_0:
578: PetscDSSetResidual(ds, 0, f0_checkerboard_0_u, f1_field_u);
579: PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_field_uu);
580: break;
581: default:
582: SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Invalid variable coefficient type %d", user->variableCoefficient);
583: }
584: switch (dim) {
585: case 2:
586: switch (user->variableCoefficient) {
587: case COEFF_BALL:
588: user->exactFuncs[0] = ball_u_2d;
589: break;
590: case COEFF_CROSS:
591: user->exactFuncs[0] = cross_u_2d;
592: break;
593: case COEFF_CHECKERBOARD_0:
594: user->exactFuncs[0] = zero;
595: break;
596: default:
597: if (L && L[0]) {
598: if (L && L[1]) {
599: user->exactFuncs[0] = xytrig_u_2d;
600: } else {
601: user->exactFuncs[0] = xtrig_u_2d;
602: }
603: } else {
604: user->exactFuncs[0] = quadratic_u_2d;
605: user->exactFields[0] = quadratic_u_field_2d;
606: }
607: }
608: if (user->bcType == NEUMANN) {
609: DMGetLabel(dm, "boundary", &label);
610: DMAddBoundary(dm, DM_BC_NATURAL, "wall", label, 1, &id, 0, 0, NULL, NULL, NULL, user, &bd);
611: PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL);
612: PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, 0, f0_bd_u, 0, NULL);
613: }
614: break;
615: case 3:
616: switch (user->variableCoefficient) {
617: case COEFF_BALL:
618: user->exactFuncs[0] = ball_u_3d;
619: break;
620: case COEFF_CROSS:
621: user->exactFuncs[0] = cross_u_3d;
622: break;
623: default:
624: user->exactFuncs[0] = quadratic_u_3d;
625: user->exactFields[0] = quadratic_u_field_3d;
626: }
627: if (user->bcType == NEUMANN) {
628: DMGetLabel(dm, "boundary", &label);
629: DMAddBoundary(dm, DM_BC_NATURAL, "wall", label, 1, &id, 0, 0, NULL, NULL, NULL, user, &bd);
630: PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL);
631: PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, 0, f0_bd_u, 0, NULL);
632: }
633: break;
634: default:
635: SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Invalid dimension %" PetscInt_FMT, dim);
636: }
637: /* Setup constants */
638: switch (user->variableCoefficient) {
639: case COEFF_CHECKERBOARD_0: {
640: PetscScalar constants[2];
642: constants[0] = user->div;
643: constants[1] = user->k;
644: PetscDSSetConstants(ds, 2, constants);
645: } break;
646: default:
647: break;
648: }
649: PetscDSSetExactSolution(ds, 0, user->exactFuncs[0], user);
650: /* Setup Boundary Conditions */
651: if (user->bcType == DIRICHLET) {
652: DMGetLabel(dm, "marker", &label);
653: if (!label) {
654: /* Right now, p4est cannot create labels immediately */
655: PetscDSAddBoundaryByName(ds, user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL, "wall", "marker", 1, &id, 0, 0, NULL, user->fieldBC ? (void (*)(void))user->exactFields[0] : (void (*)(void))user->exactFuncs[0], NULL, user, NULL);
656: } else {
657: DMAddBoundary(dm, user->fieldBC ? DM_BC_ESSENTIAL_FIELD : DM_BC_ESSENTIAL, "wall", label, 1, &id, 0, 0, NULL, user->fieldBC ? (void (*)(void))user->exactFields[0] : (void (*)(void))user->exactFuncs[0], NULL, user, NULL);
658: }
659: }
660: return 0;
661: }
663: static PetscErrorCode SetupMaterial(DM dm, DM dmAux, AppCtx *user)
664: {
665: PetscErrorCode (*matFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {nu_2d};
666: void *ctx[1];
667: Vec nu;
669: ctx[0] = user;
670: if (user->variableCoefficient == COEFF_CHECKERBOARD_0) matFuncs[0] = checkerboardCoeff;
671: DMCreateLocalVector(dmAux, &nu);
672: PetscObjectSetName((PetscObject)nu, "Coefficient");
673: DMProjectFunctionLocal(dmAux, 0.0, matFuncs, ctx, INSERT_ALL_VALUES, nu);
674: DMSetAuxiliaryVec(dm, NULL, 0, 0, nu);
675: VecDestroy(&nu);
676: return 0;
677: }
679: static PetscErrorCode SetupBC(DM dm, DM dmAux, AppCtx *user)
680: {
681: PetscErrorCode (*bcFuncs[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx);
682: Vec uexact;
683: PetscInt dim;
685: DMGetDimension(dm, &dim);
686: if (dim == 2) bcFuncs[0] = quadratic_u_2d;
687: else bcFuncs[0] = quadratic_u_3d;
688: DMCreateLocalVector(dmAux, &uexact);
689: DMProjectFunctionLocal(dmAux, 0.0, bcFuncs, NULL, INSERT_ALL_VALUES, uexact);
690: DMSetAuxiliaryVec(dm, NULL, 0, 0, uexact);
691: VecDestroy(&uexact);
692: return 0;
693: }
695: static PetscErrorCode SetupAuxDM(DM dm, PetscFE feAux, AppCtx *user)
696: {
697: DM dmAux, coordDM;
699: /* MUST call DMGetCoordinateDM() in order to get p4est setup if present */
700: DMGetCoordinateDM(dm, &coordDM);
701: if (!feAux) return 0;
702: DMClone(dm, &dmAux);
703: DMSetCoordinateDM(dmAux, coordDM);
704: DMSetField(dmAux, 0, NULL, (PetscObject)feAux);
705: DMCreateDS(dmAux);
706: if (user->fieldBC) SetupBC(dm, dmAux, user);
707: else SetupMaterial(dm, dmAux, user);
708: DMDestroy(&dmAux);
709: return 0;
710: }
712: static PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
713: {
714: DM plex, cdm = dm;
715: PetscFE fe, feAux = NULL;
716: PetscBool simplex;
717: PetscInt dim;
718: MPI_Comm comm;
721: DMGetDimension(dm, &dim);
722: DMConvert(dm, DMPLEX, &plex);
723: DMPlexIsSimplex(plex, &simplex);
724: DMDestroy(&plex);
725: PetscObjectGetComm((PetscObject)dm, &comm);
726: PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, NULL, -1, &fe);
727: PetscObjectSetName((PetscObject)fe, "potential");
728: if (user->variableCoefficient == COEFF_FIELD || user->variableCoefficient == COEFF_CHECKERBOARD_0) {
729: PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, "mat_", -1, &feAux);
730: PetscObjectSetName((PetscObject)feAux, "coefficient");
731: PetscFECopyQuadrature(fe, feAux);
732: } else if (user->fieldBC) {
733: PetscFECreateDefault(PETSC_COMM_SELF, dim, 1, simplex, "bc_", -1, &feAux);
734: PetscFECopyQuadrature(fe, feAux);
735: }
736: /* Set discretization and boundary conditions for each mesh */
737: DMSetField(dm, 0, NULL, (PetscObject)fe);
738: DMCreateDS(dm);
739: SetupProblem(dm, user);
740: while (cdm) {
741: SetupAuxDM(cdm, feAux, user);
742: if (user->bcType == DIRICHLET) {
743: PetscBool hasLabel;
745: DMHasLabel(cdm, "marker", &hasLabel);
746: if (!hasLabel) CreateBCLabel(cdm, "marker");
747: }
748: DMCopyDisc(dm, cdm);
749: DMGetCoarseDM(cdm, &cdm);
750: }
751: PetscFEDestroy(&fe);
752: PetscFEDestroy(&feAux);
753: return 0;
754: }
756: int main(int argc, char **argv)
757: {
758: DM dm; /* Problem specification */
759: SNES snes; /* nonlinear solver */
760: Vec u; /* solution vector */
761: Mat A, J; /* Jacobian matrix */
762: MatNullSpace nullSpace; /* May be necessary for Neumann conditions */
763: AppCtx user; /* user-defined work context */
764: JacActionCtx userJ; /* context for Jacobian MF action */
765: PetscReal error = 0.0; /* L_2 error in the solution */
768: PetscInitialize(&argc, &argv, NULL, help);
769: ProcessOptions(PETSC_COMM_WORLD, &user);
770: SNESCreate(PETSC_COMM_WORLD, &snes);
771: CreateMesh(PETSC_COMM_WORLD, &user, &dm);
772: SNESSetDM(snes, dm);
773: DMSetApplicationContext(dm, &user);
775: PetscMalloc2(1, &user.exactFuncs, 1, &user.exactFields);
776: SetupDiscretization(dm, &user);
778: DMCreateGlobalVector(dm, &u);
779: PetscObjectSetName((PetscObject)u, "potential");
781: DMCreateMatrix(dm, &J);
782: if (user.jacobianMF) {
783: PetscInt M, m, N, n;
785: MatGetSize(J, &M, &N);
786: MatGetLocalSize(J, &m, &n);
787: MatCreate(PETSC_COMM_WORLD, &A);
788: MatSetSizes(A, m, n, M, N);
789: MatSetType(A, MATSHELL);
790: MatSetUp(A);
791: #if 0
792: MatShellSetOperation(A, MATOP_MULT, (void (*)(void))FormJacobianAction);
793: #endif
795: userJ.dm = dm;
796: userJ.J = J;
797: userJ.user = &user;
799: DMCreateLocalVector(dm, &userJ.u);
800: if (user.fieldBC) DMProjectFieldLocal(dm, 0.0, userJ.u, user.exactFields, INSERT_BC_VALUES, userJ.u);
801: else DMProjectFunctionLocal(dm, 0.0, user.exactFuncs, NULL, INSERT_BC_VALUES, userJ.u);
802: MatShellSetContext(A, &userJ);
803: } else {
804: A = J;
805: }
807: nullSpace = NULL;
808: if (user.bcType != DIRICHLET) {
809: MatNullSpaceCreate(PetscObjectComm((PetscObject)dm), PETSC_TRUE, 0, NULL, &nullSpace);
810: MatSetNullSpace(A, nullSpace);
811: }
813: DMPlexSetSNESLocalFEM(dm, &user, &user, &user);
814: SNESSetJacobian(snes, A, J, NULL, NULL);
816: SNESSetFromOptions(snes);
818: if (user.fieldBC) DMProjectField(dm, 0.0, u, user.exactFields, INSERT_ALL_VALUES, u);
819: else DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);
820: if (user.restart) {
821: #if defined(PETSC_HAVE_HDF5)
822: PetscViewer viewer;
823: char filename[PETSC_MAX_PATH_LEN];
825: PetscOptionsGetString(NULL, NULL, "-dm_plex_filename", filename, sizeof(filename), NULL);
826: PetscViewerCreate(PETSC_COMM_WORLD, &viewer);
827: PetscViewerSetType(viewer, PETSCVIEWERHDF5);
828: PetscViewerFileSetMode(viewer, FILE_MODE_READ);
829: PetscViewerFileSetName(viewer, filename);
830: PetscViewerHDF5PushGroup(viewer, "/fields");
831: VecLoad(u, viewer);
832: PetscViewerHDF5PopGroup(viewer);
833: PetscViewerDestroy(&viewer);
834: #endif
835: }
836: if (user.showInitial) {
837: Vec lv;
838: DMGetLocalVector(dm, &lv);
839: DMGlobalToLocalBegin(dm, u, INSERT_VALUES, lv);
840: DMGlobalToLocalEnd(dm, u, INSERT_VALUES, lv);
841: DMPrintLocalVec(dm, "Local function", 1.0e-10, lv);
842: DMRestoreLocalVector(dm, &lv);
843: }
844: if (user.runType == RUN_FULL || user.runType == RUN_EXACT) {
845: PetscErrorCode (*initialGuess[1])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx) = {zero};
847: if (user.nonzInit) initialGuess[0] = ecks;
848: if (user.runType == RUN_FULL) DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
849: VecViewFromOptions(u, NULL, "-guess_vec_view");
850: SNESSolve(snes, NULL, u);
851: SNESGetSolution(snes, &u);
852: SNESGetDM(snes, &dm);
854: if (user.showSolution) {
855: PetscPrintf(PETSC_COMM_WORLD, "Solution\n");
856: VecChop(u, 3.0e-9);
857: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
858: }
859: } else if (user.runType == RUN_PERF) {
860: Vec r;
861: PetscReal res = 0.0;
863: SNESGetFunction(snes, &r, NULL, NULL);
864: SNESComputeFunction(snes, u, r);
865: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
866: VecChop(r, 1.0e-10);
867: VecNorm(r, NORM_2, &res);
868: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res);
869: } else {
870: Vec r;
871: PetscReal res = 0.0, tol = 1.0e-11;
873: /* Check discretization error */
874: SNESGetFunction(snes, &r, NULL, NULL);
875: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
876: if (!user.quiet) VecView(u, PETSC_VIEWER_STDOUT_WORLD);
877: DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
878: if (error < tol) PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < %2.1e\n", (double)tol);
879: else PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", (double)error);
880: /* Check residual */
881: SNESComputeFunction(snes, u, r);
882: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
883: VecChop(r, 1.0e-10);
884: if (!user.quiet) VecView(r, PETSC_VIEWER_STDOUT_WORLD);
885: VecNorm(r, NORM_2, &res);
886: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", (double)res);
887: /* Check Jacobian */
888: {
889: Vec b;
891: SNESComputeJacobian(snes, u, A, A);
892: VecDuplicate(u, &b);
893: VecSet(r, 0.0);
894: SNESComputeFunction(snes, r, b);
895: MatMult(A, u, r);
896: VecAXPY(r, 1.0, b);
897: PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n");
898: VecChop(r, 1.0e-10);
899: if (!user.quiet) VecView(r, PETSC_VIEWER_STDOUT_WORLD);
900: VecNorm(r, NORM_2, &res);
901: PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", (double)res);
902: /* check solver */
903: if (user.checkksp) {
904: KSP ksp;
906: if (nullSpace) MatNullSpaceRemove(nullSpace, u);
907: SNESComputeJacobian(snes, u, A, J);
908: MatMult(A, u, b);
909: SNESGetKSP(snes, &ksp);
910: KSPSetOperators(ksp, A, J);
911: KSPSolve(ksp, b, r);
912: VecAXPY(r, -1.0, u);
913: VecNorm(r, NORM_2, &res);
914: PetscPrintf(PETSC_COMM_WORLD, "KSP Error: %g\n", (double)res);
915: }
916: VecDestroy(&b);
917: }
918: }
919: VecViewFromOptions(u, NULL, "-vec_view");
920: {
921: Vec nu;
923: DMGetAuxiliaryVec(dm, NULL, 0, 0, &nu);
924: if (nu) VecViewFromOptions(nu, NULL, "-coeff_view");
925: }
927: if (user.bdIntegral) {
928: DMLabel label;
929: PetscInt id = 1;
930: PetscScalar bdInt = 0.0;
931: PetscReal exact = 3.3333333333;
933: DMGetLabel(dm, "marker", &label);
934: DMPlexComputeBdIntegral(dm, u, label, 1, &id, bd_integral_2d, &bdInt, NULL);
935: PetscPrintf(PETSC_COMM_WORLD, "Solution boundary integral: %.4g\n", (double)PetscAbsScalar(bdInt));
937: }
939: MatNullSpaceDestroy(&nullSpace);
940: if (user.jacobianMF) VecDestroy(&userJ.u);
941: if (A != J) MatDestroy(&A);
942: MatDestroy(&J);
943: VecDestroy(&u);
944: SNESDestroy(&snes);
945: DMDestroy(&dm);
946: PetscFree2(user.exactFuncs, user.exactFields);
947: PetscFree(user.kgrid);
948: PetscFinalize();
949: return 0;
950: }
952: /*TEST
953: # 2D serial P1 test 0-4
954: test:
955: suffix: 2d_p1_0
956: requires: triangle
957: args: -run_type test -bc_type dirichlet -dm_plex_interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
959: test:
960: suffix: 2d_p1_1
961: requires: triangle
962: args: -run_type test -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
964: test:
965: suffix: 2d_p1_2
966: requires: triangle
967: args: -run_type test -dm_refine_volume_limit_pre 0.0625 -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
969: test:
970: suffix: 2d_p1_neumann_0
971: requires: triangle
972: args: -dm_coord_space 0 -run_type test -bc_type neumann -dm_plex_boundary_label boundary -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
974: test:
975: suffix: 2d_p1_neumann_1
976: requires: triangle
977: args: -run_type test -dm_refine_volume_limit_pre 0.0625 -bc_type neumann -dm_plex_boundary_label boundary -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
979: # 2D serial P2 test 5-8
980: test:
981: suffix: 2d_p2_0
982: requires: triangle
983: args: -run_type test -bc_type dirichlet -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
985: test:
986: suffix: 2d_p2_1
987: requires: triangle
988: args: -run_type test -dm_refine_volume_limit_pre 0.0625 -bc_type dirichlet -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
990: test:
991: suffix: 2d_p2_neumann_0
992: requires: triangle
993: args: -dm_coord_space 0 -run_type test -bc_type neumann -dm_plex_boundary_label boundary -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
995: test:
996: suffix: 2d_p2_neumann_1
997: requires: triangle
998: args: -dm_coord_space 0 -run_type test -dm_refine_volume_limit_pre 0.0625 -bc_type neumann -dm_plex_boundary_label boundary -petscspace_degree 2 -show_initial -dm_plex_print_fem 1 -dm_view ascii::ascii_info_detail
1000: test:
1001: suffix: bd_int_0
1002: requires: triangle
1003: args: -run_type test -bc_type dirichlet -petscspace_degree 2 -bd_integral -dm_view -quiet
1005: test:
1006: suffix: bd_int_1
1007: requires: triangle
1008: args: -run_type test -dm_refine 2 -bc_type dirichlet -petscspace_degree 2 -bd_integral -dm_view -quiet
1010: # 3D serial P1 test 9-12
1011: test:
1012: suffix: 3d_p1_0
1013: requires: ctetgen
1014: args: -run_type test -dm_plex_dim 3 -bc_type dirichlet -dm_plex_interpolate 0 -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view
1016: test:
1017: suffix: 3d_p1_1
1018: requires: ctetgen
1019: args: -run_type test -dm_plex_dim 3 -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view
1021: test:
1022: suffix: 3d_p1_2
1023: requires: ctetgen
1024: args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -bc_type dirichlet -petscspace_degree 1 -show_initial -dm_plex_print_fem 1 -dm_view
1026: test:
1027: suffix: 3d_p1_neumann_0
1028: requires: ctetgen
1029: args: -run_type test -dm_plex_dim 3 -bc_type neumann -dm_plex_boundary_label boundary -petscspace_degree 1 -snes_fd -show_initial -dm_plex_print_fem 1 -dm_view
1031: # Analytic variable coefficient 13-20
1032: test:
1033: suffix: 13
1034: requires: triangle
1035: args: -run_type test -variable_coefficient analytic -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1036: test:
1037: suffix: 14
1038: requires: triangle
1039: args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient analytic -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1040: test:
1041: suffix: 15
1042: requires: triangle
1043: args: -run_type test -variable_coefficient analytic -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1044: test:
1045: suffix: 16
1046: requires: triangle
1047: args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient analytic -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1048: test:
1049: suffix: 17
1050: requires: ctetgen
1051: args: -run_type test -dm_plex_dim 3 -variable_coefficient analytic -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1053: test:
1054: suffix: 18
1055: requires: ctetgen
1056: args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient analytic -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1058: test:
1059: suffix: 19
1060: requires: ctetgen
1061: args: -run_type test -dm_plex_dim 3 -variable_coefficient analytic -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1063: test:
1064: suffix: 20
1065: requires: ctetgen
1066: args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient analytic -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1068: # P1 variable coefficient 21-28
1069: test:
1070: suffix: 21
1071: requires: triangle
1072: args: -run_type test -variable_coefficient field -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1074: test:
1075: suffix: 22
1076: requires: triangle
1077: args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient field -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1079: test:
1080: suffix: 23
1081: requires: triangle
1082: args: -run_type test -variable_coefficient field -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1084: test:
1085: suffix: 24
1086: requires: triangle
1087: args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient field -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1089: test:
1090: suffix: 25
1091: requires: ctetgen
1092: args: -run_type test -dm_plex_dim 3 -variable_coefficient field -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1094: test:
1095: suffix: 26
1096: requires: ctetgen
1097: args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient field -petscspace_degree 1 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1099: test:
1100: suffix: 27
1101: requires: ctetgen
1102: args: -run_type test -dm_plex_dim 3 -variable_coefficient field -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1104: test:
1105: suffix: 28
1106: requires: ctetgen
1107: args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient field -petscspace_degree 2 -mat_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1109: # P0 variable coefficient 29-36
1110: test:
1111: suffix: 29
1112: requires: triangle
1113: args: -run_type test -variable_coefficient field -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1115: test:
1116: suffix: 30
1117: requires: triangle
1118: args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient field -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1120: test:
1121: suffix: 31
1122: requires: triangle
1123: args: -run_type test -variable_coefficient field -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1125: test:
1126: requires: triangle
1127: suffix: 32
1128: args: -run_type test -dm_refine_volume_limit_pre 0.0625 -variable_coefficient field -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1130: test:
1131: requires: ctetgen
1132: suffix: 33
1133: args: -run_type test -dm_plex_dim 3 -variable_coefficient field -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1135: test:
1136: suffix: 34
1137: requires: ctetgen
1138: args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient field -petscspace_degree 1 -show_initial -dm_plex_print_fem 1
1140: test:
1141: suffix: 35
1142: requires: ctetgen
1143: args: -run_type test -dm_plex_dim 3 -variable_coefficient field -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1145: test:
1146: suffix: 36
1147: requires: ctetgen
1148: args: -run_type test -dm_plex_dim 3 -dm_refine_volume_limit_pre 0.0125 -variable_coefficient field -petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1150: # Full solve 39-44
1151: test:
1152: suffix: 39
1153: requires: triangle !single
1154: args: -run_type full -dm_refine_volume_limit_pre 0.015625 -petscspace_degree 2 -pc_type gamg -pc_gamg_esteig_ksp_type cg -pc_gamg_esteig_ksp_max_it 10 -snes_rtol 1.0e-6 -ksp_rtol 1.0e-7 -ksp_monitor -ksp_converged_reason -snes_monitor_short -snes_converged_reason ::ascii_info_detail
1155: test:
1156: suffix: 40
1157: requires: triangle !single
1158: args: -run_type full -dm_refine_volume_limit_pre 0.015625 -variable_coefficient nonlinear -petscspace_degree 2 -pc_type svd -ksp_rtol 1.0e-10 -snes_monitor_short -snes_converged_reason ::ascii_info_detail
1159: test:
1160: suffix: 41
1161: requires: triangle !single
1162: args: -run_type full -dm_refine_volume_limit_pre 0.03125 -variable_coefficient nonlinear -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
1163: test:
1164: suffix: 42
1165: requires: triangle !single
1166: args: -run_type full -dm_refine_volume_limit_pre 0.0625 -variable_coefficient nonlinear -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 -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
1167: test:
1168: suffix: 43
1169: requires: triangle !single
1170: nsize: 2
1171: args: -run_type full -dm_refine_volume_limit_pre 0.03125 -variable_coefficient nonlinear -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
1173: test:
1174: suffix: 44
1175: requires: triangle !single
1176: nsize: 2
1177: args: -run_type full -dm_refine_volume_limit_pre 0.0625 -variable_coefficient nonlinear -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
1179: # These tests use a loose tolerance just to exercise the PtAP operations for MATIS and multiple PCBDDC setup calls inside PCMG
1180: testset:
1181: requires: triangle !single
1182: nsize: 3
1183: args: -run_type full -petscspace_degree 1 -dm_mat_type is -pc_type mg -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
1184: test:
1185: suffix: gmg_bddc
1186: filter: sed -e "s/CONVERGED_FNORM_RELATIVE iterations 3/CONVERGED_FNORM_RELATIVE iterations 4/g"
1187: args: -mg_levels_pc_type jacobi
1188: test:
1189: filter: sed -e "s/iterations [0-4]/iterations 4/g"
1190: suffix: gmg_bddc_lev
1191: args: -mg_levels_pc_type bddc
1193: # Restarting
1194: testset:
1195: suffix: restart
1196: requires: hdf5 triangle !complex
1197: args: -run_type test -bc_type dirichlet -petscspace_degree 1
1198: test:
1199: args: -dm_view hdf5:sol.h5 -vec_view hdf5:sol.h5::append
1200: test:
1201: args: -dm_plex_filename sol.h5 -dm_plex_name box -restart
1203: # Periodicity
1204: test:
1205: suffix: periodic_0
1206: requires: triangle
1207: args: -run_type full -bc_type dirichlet -petscspace_degree 1 -snes_converged_reason ::ascii_info_detail
1209: test:
1210: requires: !complex
1211: suffix: periodic_1
1212: args: -quiet -run_type test -dm_plex_simplex 0 -dm_plex_box_faces 3,3 -dm_plex_box_bd periodic,periodic -vec_view vtk:test.vtu:vtk_vtu -petscspace_degree 1 -dm_refine 1
1214: # 2D serial P1 test with field bc
1215: test:
1216: suffix: field_bc_2d_p1_0
1217: requires: triangle
1218: args: -run_type test -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1220: test:
1221: suffix: field_bc_2d_p1_1
1222: requires: triangle
1223: args: -run_type test -dm_refine 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1225: test:
1226: suffix: field_bc_2d_p1_neumann_0
1227: requires: triangle
1228: args: -run_type test -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1230: test:
1231: suffix: field_bc_2d_p1_neumann_1
1232: requires: triangle
1233: args: -run_type test -dm_refine 1 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1235: # 3D serial P1 test with field bc
1236: test:
1237: suffix: field_bc_3d_p1_0
1238: requires: ctetgen
1239: args: -run_type test -dm_plex_dim 3 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1241: test:
1242: suffix: field_bc_3d_p1_1
1243: requires: ctetgen
1244: args: -run_type test -dm_plex_dim 3 -dm_refine 1 -bc_type dirichlet -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1246: test:
1247: suffix: field_bc_3d_p1_neumann_0
1248: requires: ctetgen
1249: args: -run_type test -dm_plex_dim 3 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1251: test:
1252: suffix: field_bc_3d_p1_neumann_1
1253: requires: ctetgen
1254: args: -run_type test -dm_plex_dim 3 -dm_refine 1 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 1 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1256: # 2D serial P2 test with field bc
1257: test:
1258: suffix: field_bc_2d_p2_0
1259: requires: triangle
1260: args: -run_type test -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1262: test:
1263: suffix: field_bc_2d_p2_1
1264: requires: triangle
1265: args: -run_type test -dm_refine 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1267: test:
1268: suffix: field_bc_2d_p2_neumann_0
1269: requires: triangle
1270: args: -run_type test -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1272: test:
1273: suffix: field_bc_2d_p2_neumann_1
1274: requires: triangle
1275: args: -run_type test -dm_refine 1 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1277: # 3D serial P2 test with field bc
1278: test:
1279: suffix: field_bc_3d_p2_0
1280: requires: ctetgen
1281: args: -run_type test -dm_plex_dim 3 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1283: test:
1284: suffix: field_bc_3d_p2_1
1285: requires: ctetgen
1286: args: -run_type test -dm_plex_dim 3 -dm_refine 1 -bc_type dirichlet -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1288: test:
1289: suffix: field_bc_3d_p2_neumann_0
1290: requires: ctetgen
1291: args: -run_type test -dm_plex_dim 3 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1293: test:
1294: suffix: field_bc_3d_p2_neumann_1
1295: requires: ctetgen
1296: args: -run_type test -dm_plex_dim 3 -dm_refine 1 -bc_type neumann -dm_plex_boundary_label boundary -field_bc -petscspace_degree 2 -bc_petscspace_degree 2 -show_initial -dm_plex_print_fem 1
1298: # Full solve simplex: Convergence
1299: test:
1300: suffix: 3d_p1_conv
1301: requires: ctetgen
1302: args: -run_type full -dm_plex_dim 3 -dm_refine 1 -bc_type dirichlet -petscspace_degree 1 \
1303: -snes_convergence_estimate -convest_num_refine 1 -pc_type lu
1305: # Full solve simplex: PCBDDC
1306: test:
1307: suffix: tri_bddc
1308: requires: triangle !single
1309: nsize: 5
1310: args: -run_type full -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -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
1312: # Full solve simplex: PCBDDC
1313: test:
1314: suffix: tri_parmetis_bddc
1315: requires: triangle !single parmetis
1316: nsize: 4
1317: args: -run_type full -petscpartitioner_type parmetis -dm_refine 2 -bc_type dirichlet -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
1319: testset:
1320: args: -run_type full -dm_plex_simplex 0 -dm_plex_box_faces 3,3 -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -petscspace_degree 2 -dm_mat_type is -pc_type bddc -ksp_type gmres -snes_monitor_short -ksp_monitor_short -snes_view -petscspace_poly_tensor -pc_bddc_corner_selection -ksp_rtol 1.e-9 -pc_bddc_use_edges 0
1321: nsize: 5
1322: output_file: output/ex12_quad_bddc.out
1323: filter: sed -e "s/aijcusparse/aij/g" -e "s/aijviennacl/aij/g" -e "s/factorization: cusparse/factorization: petsc/g"
1324: test:
1325: requires: !single
1326: suffix: quad_bddc
1327: test:
1328: requires: !single cuda
1329: suffix: quad_bddc_cuda
1330: args: -matis_localmat_type aijcusparse -pc_bddc_dirichlet_pc_factor_mat_solver_type cusparse -pc_bddc_neumann_pc_factor_mat_solver_type cusparse
1331: test:
1332: requires: !single viennacl
1333: suffix: quad_bddc_viennacl
1334: args: -matis_localmat_type aijviennacl
1336: # Full solve simplex: ASM
1337: test:
1338: suffix: tri_q2q1_asm_lu
1339: requires: triangle !single
1340: args: -run_type full -dm_refine 3 -bc_type dirichlet -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
1342: test:
1343: suffix: tri_q2q1_msm_lu
1344: requires: triangle !single
1345: args: -run_type full -dm_refine 3 -bc_type dirichlet -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
1347: test:
1348: suffix: tri_q2q1_asm_sor
1349: requires: triangle !single
1350: args: -run_type full -dm_refine 3 -bc_type dirichlet -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
1352: test:
1353: suffix: tri_q2q1_msm_sor
1354: requires: triangle !single
1355: args: -run_type full -dm_refine 3 -bc_type dirichlet -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
1357: # Full solve simplex: FAS
1358: test:
1359: suffix: fas_newton_0
1360: requires: triangle !single
1361: args: -run_type full -variable_coefficient nonlinear -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
1363: test:
1364: suffix: fas_newton_1
1365: requires: triangle !single
1366: args: -run_type full -dm_refine_hierarchy 3 -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
1367: filter: sed -e "s/total number of linear solver iterations=14/total number of linear solver iterations=15/g"
1369: test:
1370: suffix: fas_ngs_0
1371: requires: triangle !single
1372: args: -run_type full -variable_coefficient nonlinear -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
1374: # These two tests are broken because DMPlexComputeInjectorFEM() only works for regularly refined meshes
1375: test:
1376: suffix: fas_newton_coarse_0
1377: requires: pragmatic triangle
1378: TODO: broken
1379: args: -run_type full -variable_coefficient nonlinear -petscspace_degree 1 \
1380: -dm_refine 2 -dm_coarsen_hierarchy 1 -dm_plex_hash_location -dm_adaptor pragmatic \
1381: -snes_type fas -snes_fas_levels 2 -snes_converged_reason ::ascii_info_detail -snes_monitor_short -snes_view \
1382: -fas_coarse_pc_type svd -fas_coarse_ksp_rtol 1.0e-10 -fas_coarse_snes_monitor_short -fas_coarse_snes_linesearch_type basic \
1383: -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
1385: test:
1386: suffix: mg_newton_coarse_0
1387: requires: triangle pragmatic
1388: TODO: broken
1389: args: -run_type full -petscspace_degree 1 \
1390: -dm_refine 3 -dm_coarsen_hierarchy 3 -dm_plex_hash_location -dm_adaptor pragmatic \
1391: -snes_atol 1.0e-8 -snes_rtol 0.0 -snes_monitor_short -snes_converged_reason ::ascii_info_detail -snes_view \
1392: -ksp_type richardson -ksp_atol 1.0e-8 -ksp_rtol 0.0 -ksp_norm_type unpreconditioned -ksp_monitor_true_residual \
1393: -pc_type mg -pc_mg_levels 4 \
1394: -mg_levels_ksp_type gmres -mg_levels_pc_type ilu -mg_levels_ksp_max_it 10
1396: # Full solve tensor
1397: test:
1398: suffix: tensor_plex_2d
1399: args: -run_type test -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_refine_hierarchy 2
1401: test:
1402: suffix: tensor_p4est_2d
1403: requires: p4est
1404: args: -run_type test -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_forest_initial_refinement 2 -dm_forest_minimum_refinement 0 -dm_plex_convert_type p4est
1406: test:
1407: suffix: tensor_plex_3d
1408: args: -run_type test -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_plex_dim 3 -dm_refine_hierarchy 1 -dm_plex_box_faces 2,2,2
1410: test:
1411: suffix: tensor_p4est_3d
1412: requires: p4est
1413: args: -run_type test -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_forest_initial_refinement 1 -dm_forest_minimum_refinement 0 -dm_plex_dim 3 -dm_plex_convert_type p8est -dm_plex_box_faces 2,2,2
1415: test:
1416: suffix: p4est_test_q2_conformal_serial
1417: requires: p4est
1418: args: -run_type test -petscspace_degree 2 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2
1420: test:
1421: suffix: p4est_test_q2_conformal_parallel
1422: requires: p4est
1423: nsize: 7
1424: args: -run_type test -petscspace_degree 2 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type simple
1426: test:
1427: suffix: p4est_test_q2_conformal_parallel_parmetis
1428: requires: parmetis p4est
1429: nsize: 4
1430: args: -run_type test -petscspace_degree 2 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type parmetis
1432: test:
1433: suffix: p4est_test_q2_nonconformal_serial
1434: requires: p4est
1435: filter: grep -v "CG or CGNE: variant"
1436: args: -run_type test -petscspace_degree 2 -dm_plex_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
1438: test:
1439: suffix: p4est_test_q2_nonconformal_parallel
1440: requires: p4est
1441: filter: grep -v "CG or CGNE: variant"
1442: nsize: 7
1443: args: -run_type test -petscspace_degree 2 -dm_plex_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
1445: test:
1446: suffix: p4est_test_q2_nonconformal_parallel_parmetis
1447: requires: parmetis p4est
1448: nsize: 4
1449: args: -run_type test -petscspace_degree 2 -dm_plex_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
1451: test:
1452: suffix: p4est_exact_q2_conformal_serial
1453: requires: p4est !single !complex !__float128
1454: args: -run_type exact -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 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2
1456: test:
1457: suffix: p4est_exact_q2_conformal_parallel
1458: requires: p4est !single !complex !__float128
1459: nsize: 4
1460: args: -run_type exact -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 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2
1462: test:
1463: suffix: p4est_exact_q2_conformal_parallel_parmetis
1464: requires: parmetis p4est !single
1465: nsize: 4
1466: args: -run_type exact -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 -dm_plex_simplex 0 -dm_plex_convert_type p4est -dm_forest_minimum_refinement 0 -dm_forest_initial_refinement 2 -petscpartitioner_type parmetis
1468: test:
1469: suffix: p4est_exact_q2_nonconformal_serial
1470: requires: p4est
1471: args: -run_type exact -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 -dm_plex_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
1473: test:
1474: suffix: p4est_exact_q2_nonconformal_parallel
1475: requires: p4est
1476: nsize: 7
1477: args: -run_type exact -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 -dm_plex_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
1479: test:
1480: suffix: p4est_exact_q2_nonconformal_parallel_parmetis
1481: requires: parmetis p4est
1482: nsize: 4
1483: args: -run_type exact -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 -dm_plex_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
1485: test:
1486: suffix: p4est_full_q2_nonconformal_serial
1487: requires: p4est !single
1488: filter: grep -v "variant HERMITIAN"
1489: args: -run_type full -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 -dm_plex_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
1491: test:
1492: suffix: p4est_full_q2_nonconformal_parallel
1493: requires: p4est !single
1494: filter: grep -v "variant HERMITIAN"
1495: nsize: 7
1496: args: -run_type full -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 -dm_plex_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
1498: test:
1499: suffix: p4est_full_q2_nonconformal_parallel_bddcfas
1500: requires: p4est !single
1501: filter: grep -v "variant HERMITIAN"
1502: nsize: 7
1503: args: -run_type full -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 -dm_plex_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
1505: test:
1506: suffix: p4est_full_q2_nonconformal_parallel_bddc
1507: requires: p4est !single
1508: filter: grep -v "variant HERMITIAN"
1509: nsize: 7
1510: args: -run_type full -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 -dm_plex_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
1512: test:
1513: TODO: broken
1514: suffix: p4est_fas_q2_conformal_serial
1515: requires: p4est !complex !__float128
1516: args: -run_type full -variable_coefficient nonlinear -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 -dm_plex_simplex 0 -dm_refine_hierarchy 3
1518: test:
1519: TODO: broken
1520: suffix: p4est_fas_q2_nonconformal_serial
1521: requires: p4est
1522: args: -run_type full -variable_coefficient nonlinear -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 -dm_plex_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
1524: test:
1525: suffix: fas_newton_0_p4est
1526: requires: p4est !single !__float128
1527: args: -run_type full -variable_coefficient nonlinear -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 -dm_plex_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
1529: # Full solve simplicial AMR
1530: test:
1531: suffix: tri_p1_adapt_init_pragmatic
1532: requires: pragmatic
1533: args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_initial 1 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor pragmatic
1535: test:
1536: suffix: tri_p2_adapt_init_pragmatic
1537: requires: pragmatic
1538: args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_initial 1 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor pragmatic
1540: test:
1541: suffix: tri_p1_adapt_init_mmg
1542: requires: mmg
1543: args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_initial 1 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor mmg
1545: test:
1546: suffix: tri_p2_adapt_init_mmg
1547: requires: mmg
1548: args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_initial 1 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor mmg
1550: test:
1551: suffix: tri_p1_adapt_seq_pragmatic
1552: requires: pragmatic
1553: args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_sequence 2 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor pragmatic
1555: test:
1556: suffix: tri_p2_adapt_seq_pragmatic
1557: requires: pragmatic
1558: args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_sequence 2 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor pragmatic
1560: test:
1561: suffix: tri_p1_adapt_seq_mmg
1562: requires: mmg
1563: args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_sequence 2 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor mmg
1565: test:
1566: suffix: tri_p2_adapt_seq_mmg
1567: requires: mmg
1568: args: -run_type exact -dm_refine 5 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -snes_adapt_sequence 2 -adaptor_target_num 4000 -dm_plex_metric_h_max 0.5 -dm_adaptor mmg
1570: test:
1571: suffix: tri_p1_adapt_analytic_pragmatic
1572: requires: pragmatic
1573: args: -run_type exact -dm_refine 3 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -adaptor_monitor -dm_plex_metric_h_min 0.0001 -dm_plex_metric_h_max 0.05 -dm_adaptor pragmatic
1575: test:
1576: suffix: tri_p2_adapt_analytic_pragmatic
1577: requires: pragmatic
1578: args: -run_type exact -dm_refine 3 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -adaptor_monitor -dm_plex_metric_h_min 0.0001 -dm_plex_metric_h_max 0.05 -dm_adaptor pragmatic
1580: test:
1581: suffix: tri_p1_adapt_analytic_mmg
1582: requires: mmg
1583: args: -run_type exact -dm_refine 3 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -adaptor_monitor -dm_plex_metric_h_max 0.5 -dm_adaptor mmg
1585: test:
1586: suffix: tri_p2_adapt_analytic_mmg
1587: requires: mmg
1588: args: -run_type exact -dm_refine 3 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient cross -snes_adapt_initial 4 -adaptor_target_num 500 -adaptor_monitor -dm_plex_metric_h_max 0.5 -dm_adaptor mmg
1590: test:
1591: suffix: tri_p1_adapt_uniform_pragmatic
1592: requires: pragmatic tetgen
1593: nsize: 2
1594: args: -run_type full -dm_plex_box_faces 8,8,8 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 3 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor pragmatic
1595: timeoutfactor: 2
1597: test:
1598: suffix: tri_p2_adapt_uniform_pragmatic
1599: requires: pragmatic tetgen
1600: nsize: 2
1601: args: -run_type full -dm_plex_box_faces 8,8,8 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 1 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor pragmatic
1602: timeoutfactor: 1
1604: test:
1605: suffix: tri_p1_adapt_uniform_mmg
1606: requires: mmg tetgen
1607: args: -run_type full -dm_plex_box_faces 4,4,4 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 3 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor mmg
1608: timeoutfactor: 2
1610: test:
1611: suffix: tri_p2_adapt_uniform_mmg
1612: requires: mmg tetgen
1613: args: -run_type full -dm_plex_box_faces 4,4,4 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 1 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor mmg
1614: timeoutfactor: 1
1616: test:
1617: suffix: tri_p1_adapt_uniform_parmmg
1618: requires: parmmg tetgen
1619: nsize: 2
1620: args: -run_type full -dm_plex_box_faces 8,8,8 -bc_type dirichlet -petscspace_degree 1 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 3 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor parmmg
1621: timeoutfactor: 2
1623: test:
1624: suffix: tri_p2_adapt_uniform_parmmg
1625: requires: parmmg tetgen
1626: nsize: 2
1627: args: -run_type full -dm_plex_box_faces 8,8,8 -bc_type dirichlet -petscspace_degree 2 -variable_coefficient none -snes_converged_reason ::ascii_info_detail -ksp_type cg -pc_type sor -snes_adapt_sequence 1 -adaptor_target_num 400 -dm_plex_metric_h_max 0.5 -dm_plex_dim 3 -dm_adaptor parmmg
1628: timeoutfactor: 1
1630: # Full solve tensor AMR
1631: test:
1632: suffix: quad_q1_adapt_0
1633: requires: p4est
1634: args: -run_type exact -dm_plex_simplex 0 -dm_plex_convert_type p4est -bc_type dirichlet -petscspace_degree 1 -variable_coefficient ball -snes_converged_reason ::ascii_info_detail -pc_type lu -dm_forest_initial_refinement 4 -snes_adapt_initial 1 -dm_view
1635: filter: grep -v DM_
1637: test:
1638: suffix: amr_0
1639: nsize: 5
1640: args: -run_type test -petscpartitioner_type simple -dm_plex_simplex 0 -bc_type dirichlet -petscspace_degree 1 -dm_refine 1
1642: test:
1643: suffix: amr_1
1644: requires: p4est !complex
1645: args: -run_type test -dm_plex_simplex 0 -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
1647: test:
1648: suffix: p4est_solve_bddc
1649: requires: p4est !complex
1650: args: -run_type full -variable_coefficient nonlinear -nonzero_initial_guess 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 -dm_plex_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
1651: nsize: 4
1653: test:
1654: suffix: p4est_solve_fas
1655: requires: p4est
1656: args: -run_type full -variable_coefficient nonlinear -nonzero_initial_guess 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 -dm_plex_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
1657: nsize: 4
1658: TODO: identical machine two runs produce slightly different solver trackers
1660: test:
1661: suffix: p4est_convergence_test_1
1662: requires: p4est
1663: args: -quiet -run_type test -petscspace_degree 1 -dm_plex_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
1664: nsize: 4
1666: test:
1667: suffix: p4est_convergence_test_2
1668: requires: p4est
1669: args: -quiet -run_type test -petscspace_degree 1 -dm_plex_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
1671: test:
1672: suffix: p4est_convergence_test_3
1673: requires: p4est
1674: args: -quiet -run_type test -petscspace_degree 1 -dm_plex_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
1676: test:
1677: suffix: p4est_convergence_test_4
1678: requires: p4est
1679: args: -quiet -run_type test -petscspace_degree 1 -dm_plex_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
1680: timeoutfactor: 5
1682: # Serial tests with GLVis visualization
1683: test:
1684: suffix: glvis_2d_tet_p1
1685: args: -quiet -run_type test -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_plex_boundary_label marker -dm_plex_gmsh_periodic 0 -dm_coord_space 0
1686: test:
1687: suffix: glvis_2d_tet_p2
1688: args: -quiet -run_type test -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_plex_boundary_label marker -dm_plex_gmsh_periodic 0 -dm_coord_space 0
1689: test:
1690: suffix: glvis_2d_hex_p1
1691: args: -quiet -run_type test -bc_type dirichlet -petscspace_degree 1 -vec_view glvis: -dm_plex_simplex 0 -dm_refine 1 -dm_coord_space 0
1692: test:
1693: suffix: glvis_2d_hex_p2
1694: args: -quiet -run_type test -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -dm_plex_simplex 0 -dm_refine 1 -dm_coord_space 0
1695: test:
1696: suffix: glvis_2d_hex_p2_p4est
1697: requires: p4est
1698: args: -quiet -run_type test -bc_type dirichlet -petscspace_degree 2 -vec_view glvis: -dm_plex_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 -viewer_glvis_dm_plex_enable_ncmesh
1699: test:
1700: suffix: glvis_2d_tet_p0
1701: args: -run_type exact -guess_vec_view glvis: -nonzero_initial_guess 1 -dm_plex_filename ${wPETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_plex_boundary_label marker -petscspace_degree 0 -dm_coord_space 0
1702: test:
1703: suffix: glvis_2d_hex_p0
1704: args: -run_type exact -guess_vec_view glvis: -nonzero_initial_guess 1 -dm_plex_box_faces 5,7 -dm_plex_simplex 0 -petscspace_degree 0 -dm_coord_space 0
1706: # PCHPDDM tests
1707: testset:
1708: nsize: 4
1709: requires: hpddm slepc !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1710: args: -run_type test -run_test_check_ksp -quiet -petscspace_degree 1 -petscpartitioner_type simple -bc_type none -dm_plex_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
1711: test:
1712: suffix: quad_singular_hpddm
1713: args: -dm_plex_box_faces 6,7
1714: test:
1715: requires: p4est
1716: suffix: p4est_singular_2d_hpddm
1717: args: -dm_plex_convert_type p4est -dm_forest_minimum_refinement 1 -dm_forest_initial_refinement 3 -dm_forest_maximum_refinement 3
1718: test:
1719: requires: p4est
1720: suffix: p4est_nc_singular_2d_hpddm
1721: 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
1722: testset:
1723: nsize: 4
1724: requires: hpddm slepc triangle !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1725: args: -run_type full -petscpartitioner_type simple -dm_refine 2 -bc_type dirichlet -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
1726: test:
1727: args: -pc_hpddm_coarse_mat_type baij -options_left no
1728: suffix: tri_hpddm_reuse_baij
1729: test:
1730: requires: !complex
1731: suffix: tri_hpddm_reuse
1732: testset:
1733: nsize: 4
1734: requires: hpddm slepc !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1735: args: -run_type full -petscpartitioner_type simple -dm_plex_box_faces 7,5 -dm_refine 2 -dm_plex_simplex 0 -bc_type dirichlet -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
1736: test:
1737: args: -pc_hpddm_coarse_mat_type baij -options_left no
1738: suffix: quad_hpddm_reuse_baij
1739: test:
1740: requires: !complex
1741: suffix: quad_hpddm_reuse
1742: testset:
1743: nsize: 4
1744: requires: hpddm slepc !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1745: args: -run_type full -petscpartitioner_type simple -dm_plex_box_faces 7,5 -dm_refine 2 -dm_plex_simplex 0 -bc_type dirichlet -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
1746: test:
1747: args: -pc_hpddm_coarse_mat_type baij -options_left no
1748: suffix: quad_hpddm_reuse_threshold_baij
1749: test:
1750: requires: !complex
1751: suffix: quad_hpddm_reuse_threshold
1752: testset:
1753: nsize: 4
1754: requires: hpddm slepc parmetis !single defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
1755: filter: sed -e "s/linear solver iterations=17/linear solver iterations=16/g"
1756: args: -run_type full -petscpartitioner_type parmetis -dm_refine 3 -bc_type dirichlet -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 -dm_plex_filename ${PETSC_DIR}/share/petsc/datafiles/meshes/square_periodic.msh -dm_plex_boundary_label marker -pc_hpddm_levels_1_sub_pc_factor_levels 3 -variable_coefficient ball -dm_plex_gmsh_periodic 0
1757: test:
1758: args: -pc_hpddm_coarse_mat_type baij -options_left no
1759: 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"
1760: suffix: tri_parmetis_hpddm_baij
1761: test:
1762: 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"
1763: requires: !complex
1764: suffix: tri_parmetis_hpddm
1766: # 2D serial P1 tests for adaptive MG
1767: test:
1768: suffix: 2d_p1_adaptmg_0
1769: requires: triangle
1770: args: -petscpartitioner_type simple -dm_refine_hierarchy 3 -dm_plex_box_faces 4,4 -bc_type dirichlet -petscspace_degree 1 \
1771: -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
1772: -snes_max_it 1 -ksp_converged_reason \
1773: -ksp_rtol 1e-8 -pc_type mg
1774: test:
1775: suffix: 2d_p1_adaptmg_1
1776: requires: triangle bamg todo
1777: args: -petscpartitioner_type simple -dm_refine_hierarchy 3 -dm_plex_box_faces 4,4 -bc_type dirichlet -petscspace_degree 1 \
1778: -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
1779: -snes_max_it 1 -ksp_converged_reason \
1780: -ksp_rtol 1e-8 -pc_type mg -pc_mg_galerkin -pc_mg_adapt_interp_coarse_space eigenvector -pc_mg_adapt_interp_n 1 \
1781: -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
1782: test:
1783: suffix: 2d_p1_adaptmg_gdsw
1784: requires: triangle
1785: nsize: 4
1786: args: -petscpartitioner_type simple -dm_refine 3 -dm_plex_box_faces 4,4 -bc_type dirichlet -petscspace_degree 1 \
1787: -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
1788: -snes_max_it 1 -ksp_converged_reason \
1789: -ksp_rtol 1e-8 -pc_type mg -pc_mg_galerkin -pc_mg_adapt_interp_coarse_space gdsw -pc_mg_levels 2 -mg_levels_pc_type asm -dm_mat_type {{aij is}}
1791: test:
1792: suffix: 2d_p1_adaptmg_agdsw
1793: requires: triangle mumps
1794: nsize: 4
1795: args: -petscpartitioner_type simple -dm_refine 3 -dm_plex_box_faces 4,4 -bc_type dirichlet -petscspace_degree 1 \
1796: -variable_coefficient checkerboard_0 -mat_petscspace_degree 0 -div 16 -k 3 \
1797: -snes_max_it 1 -ksp_converged_reason \
1798: -ksp_rtol 1e-8 -pc_type mg -pc_mg_galerkin -pc_mg_adapt_interp_coarse_space gdsw -pc_mg_levels 2 -mg_levels_pc_type asm -dm_mat_type is -mg_levels_gdsw_tolerance 0.1 -mg_levels_gdsw_pseudo_pc_type qr
1800: TEST*/