Actual source code: ex62.c
petsc-3.10.5 2019-03-28
1: static char help[] = "Stokes Problem in 2d and 3d with simplicial finite elements.\n\
2: We solve the Stokes problem in a rectangular\n\
3: domain, using a parallel unstructured mesh (DMPLEX) to discretize it.\n\n\n";
5: /*
6: The isoviscous Stokes problem, which we discretize using the finite
7: element method on an unstructured mesh. The weak form equations are
9: < \nabla v, \nabla u + {\nabla u}^T > - < \nabla\cdot v, p > + < v, f > = 0
10: < q, \nabla\cdot u > = 0
12: Viewing:
14: To produce nice output, use
16: -dm_refine 3 -show_error -dm_view hdf5:sol1.h5 -error_vec_view hdf5:sol1.h5::append -sol_vec_view hdf5:sol1.h5::append -exact_vec_view hdf5:sol1.h5::append
18: You can get a LaTeX view of the mesh, with point numbering using
20: -dm_view :mesh.tex:ascii_latex -dm_plex_view_scale 8.0
22: The data layout can be viewed using
24: -dm_petscsection_view
26: Lots of information about the FEM assembly can be printed using
28: -dm_plex_print_fem 2
30: Field Data:
32: DMPLEX data is organized by point, and the closure operation just stacks up the
33: data from each sieve point in the closure. Thus, for a P_2-P_1 Stokes element, we
34: have
36: cl{e} = {f e_0 e_1 e_2 v_0 v_1 v_2}
37: x = [u_{e_0} v_{e_0} u_{e_1} v_{e_1} u_{e_2} v_{e_2} u_{v_0} v_{v_0} p_{v_0} u_{v_1} v_{v_1} p_{v_1} u_{v_2} v_{v_2} p_{v_2}]
39: The problem here is that we would like to loop over each field separately for
40: integration. Therefore, the closure visitor in DMPlexVecGetClosure() reorders
41: the data so that each field is contiguous
43: x' = [u_{e_0} v_{e_0} u_{e_1} v_{e_1} u_{e_2} v_{e_2} u_{v_0} v_{v_0} u_{v_1} v_{v_1} u_{v_2} v_{v_2} p_{v_0} p_{v_1} p_{v_2}]
45: Likewise, DMPlexVecSetClosure() takes data partitioned by field, and correctly
46: puts it into the Sieve ordering.
48: TODO:
49: - Update FETI test output
50: - Reorder mesh
51: - Check the q1-p0 Vanka domains are correct (I think its correct)
52: - Check scaling of iterates, right now it is bad
53: - Check the q2-q1 domains since convergence is bad
54: - Ask Patrick about domains
55: - Plot residual by fields after each smoother iterate
56: - Get Diskin checks going
57: */
59: #include <petscdmplex.h>
60: #include <petscsnes.h>
61: #include <petscds.h>
63: PetscInt spatialDim = 0;
65: typedef enum {NEUMANN, DIRICHLET, NUM_BC_TYPES} BCType;
66: const char *bcTypes[NUM_BC_TYPES+1] = {"neumann", "dirichlet", "unknown"};
67: typedef enum {RUN_FULL, RUN_TEST, NUM_RUN_TYPES} RunType;
68: const char *runTypes[NUM_RUN_TYPES+1] = {"full", "test", "unknown"};
69: typedef enum {SOL_QUADRATIC, SOL_CUBIC, SOL_TRIG, NUM_SOL_TYPES} SolType;
70: const char *solTypes[NUM_SOL_TYPES+1] = {"quadratic", "cubic", "trig", "unknown"};
72: typedef struct {
73: PetscInt debug; /* The debugging level */
74: RunType runType; /* Whether to run tests, or solve the full problem */
75: PetscLogEvent createMeshEvent;
76: PetscBool showInitial, showError;
77: /* Domain and mesh definition */
78: PetscInt dim; /* The topological mesh dimension */
79: PetscBool interpolate; /* Generate intermediate mesh elements */
80: PetscBool simplex; /* Use simplices or tensor product cells */
81: PetscReal refinementLimit; /* The largest allowable cell volume */
82: PetscBool testPartition; /* Use a fixed partitioning for testing */
83: /* Problem definition */
84: BCType bcType;
85: SolType solType;
86: PetscErrorCode (**exactFuncs)(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx);
87: } AppCtx;
89: PetscErrorCode zero_scalar(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx)
90: {
91: u[0] = 0.0;
92: return 0;
93: }
94: PetscErrorCode zero_vector(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx)
95: {
96: PetscInt d;
97: for (d = 0; d < spatialDim; ++d) u[d] = 0.0;
98: return 0;
99: }
101: /*
102: In 2D we use exact solution:
104: u = x^2 + y^2
105: v = 2 x^2 - 2xy
106: p = x + y - 1
107: f_x = f_y = 3
109: so that
111: -\Delta u + \nabla p + f = <-4, -4> + <1, 1> + <3, 3> = 0
112: \nabla \cdot u = 2x - 2x = 0
113: */
114: PetscErrorCode quadratic_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx)
115: {
116: u[0] = x[0]*x[0] + x[1]*x[1];
117: u[1] = 2.0*x[0]*x[0] - 2.0*x[0]*x[1];
118: return 0;
119: }
121: PetscErrorCode linear_p_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *p, void *ctx)
122: {
123: *p = x[0] + x[1] - 1.0;
124: return 0;
125: }
126: PetscErrorCode constant_p(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *p, void *ctx)
127: {
128: *p = 1.0;
129: return 0;
130: }
132: /*
133: In 2D we use exact solution:
135: u = x^3 + y^3
136: v = 2 x^3 - 3 x^2 y
137: p = 3/2 x^2 + 3/2 y^2 - 1
138: f_x = 6 (x + y)
139: f_y = 12 x - 3 y
141: so that
143: -\Delta u + \nabla p + f = <-6 x - 6 y, -12 x + 6 y> + <3 x, 3 y> + <6 (x + y), 12 x - 6 y> = 0
144: \nabla \cdot u = 3 x^2 - 3 x^2 = 0
145: */
146: PetscErrorCode cubic_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx)
147: {
148: u[0] = x[0]*x[0]*x[0] + x[1]*x[1]*x[1];
149: u[1] = 2.0*x[0]*x[0]*x[0] - 3.0*x[0]*x[0]*x[1];
150: return 0;
151: }
153: PetscErrorCode quadratic_p_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *p, void *ctx)
154: {
155: *p = 3.0*x[0]*x[0]/2.0 + 3.0*x[1]*x[1]/2.0 - 1.0;
156: return 0;
157: }
159: /*
160: In 2D we use exact solution:
162: u = sin(n pi x) + y^2
163: v = -sin(n pi y)
164: p = 3/2 x^2 + 3/2 y^2 - 1
165: f_x = 4 - 3x - n^2 pi^2 sin (n pi x)
166: f_y = - 3y + n^2 pi^2 sin(n pi y)
168: so that
170: -\Delta u + \nabla p + f = <n^2 pi^2 sin (n pi x) - 4, -n^2 pi^2 sin(n pi y)> + <3 x, 3 y> + <4 - 3x - n^2 pi^2 sin (n pi x), -3y + n^2 pi^2 sin(n pi y)> = 0
171: \nabla \cdot u = n pi cos(n pi x) - n pi cos(n pi y) = 0
172: */
173: PetscErrorCode trig_u_2d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx)
174: {
175: const PetscReal n = 1.0;
177: u[0] = PetscSinReal(n*PETSC_PI*x[0]) + x[1]*x[1];
178: u[1] = -PetscSinReal(n*PETSC_PI*x[1]);
179: return 0;
180: }
182: void f0_quadratic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
183: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
184: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
185: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
186: {
187: PetscInt c;
188: for (c = 0; c < dim; ++c) f0[c] = 3.0;
189: }
191: void f0_cubic_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
192: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
193: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
194: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
195: {
196: f0[0] = 3.0*x[0] + 6.0*x[1];
197: f0[1] = 12.0*x[0] - 9.0*x[1];
198: }
200: void f0_trig_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
201: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
202: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
203: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
204: {
205: const PetscReal n = 1.0;
207: f0[0] = 4.0 - 3.0*x[0] - PetscSqr(n*PETSC_PI)*PetscSinReal(n*PETSC_PI*x[0]);
208: f0[1] = -3.0*x[1] + PetscSqr(n*PETSC_PI)*PetscSinReal(n*PETSC_PI*x[1]);
209: }
211: /* gradU[comp*dim+d] = {u_x, u_y, v_x, v_y} or {u_x, u_y, u_z, v_x, v_y, v_z, w_x, w_y, w_z}
212: u[Ncomp] = {p} */
213: void f1_u(PetscInt dim, PetscInt Nf, PetscInt NfAux,
214: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
215: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
216: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
217: {
218: const PetscInt Ncomp = dim;
219: PetscInt comp, d;
221: for (comp = 0; comp < Ncomp; ++comp) {
222: for (d = 0; d < dim; ++d) {
223: /* f1[comp*dim+d] = 0.5*(gradU[comp*dim+d] + gradU[d*dim+comp]); */
224: f1[comp*dim+d] = u_x[comp*dim+d];
225: }
226: f1[comp*dim+comp] -= u[Ncomp];
227: }
228: }
230: /* gradU[comp*dim+d] = {u_x, u_y, v_x, v_y} or {u_x, u_y, u_z, v_x, v_y, v_z, w_x, w_y, w_z} */
231: void f0_p(PetscInt dim, PetscInt Nf, PetscInt NfAux,
232: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
233: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
234: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
235: {
236: PetscInt d;
237: for (d = 0, f0[0] = 0.0; d < dim; ++d) f0[0] += u_x[d*dim+d];
238: }
240: void f1_p(PetscInt dim, PetscInt Nf, PetscInt NfAux,
241: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
242: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
243: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[])
244: {
245: PetscInt d;
246: for (d = 0; d < dim; ++d) f1[d] = 0.0;
247: }
249: /* < q, \nabla\cdot u >
250: NcompI = 1, NcompJ = dim */
251: void g1_pu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
252: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
253: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
254: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[])
255: {
256: PetscInt d;
257: for (d = 0; d < dim; ++d) g1[d*dim+d] = 1.0; /* \frac{\partial\phi^{u_d}}{\partial x_d} */
258: }
260: /* -< \nabla\cdot v, p >
261: NcompI = dim, NcompJ = 1 */
262: void g2_up(PetscInt dim, PetscInt Nf, PetscInt NfAux,
263: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
264: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
265: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[])
266: {
267: PetscInt d;
268: for (d = 0; d < dim; ++d) g2[d*dim+d] = -1.0; /* \frac{\partial\psi^{u_d}}{\partial x_d} */
269: }
271: /* < \nabla v, \nabla u + {\nabla u}^T >
272: This just gives \nabla u, give the perdiagonal for the transpose */
273: void g3_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux,
274: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
275: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
276: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[])
277: {
278: const PetscInt Ncomp = dim;
279: PetscInt compI, d;
281: for (compI = 0; compI < Ncomp; ++compI) {
282: for (d = 0; d < dim; ++d) {
283: g3[((compI*Ncomp+compI)*dim+d)*dim+d] = 1.0;
284: }
285: }
286: }
288: /*
289: In 3D we use exact solution:
291: u = x^2 + y^2
292: v = y^2 + z^2
293: w = x^2 + y^2 - 2(x+y)z
294: p = x + y + z - 3/2
295: f_x = f_y = f_z = 3
297: so that
299: -\Delta u + \nabla p + f = <-4, -4, -4> + <1, 1, 1> + <3, 3, 3> = 0
300: \nabla \cdot u = 2x + 2y - 2(x + y) = 0
301: */
302: PetscErrorCode quadratic_u_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void *ctx)
303: {
304: u[0] = x[0]*x[0] + x[1]*x[1];
305: u[1] = x[1]*x[1] + x[2]*x[2];
306: u[2] = x[0]*x[0] + x[1]*x[1] - 2.0*(x[0] + x[1])*x[2];
307: return 0;
308: }
310: PetscErrorCode linear_p_3d(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *p, void *ctx)
311: {
312: *p = x[0] + x[1] + x[2] - 1.5;
313: return 0;
314: }
316: void pressure(PetscInt dim, PetscInt Nf, PetscInt NfAux,
317: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
318: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
319: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar p[])
320: {
321: p[0] = u[uOff[1]];
322: }
324: PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
325: {
326: PetscInt bc, run, sol;
330: options->debug = 0;
331: options->runType = RUN_FULL;
332: options->dim = 2;
333: options->interpolate = PETSC_FALSE;
334: options->simplex = PETSC_TRUE;
335: options->refinementLimit = 0.0;
336: options->testPartition = PETSC_FALSE;
337: options->bcType = DIRICHLET;
338: options->solType = SOL_QUADRATIC;
339: options->showInitial = PETSC_FALSE;
340: options->showError = PETSC_FALSE;
342: PetscOptionsBegin(comm, "", "Stokes Problem Options", "DMPLEX");
343: PetscOptionsInt("-debug", "The debugging level", "ex62.c", options->debug, &options->debug, NULL);
344: run = options->runType;
345: PetscOptionsEList("-run_type", "The run type", "ex62.c", runTypes, NUM_RUN_TYPES, runTypes[options->runType], &run, NULL);
346: options->runType = (RunType) run;
347: PetscOptionsInt("-dim", "The topological mesh dimension", "ex62.c", options->dim, &options->dim, NULL);
348: spatialDim = options->dim;
349: PetscOptionsBool("-interpolate", "Generate intermediate mesh elements", "ex62.c", options->interpolate, &options->interpolate, NULL);
350: PetscOptionsBool("-simplex", "Use simplices or tensor product cells", "ex62.c", options->simplex, &options->simplex, NULL);
351: PetscOptionsReal("-refinement_limit", "The largest allowable cell volume", "ex62.c", options->refinementLimit, &options->refinementLimit, NULL);
352: PetscOptionsBool("-test_partition", "Use a fixed partition for testing", "ex62.c", options->testPartition, &options->testPartition, NULL);
353: bc = options->bcType;
354: PetscOptionsEList("-bc_type","Type of boundary condition","ex62.c", bcTypes, NUM_BC_TYPES, bcTypes[options->bcType], &bc, NULL);
355: options->bcType = (BCType) bc;
356: sol = options->solType;
357: PetscOptionsEList("-sol_type", "The solution type", "ex62.c", solTypes, NUM_SOL_TYPES, solTypes[options->solType], &sol, NULL);
358: options->solType = (SolType) sol;
359: PetscOptionsBool("-show_initial", "Output the initial guess for verification", "ex62.c", options->showInitial, &options->showInitial, NULL);
360: PetscOptionsBool("-show_error", "Output the error for verification", "ex62.c", options->showError, &options->showError, NULL);
361: PetscOptionsEnd();
363: PetscLogEventRegister("CreateMesh", DM_CLASSID, &options->createMeshEvent);
364: return(0);
365: }
367: PetscErrorCode DMVecViewLocal(DM dm, Vec v)
368: {
369: Vec lv;
373: DMGetLocalVector(dm, &lv);
374: DMGlobalToLocalBegin(dm, v, INSERT_VALUES, lv);
375: DMGlobalToLocalEnd(dm, v, INSERT_VALUES, lv);
376: DMPrintLocalVec(dm, "Local function", 0.0, lv);
377: DMRestoreLocalVector(dm, &lv);
378: return(0);
379: }
381: PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
382: {
383: PetscInt dim = user->dim;
384: PetscBool interpolate = user->interpolate;
385: PetscReal refinementLimit = user->refinementLimit;
386: const PetscInt cells[3] = {3, 3, 3};
390: PetscLogEventBegin(user->createMeshEvent,0,0,0,0);
391: DMPlexCreateBoxMesh(comm, dim, user->simplex, user->simplex ? NULL : cells, NULL, NULL, NULL, interpolate, dm);
392: {
393: DM refinedMesh = NULL;
394: DM distributedMesh = NULL;
396: /* Refine mesh using a volume constraint */
397: DMPlexSetRefinementLimit(*dm, refinementLimit);
398: if (user->simplex) {DMRefine(*dm, comm, &refinedMesh);}
399: if (refinedMesh) {
400: DMDestroy(dm);
401: *dm = refinedMesh;
402: }
403: /* Setup test partitioning */
404: if (user->testPartition) {
405: PetscInt triSizes_n2[2] = {4, 4};
406: PetscInt triPoints_n2[8] = {3, 5, 6, 7, 0, 1, 2, 4};
407: PetscInt triSizes_n3[3] = {2, 3, 3};
408: PetscInt triPoints_n3[8] = {3, 5, 1, 6, 7, 0, 2, 4};
409: PetscInt triSizes_n5[5] = {1, 2, 2, 1, 2};
410: PetscInt triPoints_n5[8] = {3, 5, 6, 4, 7, 0, 1, 2};
411: PetscInt triSizes_ref_n2[2] = {8, 8};
412: PetscInt triPoints_ref_n2[16] = {1, 5, 6, 7, 10, 11, 14, 15, 0, 2, 3, 4, 8, 9, 12, 13};
413: PetscInt triSizes_ref_n3[3] = {5, 6, 5};
414: PetscInt triPoints_ref_n3[16] = {1, 7, 10, 14, 15, 2, 6, 8, 11, 12, 13, 0, 3, 4, 5, 9};
415: PetscInt triSizes_ref_n5[5] = {3, 4, 3, 3, 3};
416: PetscInt triPoints_ref_n5[16] = {1, 7, 10, 2, 11, 13, 14, 5, 6, 15, 0, 8, 9, 3, 4, 12};
417: PetscInt triSizes_ref_n5_d3[5] = {1, 1, 1, 1, 2};
418: PetscInt triPoints_ref_n5_d3[6] = {0, 1, 2, 3, 4, 5};
419: const PetscInt *sizes = NULL;
420: const PetscInt *points = NULL;
421: PetscPartitioner part;
422: PetscInt cEnd;
423: PetscMPIInt rank, size;
425: MPI_Comm_rank(comm, &rank);
426: MPI_Comm_size(comm, &size);
427: DMPlexGetHeightStratum(*dm, 0, NULL, &cEnd);
428: if (!rank) {
429: if (dim == 2 && user->simplex && size == 2 && cEnd == 8) {
430: sizes = triSizes_n2; points = triPoints_n2;
431: } else if (dim == 2 && user->simplex && size == 3 && cEnd == 8) {
432: sizes = triSizes_n3; points = triPoints_n3;
433: } else if (dim == 2 && user->simplex && size == 5 && cEnd == 8) {
434: sizes = triSizes_n5; points = triPoints_n5;
435: } else if (dim == 2 && user->simplex && size == 2 && cEnd == 16) {
436: sizes = triSizes_ref_n2; points = triPoints_ref_n2;
437: } else if (dim == 2 && user->simplex && size == 3 && cEnd == 16) {
438: sizes = triSizes_ref_n3; points = triPoints_ref_n3;
439: } else if (dim == 2 && user->simplex && size == 5 && cEnd == 16) {
440: sizes = triSizes_ref_n5; points = triPoints_ref_n5;
441: } else if (dim == 3 && user->simplex && size == 5 && cEnd == 6) {
442: sizes = triSizes_ref_n5_d3; points = triPoints_ref_n5_d3;
443: } else SETERRQ(comm, PETSC_ERR_ARG_WRONG, "No stored partition matching run parameters");
444: }
445: DMPlexGetPartitioner(*dm, &part);
446: PetscPartitionerSetType(part, PETSCPARTITIONERSHELL);
447: PetscPartitionerShellSetPartition(part, size, sizes, points);
448: } else {
449: PetscPartitioner part;
451: DMPlexGetPartitioner(*dm, &part);
452: PetscPartitionerSetFromOptions(part);
453: }
454: /* Distribute mesh over processes */
455: DMPlexDistribute(*dm, 0, NULL, &distributedMesh);
456: if (distributedMesh) {
457: DMDestroy(dm);
458: *dm = distributedMesh;
459: }
460: }
461: DMSetFromOptions(*dm);
462: DMViewFromOptions(*dm, NULL, "-dm_view");
463: PetscLogEventEnd(user->createMeshEvent,0,0,0,0);
464: return(0);
465: }
467: PetscErrorCode SetupProblem(PetscDS prob, AppCtx *user)
468: {
469: const PetscInt id = 1;
473: switch (user->solType) {
474: case SOL_QUADRATIC:
475: switch (user->dim) {
476: case 2:
477: PetscDSSetResidual(prob, 0, f0_quadratic_u, f1_u);
478: user->exactFuncs[0] = quadratic_u_2d;
479: user->exactFuncs[1] = linear_p_2d;
480: break;
481: case 3:
482: PetscDSSetResidual(prob, 0, f0_quadratic_u, f1_u);
483: user->exactFuncs[0] = quadratic_u_3d;
484: user->exactFuncs[1] = linear_p_3d;
485: break;
486: default: SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Unsupported dimension %d for quadratic solution", user->dim);
487: }
488: break;
489: case SOL_CUBIC:
490: switch (user->dim) {
491: case 2:
492: PetscDSSetResidual(prob, 0, f0_cubic_u, f1_u);
493: user->exactFuncs[0] = cubic_u_2d;
494: user->exactFuncs[1] = quadratic_p_2d;
495: break;
496: default: SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Unsupported dimension %d for quadratic solution", user->dim);
497: }
498: break;
499: case SOL_TRIG:
500: switch (user->dim) {
501: case 2:
502: PetscDSSetResidual(prob, 0, f0_trig_u, f1_u);
503: user->exactFuncs[0] = trig_u_2d;
504: user->exactFuncs[1] = quadratic_p_2d;
505: break;
506: default: SETERRQ1(PETSC_COMM_WORLD, PETSC_ERR_ARG_OUTOFRANGE, "Unsupported dimension %d for trigonometric solution", user->dim);
507: }
508: break;
509: default: SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unsupported solution type: %s (%D)", solTypes[PetscMin(user->solType, NUM_SOL_TYPES)], user->solType);
510: }
511: PetscDSSetResidual(prob, 1, f0_p, f1_p);
512: PetscDSSetJacobian(prob, 0, 0, NULL, NULL, NULL, g3_uu);
513: PetscDSSetJacobian(prob, 0, 1, NULL, NULL, g2_up, NULL);
514: PetscDSSetJacobian(prob, 1, 0, NULL, g1_pu, NULL, NULL);
516: PetscDSAddBoundary(prob, user->bcType == DIRICHLET ? DM_BC_ESSENTIAL : DM_BC_NATURAL, "wall", user->bcType == NEUMANN ? "boundary" : "marker", 0, 0, NULL, (void (*)(void)) user->exactFuncs[0], 1, &id, user);
517: PetscDSSetExactSolution(prob, 0, user->exactFuncs[0]);
518: PetscDSSetExactSolution(prob, 1, user->exactFuncs[1]);
519: PetscDSSetFromOptions(prob);
520: return(0);
521: }
523: PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
524: {
525: DM cdm = dm;
526: const PetscInt dim = user->dim;
527: PetscFE fe[2];
528: PetscQuadrature q;
529: PetscDS prob;
530: MPI_Comm comm;
531: PetscErrorCode ierr;
534: /* Create finite element */
535: PetscObjectGetComm((PetscObject) dm, &comm);
536: PetscFECreateDefault(comm, dim, dim, user->simplex, "vel_", PETSC_DEFAULT, &fe[0]);
537: PetscObjectSetName((PetscObject) fe[0], "velocity");
538: PetscFEGetQuadrature(fe[0], &q);
539: PetscFECreateDefault(comm, dim, 1, user->simplex, "pres_", PETSC_DEFAULT, &fe[1]);
540: PetscFESetQuadrature(fe[1], q);
541: PetscObjectSetName((PetscObject) fe[1], "pressure");
542: /* Set discretization and boundary conditions for each mesh */
543: DMGetDS(dm, &prob);
544: PetscDSSetDiscretization(prob, 0, (PetscObject) fe[0]);
545: PetscDSSetDiscretization(prob, 1, (PetscObject) fe[1]);
546: SetupProblem(prob, user);
547: while (cdm) {
548: DMSetDS(cdm, prob);
549: DMGetCoarseDM(cdm, &cdm);
550: }
551: PetscFEDestroy(&fe[0]);
552: PetscFEDestroy(&fe[1]);
553: return(0);
554: }
556: static PetscErrorCode CreatePressureNullSpace(DM dm, PetscInt dummy, MatNullSpace *nullspace)
557: {
558: Vec vec;
559: PetscErrorCode (*funcs[2])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void* ctx) = {zero_vector, constant_p};
560: PetscErrorCode ierr;
563: DMCreateGlobalVector(dm, &vec);
564: DMProjectFunction(dm, 0.0, funcs, NULL, INSERT_ALL_VALUES, vec);
565: VecNormalize(vec, NULL);
566: PetscObjectSetName((PetscObject) vec, "Pressure Null Space");
567: VecViewFromOptions(vec, NULL, "-pressure_nullspace_view");
568: MatNullSpaceCreate(PetscObjectComm((PetscObject)dm), PETSC_FALSE, 1, &vec, nullspace);
569: VecDestroy(&vec);
570: /* New style for field null spaces */
571: {
572: PetscObject pressure;
573: MatNullSpace nullspacePres;
575: DMGetField(dm, 1, &pressure);
576: MatNullSpaceCreate(PetscObjectComm(pressure), PETSC_TRUE, 0, NULL, &nullspacePres);
577: PetscObjectCompose(pressure, "nullspace", (PetscObject) nullspacePres);
578: MatNullSpaceDestroy(&nullspacePres);
579: }
580: return(0);
581: }
583: /* Add a vector in the nullspace to make the continuum integral 0.
585: If int(u) = a and int(n) = b, then int(u - a/b n) = a - a/b b = 0
586: */
587: static PetscErrorCode CorrectDiscretePressure(DM dm, MatNullSpace nullspace, Vec u, AppCtx *user)
588: {
589: PetscDS prob;
590: const Vec *nullvecs;
591: PetscScalar pintd, intc[2], intn[2];
592: MPI_Comm comm;
596: PetscObjectGetComm((PetscObject) dm, &comm);
597: DMGetDS(dm, &prob);
598: PetscDSSetObjective(prob, 1, pressure);
599: MatNullSpaceGetVecs(nullspace, NULL, NULL, &nullvecs);
600: VecDot(nullvecs[0], u, &pintd);
601: if (PetscAbsScalar(pintd) > 1.0e-10) SETERRQ1(comm, PETSC_ERR_ARG_WRONG, "Discrete integral of pressure: %g\n", (double) PetscRealPart(pintd));
602: DMPlexComputeIntegralFEM(dm, nullvecs[0], intn, user);
603: DMPlexComputeIntegralFEM(dm, u, intc, user);
604: VecAXPY(u, -intc[1]/intn[1], nullvecs[0]);
605: DMPlexComputeIntegralFEM(dm, u, intc, user);
606: if (PetscAbsScalar(intc[1]) > 1.0e-10) SETERRQ1(comm, PETSC_ERR_ARG_WRONG, "Continuum integral of pressure after correction: %g\n", (double) PetscRealPart(intc[1]));
607: return(0);
608: }
610: static PetscErrorCode SNESConvergenceCorrectPressure(SNES snes, PetscInt it, PetscReal xnorm, PetscReal gnorm, PetscReal f, SNESConvergedReason *reason, void *user)
611: {
615: SNESConvergedDefault(snes, it, xnorm, gnorm, f, reason, user);
616: if (*reason > 0) {
617: DM dm;
618: Mat J;
619: Vec u;
620: MatNullSpace nullspace;
622: SNESGetDM(snes, &dm);
623: SNESGetSolution(snes, &u);
624: SNESGetJacobian(snes, &J, NULL, NULL, NULL);
625: MatGetNullSpace(J, &nullspace);
626: CorrectDiscretePressure(dm, nullspace, u, (AppCtx *) user);
627: }
628: return(0);
629: }
631: int main(int argc, char **argv)
632: {
633: SNES snes; /* nonlinear solver */
634: DM dm; /* problem definition */
635: Vec u, r; /* solution and residual */
636: AppCtx user; /* user-defined work context */
637: PetscReal error = 0.0; /* L_2 error in the solution */
638: PetscReal ferrors[2];
641: PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr;
642: ProcessOptions(PETSC_COMM_WORLD, &user);
643: SNESCreate(PETSC_COMM_WORLD, &snes);
644: CreateMesh(PETSC_COMM_WORLD, &user, &dm);
645: SNESSetDM(snes, dm);
646: DMSetApplicationContext(dm, &user);
648: PetscMalloc(2 * sizeof(void (*)(const PetscReal[], PetscScalar *, void *)), &user.exactFuncs);
649: SetupDiscretization(dm, &user);
650: DMPlexCreateClosureIndex(dm, NULL);
652: DMCreateGlobalVector(dm, &u);
653: VecDuplicate(u, &r);
655: DMSetNullSpaceConstructor(dm, 2, CreatePressureNullSpace);
656: SNESSetConvergenceTest(snes, SNESConvergenceCorrectPressure, &user, NULL);
658: DMPlexSetSNESLocalFEM(dm,&user,&user,&user);
660: SNESSetFromOptions(snes);
662: DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, u);
663: PetscObjectSetName((PetscObject) u, "Exact Solution");
664: VecViewFromOptions(u, NULL, "-exact_vec_view");
665: PetscObjectSetName((PetscObject) u, "Solution");
666: if (user.showInitial) {DMVecViewLocal(dm, u);}
667: PetscObjectSetName((PetscObject) u, "Solution");
668: if (user.runType == RUN_FULL) {
669: PetscErrorCode (*initialGuess[2])(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nf, PetscScalar *u, void* ctx) = {zero_vector, zero_scalar};
671: DMProjectFunction(dm, 0.0, initialGuess, NULL, INSERT_VALUES, u);
672: if (user.debug) {
673: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
674: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
675: }
676: SNESSolve(snes, NULL, u);
677: DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
678: DMComputeL2FieldDiff(dm, 0.0, user.exactFuncs, NULL, u, ferrors);
679: PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g [%g, %g]\n", error, ferrors[0], ferrors[1]);
680: if (user.showError) {
681: Vec r;
683: DMGetGlobalVector(dm, &r);
684: DMProjectFunction(dm, 0.0, user.exactFuncs, NULL, INSERT_ALL_VALUES, r);
685: VecAXPY(r, -1.0, u);
686: PetscObjectSetName((PetscObject) r, "Solution Error");
687: VecViewFromOptions(r, NULL, "-error_vec_view");
688: DMRestoreGlobalVector(dm, &r);
689: }
690: } else {
691: PetscReal res = 0.0;
693: /* Check discretization error */
694: PetscPrintf(PETSC_COMM_WORLD, "Initial guess\n");
695: VecView(u, PETSC_VIEWER_STDOUT_WORLD);
696: DMComputeL2Diff(dm, 0.0, user.exactFuncs, NULL, u, &error);
697: if (error >= 1.0e-11) {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: %g\n", error);}
698: else {PetscPrintf(PETSC_COMM_WORLD, "L_2 Error: < 1.0e-11\n");}
699: /* Check residual */
700: SNESComputeFunction(snes, u, r);
701: PetscPrintf(PETSC_COMM_WORLD, "Initial Residual\n");
702: VecChop(r, 1.0e-10);
703: VecView(r, PETSC_VIEWER_STDOUT_WORLD);
704: VecNorm(r, NORM_2, &res);
705: PetscPrintf(PETSC_COMM_WORLD, "L_2 Residual: %g\n", res);
706: /* Check Jacobian */
707: {
708: Mat J, M;
709: MatNullSpace nullspace;
710: Vec b;
711: PetscBool isNull;
713: SNESSetUpMatrices(snes);
714: SNESGetJacobian(snes, &J, &M, NULL, NULL);
715: SNESComputeJacobian(snes, u, J, M);
716: MatGetNullSpace(J, &nullspace);
717: MatNullSpaceTest(nullspace, J, &isNull);
718: VecDuplicate(u, &b);
719: VecSet(r, 0.0);
720: SNESComputeFunction(snes, r, b);
721: MatMult(J, u, r);
722: VecAXPY(r, 1.0, b);
723: VecDestroy(&b);
724: PetscPrintf(PETSC_COMM_WORLD, "Au - b = Au + F(0)\n");
725: VecChop(r, 1.0e-10);
726: VecView(r, PETSC_VIEWER_STDOUT_WORLD);
727: VecNorm(r, NORM_2, &res);
728: PetscPrintf(PETSC_COMM_WORLD, "Linear L_2 Residual: %g\n", res);
729: }
730: }
731: VecViewFromOptions(u, NULL, "-sol_vec_view");
733: VecDestroy(&u);
734: VecDestroy(&r);
735: SNESDestroy(&snes);
736: DMDestroy(&dm);
737: PetscFree(user.exactFuncs);
738: PetscFinalize();
739: return ierr;
740: }
742: /*TEST
744: # 2D serial P1 tests 0-3
745: test:
746: suffix: 0
747: requires: triangle
748: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
749: test:
750: suffix: 1
751: requires: triangle
752: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
753: test:
754: suffix: 2
755: requires: triangle
756: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
757: test:
758: suffix: 3
759: requires: triangle
760: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
761: # 2D serial P2 tests 4-5
762: test:
763: suffix: 4
764: requires: triangle
765: args: -run_type test -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
766: test:
767: suffix: 5
768: requires: triangle
769: args: -run_type test -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
770: # 2D serial P3 tests
771: test:
772: suffix: 2d_p3_0
773: requires: triangle
774: args: -run_type test -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 3 -pres_petscspace_degree 2
775: test:
776: suffix: 2d_p3_1
777: requires: triangle !single
778: args: -run_type full -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 3 -pres_petscspace_degree 2
779: # Parallel tests 6-17
780: test:
781: suffix: 6
782: requires: triangle
783: nsize: 2
784: args: -run_type test -refinement_limit 0.0 -test_partition -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
785: test:
786: suffix: 7
787: requires: triangle
788: nsize: 3
789: args: -run_type test -refinement_limit 0.0 -test_partition -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
790: test:
791: suffix: 8
792: requires: triangle
793: nsize: 5
794: args: -run_type test -refinement_limit 0.0 -test_partition -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
795: test:
796: suffix: 9
797: requires: triangle
798: nsize: 2
799: args: -run_type test -refinement_limit 0.0 -test_partition -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
800: test:
801: suffix: 10
802: requires: triangle
803: nsize: 3
804: args: -run_type test -refinement_limit 0.0 -test_partition -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
805: test:
806: suffix: 11
807: requires: triangle
808: nsize: 5
809: args: -run_type test -refinement_limit 0.0 -test_partition -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
810: test:
811: suffix: 12
812: requires: triangle
813: nsize: 2
814: args: -run_type test -refinement_limit 0.0625 -test_partition -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
815: test:
816: suffix: 13
817: requires: triangle
818: nsize: 3
819: args: -run_type test -refinement_limit 0.0625 -test_partition -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
820: test:
821: suffix: 14
822: requires: triangle
823: nsize: 5
824: args: -run_type test -refinement_limit 0.0625 -test_partition -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
825: test:
826: suffix: 15
827: requires: triangle
828: nsize: 2
829: args: -run_type test -refinement_limit 0.0625 -test_partition -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
830: test:
831: suffix: 16
832: requires: triangle
833: nsize: 3
834: args: -run_type test -refinement_limit 0.0625 -test_partition -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
835: test:
836: suffix: 17
837: requires: triangle
838: nsize: 5
839: args: -run_type test -refinement_limit 0.0625 -test_partition -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -dm_plex_print_fem 1
840: # 3D serial P1 tests 43-46
841: test:
842: suffix: 43
843: requires: ctetgen
844: args: -run_type test -dim 3 -refinement_limit 0.0 -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
845: test:
846: suffix: 44
847: requires: ctetgen
848: args: -run_type test -dim 3 -refinement_limit 0.0 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
849: test:
850: suffix: 45
851: requires: ctetgen
852: args: -run_type test -dim 3 -refinement_limit 0.0125 -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
853: test:
854: suffix: 46
855: requires: ctetgen
856: args: -run_type test -dim 3 -refinement_limit 0.0125 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -show_initial -dm_plex_print_fem 1
857: # Full solutions 18-29
858: test:
859: suffix: 18
860: requires: triangle !single
861: filter: sed -e "s/total number of linear solver iterations=11/total number of linear solver iterations=12/g"
862: args: -run_type full -refinement_limit 0.0625 -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
863: test:
864: suffix: 19
865: requires: triangle !single
866: nsize: 2
867: filter: sed -e "s/total number of linear solver iterations=11/total number of linear solver iterations=12/g"
868: args: -run_type full -petscpartitioner_type simple -refinement_limit 0.0625 -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
869: test:
870: suffix: 20
871: requires: triangle !single
872: filter: sed -e "s/total number of linear solver iterations=11/total number of linear solver iterations=12/g"
873: nsize: 3
874: args: -run_type full -petscpartitioner_type simple -refinement_limit 0.0625 -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
875: test:
876: suffix: 20_parmetis
877: requires: parmetis triangle !single
878: filter: sed -e "s/total number of linear solver iterations=11/total number of linear solver iterations=12/g"
879: nsize: 3
880: args: -run_type full -petscpartitioner_type parmetis -refinement_limit 0.0625 -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
881: test:
882: suffix: 21
883: requires: triangle !single
884: filter: sed -e "s/total number of linear solver iterations=11/total number of linear solver iterations=12/g"
885: nsize: 5
886: args: -run_type full -petscpartitioner_type simple -refinement_limit 0.0625 -bc_type dirichlet -interpolate 0 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
887: test:
888: suffix: 22
889: requires: triangle !single
890: filter: sed -e "s/total number of linear solver iterations=11/total number of linear solver iterations=12/g"
891: args: -run_type full -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
892: test:
893: suffix: 23
894: requires: triangle !single
895: nsize: 2
896: filter: sed -e "s/total number of linear solver iterations=11/total number of linear solver iterations=12/g"
897: args: -run_type full -petscpartitioner_type simple -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
898: test:
899: suffix: 24
900: requires: triangle !single
901: nsize: 3
902: filter: sed -e "s/total number of linear solver iterations=11/total number of linear solver iterations=12/g"
903: args: -run_type full -petscpartitioner_type simple -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
904: test:
905: suffix: 25
906: requires: triangle !single
907: filter: sed -e "s/total number of linear solver iterations=11/total number of linear solver iterations=12/g"
908: nsize: 5
909: args: -run_type full -petscpartitioner_type simple -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
910: test:
911: suffix: 26
912: requires: triangle !single
913: args: -run_type full -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
914: test:
915: suffix: 27
916: requires: triangle !single
917: nsize: 2
918: args: -run_type full -petscpartitioner_type simple -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
919: test:
920: suffix: 28
921: requires: triangle !single
922: nsize: 3
923: args: -run_type full -petscpartitioner_type simple -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
924: test:
925: suffix: 29
926: requires: triangle !single
927: nsize: 5
928: args: -run_type full -petscpartitioner_type simple -refinement_limit 0.0625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
929: # Full solutions with quads
930: # FULL Schur with LU/Jacobi
931: test:
932: suffix: quad_q2q1_full
933: requires: !single
934: args: -run_type full -simplex 0 -refinement_limit 0.00625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_velocity_ksp_type gmres -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
935: test:
936: suffix: quad_q2p1_full
937: requires: !single
938: args: -run_type full -simplex 0 -refinement_limit 0.00625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -pres_petscspace_poly_tensor 0 -pres_petscdualspace_lagrange_continuity 0 -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_velocity_ksp_type gmres -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
939: # Stokes preconditioners 30-36
940: # Jacobi
941: test:
942: suffix: 30
943: requires: triangle !single
944: filter: sed -e "s/total number of linear solver iterations=756/total number of linear solver iterations=757/g" -e "s/total number of linear solver iterations=758/total number of linear solver iterations=757/g"
945: args: -run_type full -refinement_limit 0.00625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -ksp_gmres_restart 100 -pc_type jacobi -ksp_rtol 1.0e-9 -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
946: # Block diagonal \begin{pmatrix} A & 0 \\ 0 & I \end{pmatrix}
947: test:
948: suffix: 31
949: requires: triangle !single
950: args: -run_type full -refinement_limit 0.00625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -ksp_type fgmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-4 -pc_type fieldsplit -pc_fieldsplit_type additive -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
951: # Block triangular \begin{pmatrix} A & B \\ 0 & I \end{pmatrix}
952: test:
953: suffix: 32
954: requires: triangle !single
955: args: -run_type full -refinement_limit 0.00625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -ksp_type fgmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type fieldsplit -pc_fieldsplit_type multiplicative -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
956: # Diagonal Schur complement \begin{pmatrix} A & 0 \\ 0 & S \end{pmatrix}
957: test:
958: suffix: 33
959: requires: triangle !single
960: args: -run_type full -refinement_limit 0.00625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -ksp_type fgmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type diag -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_velocity_ksp_type gmres -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
961: # Upper triangular Schur complement \begin{pmatrix} A & B \\ 0 & S \end{pmatrix}
962: test:
963: suffix: 34
964: requires: triangle !single
965: args: -run_type full -refinement_limit 0.00625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -ksp_type fgmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type upper -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_velocity_ksp_type gmres -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
966: # Lower triangular Schur complement \begin{pmatrix} A & B \\ 0 & S \end{pmatrix}
967: test:
968: suffix: 35
969: requires: triangle !single
970: args: -run_type full -refinement_limit 0.00625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -ksp_type fgmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type lower -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_velocity_ksp_type gmres -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
971: # Full Schur complement \begin{pmatrix} I & 0 \\ B^T A^{-1} & I \end{pmatrix} \begin{pmatrix} A & 0 \\ 0 & S \end{pmatrix} \begin{pmatrix} I & A^{-1} B \\ 0 & I \end{pmatrix}
972: test:
973: suffix: 36
974: requires: triangle !single
975: args: -run_type full -refinement_limit 0.00625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -ksp_type fgmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_velocity_ksp_type gmres -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
976: # SIMPLE \begin{pmatrix} I & 0 \\ B^T A^{-1} & I \end{pmatrix} \begin{pmatrix} A & 0 \\ 0 & B^T diag(A)^{-1} B \end{pmatrix} \begin{pmatrix} I & diag(A)^{-1} B \\ 0 & I \end{pmatrix}
977: test:
978: suffix: pc_simple
979: requires: triangle !single
980: args: -run_type full -refinement_limit 0.00625 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -ksp_type fgmres -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_velocity_ksp_type gmres -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi -fieldsplit_pressure_inner_ksp_type preonly -fieldsplit_pressure_inner_pc_type jacobi -fieldsplit_pressure_upper_ksp_type preonly -fieldsplit_pressure_upper_pc_type jacobi -snes_error_if_not_converged -ksp_error_if_not_converged -snes_view
981: # SIMPLEC \begin{pmatrix} I & 0 \\ B^T A^{-1} & I \end{pmatrix} \begin{pmatrix} A & 0 \\ 0 & B^T rowsum(A)^{-1} B \end{pmatrix} \begin{pmatrix} I & rowsum(A)^{-1} B \\ 0 & I \end{pmatrix}
982: test:
983: suffix: pc_simplec
984: requires: triangle
985: args: -run_type full -dm_refine 3 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -ksp_type fgmres -ksp_max_it 5 -ksp_gmres_restart 100 -ksp_rtol 1.0e-9 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_ksp_max_it 10 -fieldsplit_velocity_ksp_type gmres -fieldsplit_velocity_pc_type lu -fieldsplit_pressure_pc_type jacobi -fieldsplit_pressure_inner_ksp_type preonly -fieldsplit_pressure_inner_pc_type jacobi -fieldsplit_pressure_inner_pc_jacobi_type rowsum -fieldsplit_pressure_upper_ksp_type preonly -fieldsplit_pressure_upper_pc_type jacobi -fieldsplit_pressure_upper_pc_jacobi_type rowsum -snes_converged_reason -ksp_converged_reason -snes_view
986: # FETI-DP solvers
987: test:
988: suffix: fetidp_2d_tri
989: requires: triangle mumps
990: filter: grep -v "variant HERMITIAN"
991: nsize: 5
992: args: -run_type full -dm_refine 2 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -snes_view -snes_error_if_not_converged -dm_mat_type is -ksp_type fetidp -ksp_rtol 1.0e-8 -ksp_fetidp_saddlepoint -fetidp_ksp_type cg -fetidp_fieldsplit_p_ksp_max_it 1 -fetidp_fieldsplit_p_ksp_type richardson -fetidp_fieldsplit_p_ksp_richardson_scale 200 -fetidp_fieldsplit_p_pc_type none -ksp_fetidp_saddlepoint_flip 1 -fetidp_bddc_pc_bddc_vertex_size 2 -fetidp_bddc_pc_bddc_dirichlet_pc_factor_mat_solver_type mumps -fetidp_bddc_pc_bddc_neumann_pc_factor_mat_solver_type mumps -petscpartitioner_type simple -fetidp_fieldsplit_lag_ksp_type preonly
993: test:
994: suffix: fetidp_3d_tet
995: requires: ctetgen suitesparse
996: filter: grep -v "variant HERMITIAN" | sed -e "s/linear solver iterations=10[0-9]/linear solver iterations=100/g" | sed -e "s/linear solver iterations=9[0-9]/linear solver iterations=100/g"
997: nsize: 5
998: args: -run_type full -dm_refine 2 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -snes_view -snes_error_if_not_converged -dm_mat_type is -ksp_type fetidp -ksp_rtol 1.0e-8 -ksp_fetidp_saddlepoint -fetidp_ksp_type cg -fetidp_fieldsplit_p_ksp_max_it 1 -fetidp_fieldsplit_p_ksp_type richardson -fetidp_fieldsplit_p_ksp_richardson_scale 1000 -fetidp_fieldsplit_p_pc_type none -ksp_fetidp_saddlepoint_flip 1 -fetidp_bddc_pc_bddc_vertex_size 3 -fetidp_bddc_pc_bddc_use_deluxe_scaling -fetidp_bddc_pc_bddc_benign_trick -fetidp_bddc_pc_bddc_deluxe_singlemat -dim 3 -fetidp_pc_discrete_harmonic -fetidp_harmonic_pc_factor_mat_solver_type cholmod -fetidp_harmonic_pc_type cholesky -fetidp_bddelta_pc_factor_mat_solver_type umfpack -fetidp_fieldsplit_lag_ksp_type preonly -test_partition
1000: test:
1001: suffix: fetidp_2d_quad
1002: requires: mumps
1003: filter: grep -v "variant HERMITIAN"
1004: nsize: 5
1005: args: -run_type full -dm_refine 2 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -snes_view -snes_error_if_not_converged -dm_mat_type is -ksp_type fetidp -ksp_rtol 1.0e-8 -ksp_fetidp_saddlepoint -fetidp_ksp_type cg -fetidp_fieldsplit_p_ksp_max_it 1 -fetidp_fieldsplit_p_ksp_type richardson -fetidp_fieldsplit_p_ksp_richardson_scale 200 -fetidp_fieldsplit_p_pc_type none -ksp_fetidp_saddlepoint_flip 1 -fetidp_bddc_pc_bddc_vertex_size 2 -fetidp_bddc_pc_bddc_dirichlet_pc_factor_mat_solver_type mumps -fetidp_bddc_pc_bddc_neumann_pc_factor_mat_solver_type mumps -simplex 0 -petscpartitioner_type simple -fetidp_fieldsplit_lag_ksp_type preonly
1006: test:
1007: suffix: fetidp_3d_hex
1008: requires: suitesparse
1009: filter: grep -v "variant HERMITIAN" | sed -e "s/linear solver iterations=7[0-9]/linear solver iterations=71/g"
1010: nsize: 5
1011: args: -run_type full -dm_refine 1 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -snes_view -snes_error_if_not_converged -dm_mat_type is -ksp_type fetidp -ksp_rtol 1.0e-8 -ksp_fetidp_saddlepoint -fetidp_ksp_type cg -fetidp_fieldsplit_p_ksp_max_it 1 -fetidp_fieldsplit_p_ksp_type richardson -fetidp_fieldsplit_p_ksp_richardson_scale 2000 -fetidp_fieldsplit_p_pc_type none -ksp_fetidp_saddlepoint_flip 1 -fetidp_bddc_pc_bddc_vertex_size 3 -dim 3 -simplex 0 -fetidp_pc_discrete_harmonic -fetidp_harmonic_pc_factor_mat_solver_type cholmod -fetidp_harmonic_pc_type cholesky -petscpartitioner_type simple -fetidp_fieldsplit_lag_ksp_type preonly -fetidp_bddc_pc_bddc_dirichlet_pc_factor_mat_solver_type umfpack -fetidp_bddc_pc_bddc_neumann_pc_factor_mat_solver_type umfpack
1012: # Convergence
1013: test:
1014: suffix: 2d_quad_q1_p0_conv
1015: requires: !single
1016: args: -run_type full -bc_type dirichlet -simplex 0 -interpolate 1 -dm_refine 0 -vel_petscspace_degree 1 -pres_petscspace_degree 0 \
1017: -snes_convergence_estimate -convest_num_refine 3 -snes_error_if_not_converged \
1018: -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \
1019: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \
1020: -fieldsplit_velocity_pc_type lu \
1021: -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi
1022: test:
1023: suffix: 2d_tri_p2_p1_conv
1024: requires: triangle !single
1025: args: -run_type full -sol_type cubic -bc_type dirichlet -interpolate 1 -dm_refine 0 \
1026: -vel_petscspace_degree 2 -pres_petscspace_degree 1 \
1027: -snes_convergence_estimate -convest_num_refine 3 -snes_error_if_not_converged \
1028: -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \
1029: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \
1030: -fieldsplit_velocity_pc_type lu \
1031: -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi
1032: test:
1033: suffix: 2d_quad_q2_q1_conv
1034: requires: !single
1035: args: -run_type full -sol_type cubic -bc_type dirichlet -simplex 0 -interpolate 1 -dm_refine 0 \
1036: -vel_petscspace_degree 2 -pres_petscspace_degree 1 \
1037: -snes_convergence_estimate -convest_num_refine 3 -snes_error_if_not_converged \
1038: -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \
1039: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \
1040: -fieldsplit_velocity_pc_type lu \
1041: -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi
1042: test:
1043: suffix: 2d_quad_q2_p1_conv
1044: requires: !single
1045: args: -run_type full -sol_type cubic -bc_type dirichlet -simplex 0 -interpolate 1 -dm_refine 0 \
1046: -vel_petscspace_degree 2 -pres_petscspace_degree 1 -pres_petscspace_poly_tensor 0 -pres_petscdualspace_lagrange_continuity 0 \
1047: -snes_convergence_estimate -convest_num_refine 3 -snes_error_if_not_converged \
1048: -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \
1049: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \
1050: -fieldsplit_velocity_pc_type lu \
1051: -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi
1052: # Vanka solver
1053: test:
1054: suffix: 2d_quad_q1_p0_vanka_add
1055: requires: double !complex
1056: filter: sed -e "s/linear solver iterations=52/linear solver iterations=49/g" -e "s/Linear solve converged due to CONVERGED_RTOL iterations 52/Linear solve converged due to CONVERGED_RTOL iterations 49/g"
1057: args: -run_type full -bc_type dirichlet -simplex 0 -dm_refine 1 -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 0 -petscds_jac_pre 0 \
1058: -snes_rtol 1.0e-4 -snes_error_if_not_converged -snes_view -snes_monitor -snes_converged_reason \
1059: -ksp_type gmres -ksp_rtol 1.0e-5 -ksp_error_if_not_converged -ksp_converged_reason \
1060: -pc_type patch -pc_patch_partition_of_unity 0 -pc_patch_construct_codim 0 -pc_patch_construct_type vanka \
1061: -sub_ksp_type preonly -sub_pc_type lu
1062: test:
1063: suffix: 2d_quad_q1_p0_vanka_add_unity
1064: requires: double !complex
1065: filter: sed -e "s/linear solver iterations=46/linear solver iterations=45/g" -e "s/Linear solve converged due to CONVERGED_RTOL iterations 46/Linear solve converged due to CONVERGED_RTOL iterations 45/g"
1066: args: -run_type full -bc_type dirichlet -simplex 0 -dm_refine 1 -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 0 -petscds_jac_pre 0 \
1067: -snes_rtol 1.0e-4 -snes_error_if_not_converged -snes_view -snes_monitor -snes_converged_reason \
1068: -ksp_type gmres -ksp_rtol 1.0e-5 -ksp_error_if_not_converged -ksp_converged_reason \
1069: -pc_type patch -pc_patch_partition_of_unity 1 -pc_patch_construct_codim 0 -pc_patch_construct_type vanka \
1070: -sub_ksp_type preonly -sub_pc_type lu
1071: test:
1072: suffix: 2d_quad_q2_q1_vanka_add
1073: requires: double !complex
1074: filter: sed -e "s/linear solver iterations=[4-9][0-9][0-9]/linear solver iterations=489/g"
1075: args: -run_type full -bc_type dirichlet -simplex 0 -dm_refine 0 -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -petscds_jac_pre 0 \
1076: -snes_rtol 1.0e-4 -snes_error_if_not_converged -snes_view -snes_monitor -snes_converged_reason \
1077: -ksp_type gmres -ksp_rtol 1.0e-5 -ksp_error_if_not_converged \
1078: -pc_type patch -pc_patch_partition_of_unity 0 -pc_patch_construct_dim 0 -pc_patch_construct_type vanka \
1079: -sub_ksp_type preonly -sub_pc_type lu
1080: test:
1081: suffix: 2d_quad_q2_q1_vanka_add_unity
1082: requires: double !complex
1083: filter: sed -e "s/linear solver iterations=[4-9][0-9][0-9]/linear solver iterations=795/g"
1084: args: -run_type full -bc_type dirichlet -simplex 0 -dm_refine 0 -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -petscds_jac_pre 0 \
1085: -snes_rtol 1.0e-4 -snes_error_if_not_converged -snes_view -snes_monitor -snes_converged_reason \
1086: -ksp_type gmres -ksp_rtol 1.0e-5 -ksp_error_if_not_converged \
1087: -pc_type patch -pc_patch_partition_of_unity 1 -pc_patch_construct_dim 0 -pc_patch_construct_type vanka \
1088: -sub_ksp_type preonly -sub_pc_type lu
1089: # Vanka smoother
1090: test:
1091: suffix: 2d_quad_q1_p0_gmg_vanka_add
1092: requires: double !complex long_runtime
1093: args: -run_type full -bc_type dirichlet -simplex 0 -dm_refine_hierarchy 3 -interpolate 1 -vel_petscspace_degree 1 -pres_petscspace_degree 0 -petscds_jac_pre 0 \
1094: -snes_rtol 1.0e-4 -snes_error_if_not_converged -snes_view -snes_monitor -snes_converged_reason \
1095: -ksp_type gmres -ksp_rtol 1.0e-5 -ksp_error_if_not_converged -ksp_monitor_true_residual \
1096: -pc_type mg -pc_mg_levels 3 \
1097: -mg_levels_ksp_type gmres -mg_levels_ksp_max_it 30 -mg_levels_ksp_monitor_true_residual_no \
1098: -mg_levels_pc_type patch -mg_levels_pc_patch_partition_of_unity 0 -mg_levels_pc_patch_construct_codim 0 -mg_levels_pc_patch_construct_type vanka \
1099: -mg_levels_sub_ksp_type preonly -mg_levels_sub_pc_type lu \
1100: -mg_coarse_pc_type svd
1102: test:
1103: requires: !single
1104: suffix: bddc_quad
1105: nsize: 5
1106: args: -run_type full -dm_refine 2 -bc_type dirichlet -interpolate 1 -vel_petscspace_degree 2 -pres_petscspace_degree 1 -snes_view -snes_error_if_not_converged -dm_mat_type is -ksp_type gmres -ksp_rtol 1.e-8 -pc_type bddc -pc_bddc_corner_selection -pc_bddc_dirichlet_pc_type svd -pc_bddc_neumann_pc_type svd -pc_bddc_coarse_redundant_pc_type svd -simplex 0 -petscpartitioner_type simple -ksp_monitor_short -pc_bddc_symmetric 0
1108: TEST*/