Actual source code: ex71.c
1: static char help[] = "Poiseuille Flow in 2d and 3d channels with finite elements.\n\
2: We solve the Poiseuille flow problem in a rectangular\n\
3: domain, using a parallel unstructured mesh (DMPLEX) to discretize it.\n\n\n";
5: /*F
6: A Poiseuille flow is a steady-state isoviscous Stokes flow in a pipe of constant cross-section. We discretize using the
7: finite element method on an unstructured mesh. The weak form equations are
8: \begin{align*}
9: < \nabla v, \nu (\nabla u + {\nabla u}^T) > - < \nabla\cdot v, p > + < v, \Delta \hat n >_{\Gamma_o} = 0
10: < q, \nabla\cdot u > = 0
11: \end{align*}
12: where $\nu$ is the kinematic viscosity, $\Delta$ is the pressure drop per unit length, assuming that pressure is 0 on
13: the left edge, and $\Gamma_o$ is the outlet boundary at the right edge of the pipe. The normal velocity will be zero at
14: the wall, but we will allow a fixed tangential velocity $u_0$.
16: In order to test our global to local basis transformation, we will allow the pipe to be at an angle $\alpha$ to the
17: coordinate axes.
19: For visualization, use
21: -dm_view hdf5:$PWD/sol.h5 -sol_vec_view hdf5:$PWD/sol.h5::append -exact_vec_view hdf5:$PWD/sol.h5::append
22: F*/
24: #include <petscdmplex.h>
25: #include <petscsnes.h>
26: #include <petscds.h>
27: #include <petscbag.h>
29: typedef struct {
30: PetscReal Delta; /* Pressure drop per unit length */
31: PetscReal nu; /* Kinematic viscosity */
32: PetscReal u_0; /* Tangential velocity at the wall */
33: PetscReal alpha; /* Angle of pipe wall to x-axis */
34: } Parameter;
36: typedef struct {
37: PetscBag bag; /* Holds problem parameters */
38: } AppCtx;
40: /*
41: In 2D, plane Poiseuille flow has exact solution:
43: u = \Delta/(2 \nu) y (1 - y) + u_0
44: v = 0
45: p = -\Delta x
46: f = 0
48: so that
50: -\nu \Delta u + \nabla p + f = <\Delta, 0> + <-\Delta, 0> + <0, 0> = 0
51: \nabla \cdot u = 0 + 0 = 0
53: In 3D we use exact solution:
55: u = \Delta/(4 \nu) (y (1 - y) + z (1 - z)) + u_0
56: v = 0
57: w = 0
58: p = -\Delta x
59: f = 0
61: so that
63: -\nu \Delta u + \nabla p + f = <Delta, 0, 0> + <-Delta, 0, 0> + <0, 0, 0> = 0
64: \nabla \cdot u = 0 + 0 + 0 = 0
66: Note that these functions use coordinates X in the global (rotated) frame
67: */
68: PetscErrorCode quadratic_u(PetscInt dim, PetscReal time, const PetscReal X[], PetscInt Nf, PetscScalar *u, PetscCtx ctx)
69: {
70: Parameter *param = (Parameter *)ctx;
71: PetscReal Delta = param->Delta;
72: PetscReal nu = param->nu;
73: PetscReal u_0 = param->u_0;
74: PetscReal fac = (PetscReal)(dim - 1);
76: u[0] = u_0;
77: for (PetscInt d = 1; d < dim; ++d) u[0] += Delta / (fac * 2.0 * nu) * X[d] * (1.0 - X[d]);
78: for (PetscInt d = 1; d < dim; ++d) u[d] = 0.0;
79: return PETSC_SUCCESS;
80: }
82: PetscErrorCode linear_p(PetscInt dim, PetscReal time, const PetscReal X[], PetscInt Nf, PetscScalar *p, PetscCtx ctx)
83: {
84: Parameter *param = (Parameter *)ctx;
85: PetscReal Delta = param->Delta;
87: p[0] = -Delta * X[0];
88: return PETSC_SUCCESS;
89: }
91: PetscErrorCode wall_velocity(PetscInt dim, PetscReal time, const PetscReal X[], PetscInt Nf, PetscScalar *u, PetscCtx ctx)
92: {
93: Parameter *param = (Parameter *)ctx;
94: PetscReal u_0 = param->u_0;
96: u[0] = u_0;
97: for (PetscInt d = 1; d < dim; ++d) u[d] = 0.0;
98: return PETSC_SUCCESS;
99: }
101: /* 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}
102: u[Ncomp] = {p} */
103: 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[])
104: {
105: const PetscReal nu = PetscRealPart(constants[1]);
106: const PetscInt Nc = dim;
108: for (PetscInt c = 0; c < Nc; ++c) {
109: for (PetscInt d = 0; d < dim; ++d) {
110: /* f1[c*dim+d] = 0.5*nu*(u_x[c*dim+d] + u_x[d*dim+c]); */
111: f1[c * dim + d] = nu * u_x[c * dim + d];
112: }
113: f1[c * dim + c] -= u[uOff[1]];
114: }
115: }
117: /* 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} */
118: void f0_p(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[])
119: {
120: f0[0] = 0.0;
121: for (PetscInt d = 0; d < dim; ++d) f0[0] += u_x[d * dim + d];
122: }
124: /* Residual functions are in reference coordinates */
125: 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[])
126: {
127: const PetscReal Delta = PetscRealPart(constants[0]);
128: PetscReal alpha = PetscRealPart(constants[3]);
129: PetscReal X = PetscCosReal(alpha) * x[0] + PetscSinReal(alpha) * x[1];
131: for (PetscInt d = 0; d < dim; ++d) f0[d] = -Delta * X * n[d];
132: }
134: /* < q, \nabla\cdot u >
135: NcompI = 1, NcompJ = dim */
136: void g1_pu(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 g1[])
137: {
138: for (PetscInt d = 0; d < dim; ++d) g1[d * dim + d] = 1.0; /* \frac{\partial\phi^{u_d}}{\partial x_d} */
139: }
141: /* -< \nabla\cdot v, p >
142: NcompI = dim, NcompJ = 1 */
143: void g2_up(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 g2[])
144: {
145: for (PetscInt d = 0; d < dim; ++d) g2[d * dim + d] = -1.0; /* \frac{\partial\psi^{u_d}}{\partial x_d} */
146: }
148: /* < \nabla v, \nabla u + {\nabla u}^T >
149: This just gives \nabla u, give the perdiagonal for the transpose */
150: 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[])
151: {
152: const PetscReal nu = PetscRealPart(constants[1]);
153: const PetscInt Nc = dim;
155: for (PetscInt c = 0; c < Nc; ++c) {
156: for (PetscInt d = 0; d < dim; ++d) g3[((c * Nc + c) * dim + d) * dim + d] = nu;
157: }
158: }
160: static PetscErrorCode SetupParameters(AppCtx *user)
161: {
162: PetscBag bag;
163: Parameter *p;
165: PetscFunctionBeginUser;
166: /* setup PETSc parameter bag */
167: PetscCall(PetscBagGetData(user->bag, &p));
168: PetscCall(PetscBagSetName(user->bag, "par", "Poiseuille flow parameters"));
169: bag = user->bag;
170: PetscCall(PetscBagRegisterReal(bag, &p->Delta, 1.0, "Delta", "Pressure drop per unit length"));
171: PetscCall(PetscBagRegisterReal(bag, &p->nu, 1.0, "nu", "Kinematic viscosity"));
172: PetscCall(PetscBagRegisterReal(bag, &p->u_0, 0.0, "u_0", "Tangential velocity at the wall"));
173: PetscCall(PetscBagRegisterReal(bag, &p->alpha, 0.0, "alpha", "Angle of pipe wall to x-axis"));
174: PetscFunctionReturn(PETSC_SUCCESS);
175: }
177: PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm)
178: {
179: PetscFunctionBeginUser;
180: PetscCall(DMCreate(comm, dm));
181: PetscCall(DMSetType(*dm, DMPLEX));
182: PetscCall(DMSetFromOptions(*dm));
183: {
184: Parameter *param;
185: Vec coordinates;
186: PetscScalar *coords;
187: PetscReal alpha;
188: PetscInt cdim, N, bs, i;
190: PetscCall(DMGetCoordinateDim(*dm, &cdim));
191: PetscCall(DMGetCoordinates(*dm, &coordinates));
192: PetscCall(VecGetLocalSize(coordinates, &N));
193: PetscCall(VecGetBlockSize(coordinates, &bs));
194: PetscCheck(bs == cdim, comm, PETSC_ERR_ARG_WRONG, "Invalid coordinate blocksize %" PetscInt_FMT " != embedding dimension %" PetscInt_FMT, bs, cdim);
195: PetscCall(VecGetArray(coordinates, &coords));
196: PetscCall(PetscBagGetData(user->bag, ¶m));
197: alpha = param->alpha;
198: for (i = 0; i < N; i += cdim) {
199: PetscScalar x = coords[i + 0];
200: PetscScalar y = coords[i + 1];
202: coords[i + 0] = PetscCosReal(alpha) * x - PetscSinReal(alpha) * y;
203: coords[i + 1] = PetscSinReal(alpha) * x + PetscCosReal(alpha) * y;
204: }
205: PetscCall(VecRestoreArray(coordinates, &coords));
206: PetscCall(DMSetCoordinates(*dm, coordinates));
207: }
208: PetscCall(DMViewFromOptions(*dm, NULL, "-dm_view"));
209: PetscFunctionReturn(PETSC_SUCCESS);
210: }
212: PetscErrorCode SetupProblem(DM dm, AppCtx *user)
213: {
214: PetscDS ds;
215: PetscWeakForm wf;
216: DMLabel label;
217: Parameter *ctx;
218: PetscInt id, bd;
220: PetscFunctionBeginUser;
221: PetscCall(PetscBagGetData(user->bag, &ctx));
222: PetscCall(DMGetDS(dm, &ds));
223: PetscCall(PetscDSSetResidual(ds, 0, NULL, f1_u));
224: PetscCall(PetscDSSetResidual(ds, 1, f0_p, NULL));
225: PetscCall(PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu));
226: PetscCall(PetscDSSetJacobian(ds, 0, 1, NULL, NULL, g2_up, NULL));
227: PetscCall(PetscDSSetJacobian(ds, 1, 0, NULL, g1_pu, NULL, NULL));
229: id = 2;
230: PetscCall(DMGetLabel(dm, "marker", &label));
231: PetscCall(DMAddBoundary(dm, DM_BC_NATURAL, "right wall", label, 1, &id, 0, 0, NULL, NULL, NULL, ctx, &bd));
232: PetscCall(PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL));
233: PetscCall(PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, 0, f0_bd_u, 0, NULL));
234: /* Setup constants */
235: {
236: Parameter *param;
237: PetscScalar constants[4];
239: PetscCall(PetscBagGetData(user->bag, ¶m));
241: constants[0] = param->Delta;
242: constants[1] = param->nu;
243: constants[2] = param->u_0;
244: constants[3] = param->alpha;
245: PetscCall(PetscDSSetConstants(ds, 4, constants));
246: }
247: /* Setup Boundary Conditions */
248: id = 3;
249: PetscCall(DMAddBoundary(dm, DM_BC_ESSENTIAL, "top wall", label, 1, &id, 0, 0, NULL, (PetscVoidFn *)wall_velocity, NULL, ctx, NULL));
250: id = 1;
251: PetscCall(DMAddBoundary(dm, DM_BC_ESSENTIAL, "bottom wall", label, 1, &id, 0, 0, NULL, (PetscVoidFn *)wall_velocity, NULL, ctx, NULL));
252: /* Setup exact solution */
253: PetscCall(PetscDSSetExactSolution(ds, 0, quadratic_u, ctx));
254: PetscCall(PetscDSSetExactSolution(ds, 1, linear_p, ctx));
255: PetscFunctionReturn(PETSC_SUCCESS);
256: }
258: PetscErrorCode SetupDiscretization(DM dm, AppCtx *user)
259: {
260: DM cdm = dm;
261: PetscFE fe[2];
262: Parameter *param;
263: PetscBool simplex;
264: PetscInt dim;
265: MPI_Comm comm;
267: PetscFunctionBeginUser;
268: PetscCall(DMGetDimension(dm, &dim));
269: PetscCall(DMPlexIsSimplex(dm, &simplex));
270: PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
271: PetscCall(PetscFECreateDefault(comm, dim, dim, simplex, "vel_", PETSC_DEFAULT, &fe[0]));
272: PetscCall(PetscObjectSetName((PetscObject)fe[0], "velocity"));
273: PetscCall(PetscFECreateDefault(comm, dim, 1, simplex, "pres_", PETSC_DEFAULT, &fe[1]));
274: PetscCall(PetscFECopyQuadrature(fe[0], fe[1]));
275: PetscCall(PetscObjectSetName((PetscObject)fe[1], "pressure"));
276: /* Set discretization and boundary conditions for each mesh */
277: PetscCall(DMSetField(dm, 0, NULL, (PetscObject)fe[0]));
278: PetscCall(DMSetField(dm, 1, NULL, (PetscObject)fe[1]));
279: PetscCall(DMCreateDS(dm));
280: PetscCall(SetupProblem(dm, user));
281: PetscCall(PetscBagGetData(user->bag, ¶m));
282: while (cdm) {
283: PetscCall(DMCopyDisc(dm, cdm));
284: PetscCall(DMPlexCreateBasisRotation(cdm, param->alpha, 0.0, 0.0));
285: PetscCall(DMGetCoarseDM(cdm, &cdm));
286: }
287: PetscCall(PetscFEDestroy(&fe[0]));
288: PetscCall(PetscFEDestroy(&fe[1]));
289: PetscFunctionReturn(PETSC_SUCCESS);
290: }
292: int main(int argc, char **argv)
293: {
294: SNES snes; /* nonlinear solver */
295: DM dm; /* problem definition */
296: Vec u, r; /* solution and residual */
297: AppCtx user; /* user-defined work context */
299: PetscFunctionBeginUser;
300: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
301: PetscCall(PetscBagCreate(PETSC_COMM_WORLD, sizeof(Parameter), &user.bag));
302: PetscCall(SetupParameters(&user));
303: PetscCall(PetscBagSetFromOptions(user.bag));
304: PetscCall(SNESCreate(PETSC_COMM_WORLD, &snes));
305: PetscCall(CreateMesh(PETSC_COMM_WORLD, &user, &dm));
306: PetscCall(SNESSetDM(snes, dm));
307: PetscCall(DMSetApplicationContext(dm, &user));
308: /* Setup problem */
309: PetscCall(SetupDiscretization(dm, &user));
310: PetscCall(DMPlexCreateClosureIndex(dm, NULL));
312: PetscCall(DMCreateGlobalVector(dm, &u));
313: PetscCall(VecDuplicate(u, &r));
315: PetscCall(DMPlexSetSNESLocalFEM(dm, PETSC_FALSE, &user));
317: PetscCall(SNESSetFromOptions(snes));
319: {
320: PetscDS ds;
321: PetscSimplePointFn *exactFuncs[2];
322: void *ctxs[2];
324: PetscCall(DMGetDS(dm, &ds));
325: PetscCall(PetscDSGetExactSolution(ds, 0, &exactFuncs[0], &ctxs[0]));
326: PetscCall(PetscDSGetExactSolution(ds, 1, &exactFuncs[1], &ctxs[1]));
327: PetscCall(DMProjectFunction(dm, 0.0, exactFuncs, ctxs, INSERT_ALL_VALUES, u));
328: PetscCall(PetscObjectSetName((PetscObject)u, "Exact Solution"));
329: PetscCall(VecViewFromOptions(u, NULL, "-exact_vec_view"));
330: }
331: PetscCall(DMSNESCheckFromOptions(snes, u));
332: PetscCall(VecSet(u, 0.0));
333: PetscCall(PetscObjectSetName((PetscObject)u, "Solution"));
334: PetscCall(SNESSolve(snes, NULL, u));
335: PetscCall(VecViewFromOptions(u, NULL, "-sol_vec_view"));
337: PetscCall(VecDestroy(&u));
338: PetscCall(VecDestroy(&r));
339: PetscCall(DMDestroy(&dm));
340: PetscCall(SNESDestroy(&snes));
341: PetscCall(PetscBagDestroy(&user.bag));
342: PetscCall(PetscFinalize());
343: return 0;
344: }
346: /*TEST
348: # Convergence
349: test:
350: suffix: 2d_quad_q1_p0_conv
351: requires: !single
352: args: -dm_plex_simplex 0 -dm_plex_separate_marker -dm_refine 1 \
353: -vel_petscspace_degree 1 -pres_petscspace_degree 0 \
354: -snes_convergence_estimate -convest_num_refine 2 -snes_error_if_not_converged \
355: -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \
356: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \
357: -fieldsplit_velocity_pc_type lu \
358: -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi
359: test:
360: suffix: 2d_quad_q1_p0_conv_u0
361: requires: !single
362: args: -dm_plex_simplex 0 -dm_plex_separate_marker -dm_refine 1 -u_0 0.125 \
363: -vel_petscspace_degree 1 -pres_petscspace_degree 0 \
364: -snes_convergence_estimate -convest_num_refine 2 -snes_error_if_not_converged \
365: -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \
366: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \
367: -fieldsplit_velocity_pc_type lu \
368: -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi
369: test:
370: suffix: 2d_quad_q1_p0_conv_u0_alpha
371: requires: !single
372: args: -dm_plex_simplex 0 -dm_plex_separate_marker -dm_refine 1 -u_0 0.125 -alpha 0.3927 \
373: -vel_petscspace_degree 1 -pres_petscspace_degree 0 \
374: -snes_convergence_estimate -convest_num_refine 2 -snes_error_if_not_converged \
375: -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \
376: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \
377: -fieldsplit_velocity_pc_type lu \
378: -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi
379: test:
380: suffix: 2d_quad_q1_p0_conv_gmg_vanka
381: requires: !single long_runtime
382: args: -dm_plex_simplex 0 -dm_plex_separate_marker -dm_plex_box_faces 2,2 -dm_refine_hierarchy 1 \
383: -vel_petscspace_degree 1 -pres_petscspace_degree 0 \
384: -snes_convergence_estimate -convest_num_refine 1 -snes_error_if_not_converged \
385: -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \
386: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \
387: -fieldsplit_velocity_pc_type mg \
388: -fieldsplit_velocity_mg_levels_pc_type patch -fieldsplit_velocity_mg_levels_pc_patch_exclude_subspaces 1 \
389: -fieldsplit_velocity_mg_levels_pc_patch_construct_codim 0 -fieldsplit_velocity_mg_levels_pc_patch_construct_type vanka \
390: -fieldsplit_pressure_ksp_rtol 1e-5 -fieldsplit_pressure_pc_type jacobi
391: test:
392: suffix: 2d_tri_p2_p1_conv
393: requires: triangle !single
394: args: -dm_plex_separate_marker -dm_refine 1 \
395: -vel_petscspace_degree 2 -pres_petscspace_degree 1 \
396: -dmsnes_check .001 -snes_error_if_not_converged \
397: -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \
398: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \
399: -fieldsplit_velocity_pc_type lu \
400: -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi
401: test:
402: suffix: 2d_tri_p2_p1_conv_u0_alpha
403: requires: triangle !single
404: args: -dm_plex_separate_marker -dm_refine 0 -u_0 0.125 -alpha 0.3927 \
405: -vel_petscspace_degree 2 -pres_petscspace_degree 1 \
406: -dmsnes_check .001 -snes_error_if_not_converged \
407: -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \
408: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \
409: -fieldsplit_velocity_pc_type lu \
410: -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi
411: test:
412: suffix: 2d_tri_p2_p1_conv_gmg_vcycle
413: TODO: broken (requires subDMs hooks)
414: requires: triangle !single
415: args: -dm_plex_separate_marker -dm_plex_box_faces 2,2 -dm_refine_hierarchy 1 \
416: -vel_petscspace_degree 2 -pres_petscspace_degree 1 \
417: -dmsnes_check .001 -snes_error_if_not_converged \
418: -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \
419: -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \
420: -fieldsplit_velocity_pc_type mg \
421: -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi
422: TEST*/