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"; 4 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$. 15 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. 18 19 For visualization, use 20 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*/ 23 24 #include <petscdmplex.h> 25 #include <petscsnes.h> 26 #include <petscds.h> 27 #include <petscbag.h> 28 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; 35 36 typedef struct { 37 PetscBag bag; /* Holds problem parameters */ 38 } AppCtx; 39 40 /* 41 In 2D, plane Poiseuille flow has exact solution: 42 43 u = \Delta/(2 \nu) y (1 - y) + u_0 44 v = 0 45 p = -\Delta x 46 f = 0 47 48 so that 49 50 -\nu \Delta u + \nabla p + f = <\Delta, 0> + <-\Delta, 0> + <0, 0> = 0 51 \nabla \cdot u = 0 + 0 = 0 52 53 In 3D we use exact solution: 54 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 60 61 so that 62 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 65 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, void *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); 75 PetscInt d; 76 77 u[0] = u_0; 78 for (d = 1; d < dim; ++d) u[0] += Delta/(fac * 2.0*nu) * X[d] * (1.0 - X[d]); 79 for (d = 1; d < dim; ++d) u[d] = 0.0; 80 return 0; 81 } 82 83 PetscErrorCode linear_p(PetscInt dim, PetscReal time, const PetscReal X[], PetscInt Nf, PetscScalar *p, void *ctx) 84 { 85 Parameter *param = (Parameter *) ctx; 86 PetscReal Delta = param->Delta; 87 88 p[0] = -Delta * X[0]; 89 return 0; 90 } 91 92 PetscErrorCode wall_velocity(PetscInt dim, PetscReal time, const PetscReal X[], PetscInt Nf, PetscScalar *u, void *ctx) 93 { 94 Parameter *param = (Parameter *) ctx; 95 PetscReal u_0 = param->u_0; 96 PetscInt d; 97 98 u[0] = u_0; 99 for (d = 1; d < dim; ++d) u[d] = 0.0; 100 return 0; 101 } 102 103 /* 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} 104 u[Ncomp] = {p} */ 105 void f1_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, 106 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 107 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 108 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]) 109 { 110 const PetscReal nu = PetscRealPart(constants[1]); 111 const PetscInt Nc = dim; 112 PetscInt c, d; 113 114 for (c = 0; c < Nc; ++c) { 115 for (d = 0; d < dim; ++d) { 116 /* f1[c*dim+d] = 0.5*nu*(u_x[c*dim+d] + u_x[d*dim+c]); */ 117 f1[c*dim+d] = nu*u_x[c*dim+d]; 118 } 119 f1[c*dim+c] -= u[uOff[1]]; 120 } 121 } 122 123 /* 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} */ 124 void f0_p(PetscInt dim, PetscInt Nf, PetscInt NfAux, 125 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 126 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 127 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]) 128 { 129 PetscInt d; 130 for (d = 0, f0[0] = 0.0; d < dim; ++d) f0[0] += u_x[d*dim+d]; 131 } 132 133 /* Residual functions are in reference coordinates */ 134 static void f0_bd_u(PetscInt dim, PetscInt Nf, PetscInt NfAux, 135 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 136 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 137 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]) 138 { 139 const PetscReal Delta = PetscRealPart(constants[0]); 140 PetscReal alpha = PetscRealPart(constants[3]); 141 PetscReal X = PetscCosReal(alpha)*x[0] + PetscSinReal(alpha)*x[1]; 142 PetscInt d; 143 144 for (d = 0; d < dim; ++d) { 145 f0[d] = -Delta * X * n[d]; 146 } 147 } 148 149 /* < q, \nabla\cdot u > 150 NcompI = 1, NcompJ = dim */ 151 void g1_pu(PetscInt dim, PetscInt Nf, PetscInt NfAux, 152 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 153 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 154 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]) 155 { 156 PetscInt d; 157 for (d = 0; d < dim; ++d) g1[d*dim+d] = 1.0; /* \frac{\partial\phi^{u_d}}{\partial x_d} */ 158 } 159 160 /* -< \nabla\cdot v, p > 161 NcompI = dim, NcompJ = 1 */ 162 void g2_up(PetscInt dim, PetscInt Nf, PetscInt NfAux, 163 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 164 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 165 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]) 166 { 167 PetscInt d; 168 for (d = 0; d < dim; ++d) g2[d*dim+d] = -1.0; /* \frac{\partial\psi^{u_d}}{\partial x_d} */ 169 } 170 171 /* < \nabla v, \nabla u + {\nabla u}^T > 172 This just gives \nabla u, give the perdiagonal for the transpose */ 173 void g3_uu(PetscInt dim, PetscInt Nf, PetscInt NfAux, 174 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], 175 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], 176 PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]) 177 { 178 const PetscReal nu = PetscRealPart(constants[1]); 179 const PetscInt Nc = dim; 180 PetscInt c, d; 181 182 for (c = 0; c < Nc; ++c) { 183 for (d = 0; d < dim; ++d) { 184 g3[((c*Nc+c)*dim+d)*dim+d] = nu; 185 } 186 } 187 } 188 189 static PetscErrorCode SetupParameters(AppCtx *user) 190 { 191 PetscBag bag; 192 Parameter *p; 193 PetscErrorCode ierr; 194 195 PetscFunctionBeginUser; 196 /* setup PETSc parameter bag */ 197 ierr = PetscBagGetData(user->bag, (void **) &p);CHKERRQ(ierr); 198 ierr = PetscBagSetName(user->bag, "par", "Poiseuille flow parameters");CHKERRQ(ierr); 199 bag = user->bag; 200 ierr = PetscBagRegisterReal(bag, &p->Delta, 1.0, "Delta", "Pressure drop per unit length");CHKERRQ(ierr); 201 ierr = PetscBagRegisterReal(bag, &p->nu, 1.0, "nu", "Kinematic viscosity");CHKERRQ(ierr); 202 ierr = PetscBagRegisterReal(bag, &p->u_0, 0.0, "u_0", "Tangential velocity at the wall");CHKERRQ(ierr); 203 ierr = PetscBagRegisterReal(bag, &p->alpha, 0.0, "alpha", "Angle of pipe wall to x-axis");CHKERRQ(ierr); 204 PetscFunctionReturn(0); 205 } 206 207 PetscErrorCode CreateMesh(MPI_Comm comm, AppCtx *user, DM *dm) 208 { 209 PetscErrorCode ierr; 210 211 PetscFunctionBeginUser; 212 ierr = DMCreate(comm, dm);CHKERRQ(ierr); 213 ierr = DMSetType(*dm, DMPLEX);CHKERRQ(ierr); 214 ierr = DMSetFromOptions(*dm);CHKERRQ(ierr); 215 { 216 Parameter *param; 217 Vec coordinates; 218 PetscScalar *coords; 219 PetscReal alpha; 220 PetscInt cdim, N, bs, i; 221 222 ierr = DMGetCoordinateDim(*dm, &cdim);CHKERRQ(ierr); 223 ierr = DMGetCoordinates(*dm, &coordinates);CHKERRQ(ierr); 224 ierr = VecGetLocalSize(coordinates, &N);CHKERRQ(ierr); 225 ierr = VecGetBlockSize(coordinates, &bs);CHKERRQ(ierr); 226 if (bs != cdim) SETERRQ2(comm, PETSC_ERR_ARG_WRONG, "Invalid coordinate blocksize %D != embedding dimension %D", bs, cdim); 227 ierr = VecGetArray(coordinates, &coords);CHKERRQ(ierr); 228 ierr = PetscBagGetData(user->bag, (void **) ¶m);CHKERRQ(ierr); 229 alpha = param->alpha; 230 for (i = 0; i < N; i += cdim) { 231 PetscScalar x = coords[i+0]; 232 PetscScalar y = coords[i+1]; 233 234 coords[i+0] = PetscCosReal(alpha)*x - PetscSinReal(alpha)*y; 235 coords[i+1] = PetscSinReal(alpha)*x + PetscCosReal(alpha)*y; 236 } 237 ierr = VecRestoreArray(coordinates, &coords);CHKERRQ(ierr); 238 ierr = DMSetCoordinates(*dm, coordinates);CHKERRQ(ierr); 239 } 240 ierr = DMViewFromOptions(*dm, NULL, "-dm_view");CHKERRQ(ierr); 241 PetscFunctionReturn(0); 242 } 243 244 PetscErrorCode SetupProblem(DM dm, AppCtx *user) 245 { 246 PetscDS ds; 247 PetscWeakForm wf; 248 DMLabel label; 249 Parameter *ctx; 250 PetscInt id, bd; 251 PetscErrorCode ierr; 252 253 PetscFunctionBeginUser; 254 ierr = PetscBagGetData(user->bag, (void **) &ctx);CHKERRQ(ierr); 255 ierr = DMGetDS(dm, &ds);CHKERRQ(ierr); 256 ierr = PetscDSSetResidual(ds, 0, NULL, f1_u);CHKERRQ(ierr); 257 ierr = PetscDSSetResidual(ds, 1, f0_p, NULL);CHKERRQ(ierr); 258 ierr = PetscDSSetJacobian(ds, 0, 0, NULL, NULL, NULL, g3_uu);CHKERRQ(ierr); 259 ierr = PetscDSSetJacobian(ds, 0, 1, NULL, NULL, g2_up, NULL);CHKERRQ(ierr); 260 ierr = PetscDSSetJacobian(ds, 1, 0, NULL, g1_pu, NULL, NULL);CHKERRQ(ierr); 261 262 id = 2; 263 ierr = DMGetLabel(dm, "marker", &label);CHKERRQ(ierr); 264 ierr = DMAddBoundary(dm, DM_BC_NATURAL, "right wall", label, 1, &id, 0, 0, NULL, NULL, NULL, ctx, &bd);CHKERRQ(ierr); 265 ierr = PetscDSGetBoundary(ds, bd, &wf, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL);CHKERRQ(ierr); 266 ierr = PetscWeakFormSetIndexBdResidual(wf, label, id, 0, 0, f0_bd_u, 0, NULL);CHKERRQ(ierr); 267 /* Setup constants */ 268 { 269 Parameter *param; 270 PetscScalar constants[4]; 271 272 ierr = PetscBagGetData(user->bag, (void **) ¶m);CHKERRQ(ierr); 273 274 constants[0] = param->Delta; 275 constants[1] = param->nu; 276 constants[2] = param->u_0; 277 constants[3] = param->alpha; 278 ierr = PetscDSSetConstants(ds, 4, constants);CHKERRQ(ierr); 279 } 280 /* Setup Boundary Conditions */ 281 id = 3; 282 ierr = DMAddBoundary(dm, DM_BC_ESSENTIAL, "top wall", label, 1, &id, 0, 0, NULL, (void (*)(void)) wall_velocity, NULL, ctx, NULL);CHKERRQ(ierr); 283 id = 1; 284 ierr = DMAddBoundary(dm, DM_BC_ESSENTIAL, "bottom wall", label, 1, &id, 0, 0, NULL, (void (*)(void)) wall_velocity, NULL, ctx, NULL);CHKERRQ(ierr); 285 /* Setup exact solution */ 286 ierr = PetscDSSetExactSolution(ds, 0, quadratic_u, ctx);CHKERRQ(ierr); 287 ierr = PetscDSSetExactSolution(ds, 1, linear_p, ctx);CHKERRQ(ierr); 288 PetscFunctionReturn(0); 289 } 290 291 PetscErrorCode SetupDiscretization(DM dm, AppCtx *user) 292 { 293 DM cdm = dm; 294 PetscFE fe[2]; 295 Parameter *param; 296 PetscBool simplex; 297 PetscInt dim; 298 MPI_Comm comm; 299 PetscErrorCode ierr; 300 301 PetscFunctionBeginUser; 302 ierr = DMGetDimension(dm, &dim);CHKERRQ(ierr); 303 ierr = DMPlexIsSimplex(dm, &simplex);CHKERRQ(ierr); 304 ierr = PetscObjectGetComm((PetscObject) dm, &comm);CHKERRQ(ierr); 305 ierr = PetscFECreateDefault(comm, dim, dim, simplex, "vel_", PETSC_DEFAULT, &fe[0]);CHKERRQ(ierr); 306 ierr = PetscObjectSetName((PetscObject) fe[0], "velocity");CHKERRQ(ierr); 307 ierr = PetscFECreateDefault(comm, dim, 1, simplex, "pres_", PETSC_DEFAULT, &fe[1]);CHKERRQ(ierr); 308 ierr = PetscFECopyQuadrature(fe[0], fe[1]);CHKERRQ(ierr); 309 ierr = PetscObjectSetName((PetscObject) fe[1], "pressure");CHKERRQ(ierr); 310 /* Set discretization and boundary conditions for each mesh */ 311 ierr = DMSetField(dm, 0, NULL, (PetscObject) fe[0]);CHKERRQ(ierr); 312 ierr = DMSetField(dm, 1, NULL, (PetscObject) fe[1]);CHKERRQ(ierr); 313 ierr = DMCreateDS(dm);CHKERRQ(ierr); 314 ierr = SetupProblem(dm, user);CHKERRQ(ierr); 315 ierr = PetscBagGetData(user->bag, (void **) ¶m);CHKERRQ(ierr); 316 while (cdm) { 317 ierr = DMCopyDisc(dm, cdm);CHKERRQ(ierr); 318 ierr = DMPlexCreateBasisRotation(cdm, param->alpha, 0.0, 0.0);CHKERRQ(ierr); 319 ierr = DMGetCoarseDM(cdm, &cdm);CHKERRQ(ierr); 320 } 321 ierr = PetscFEDestroy(&fe[0]);CHKERRQ(ierr); 322 ierr = PetscFEDestroy(&fe[1]);CHKERRQ(ierr); 323 PetscFunctionReturn(0); 324 } 325 326 int main(int argc, char **argv) 327 { 328 SNES snes; /* nonlinear solver */ 329 DM dm; /* problem definition */ 330 Vec u, r; /* solution and residual */ 331 AppCtx user; /* user-defined work context */ 332 PetscErrorCode ierr; 333 334 ierr = PetscInitialize(&argc, &argv, NULL,help);if (ierr) return ierr; 335 ierr = PetscBagCreate(PETSC_COMM_WORLD, sizeof(Parameter), &user.bag);CHKERRQ(ierr); 336 ierr = SetupParameters(&user);CHKERRQ(ierr); 337 ierr = PetscBagSetFromOptions(user.bag);CHKERRQ(ierr); 338 ierr = SNESCreate(PETSC_COMM_WORLD, &snes);CHKERRQ(ierr); 339 ierr = CreateMesh(PETSC_COMM_WORLD, &user, &dm);CHKERRQ(ierr); 340 ierr = SNESSetDM(snes, dm);CHKERRQ(ierr); 341 ierr = DMSetApplicationContext(dm, &user);CHKERRQ(ierr); 342 /* Setup problem */ 343 ierr = SetupDiscretization(dm, &user);CHKERRQ(ierr); 344 ierr = DMPlexCreateClosureIndex(dm, NULL);CHKERRQ(ierr); 345 346 ierr = DMCreateGlobalVector(dm, &u);CHKERRQ(ierr); 347 ierr = VecDuplicate(u, &r);CHKERRQ(ierr); 348 349 ierr = DMPlexSetSNESLocalFEM(dm,&user,&user,&user);CHKERRQ(ierr); 350 351 ierr = SNESSetFromOptions(snes);CHKERRQ(ierr); 352 353 { 354 PetscDS ds; 355 PetscSimplePointFunc exactFuncs[2]; 356 void *ctxs[2]; 357 358 ierr = DMGetDS(dm, &ds);CHKERRQ(ierr); 359 ierr = PetscDSGetExactSolution(ds, 0, &exactFuncs[0], &ctxs[0]);CHKERRQ(ierr); 360 ierr = PetscDSGetExactSolution(ds, 1, &exactFuncs[1], &ctxs[1]);CHKERRQ(ierr); 361 ierr = DMProjectFunction(dm, 0.0, exactFuncs, ctxs, INSERT_ALL_VALUES, u);CHKERRQ(ierr); 362 ierr = PetscObjectSetName((PetscObject) u, "Exact Solution");CHKERRQ(ierr); 363 ierr = VecViewFromOptions(u, NULL, "-exact_vec_view");CHKERRQ(ierr); 364 } 365 ierr = DMSNESCheckFromOptions(snes, u);CHKERRQ(ierr); 366 ierr = VecSet(u, 0.0);CHKERRQ(ierr); 367 ierr = PetscObjectSetName((PetscObject) u, "Solution");CHKERRQ(ierr); 368 ierr = SNESSolve(snes, NULL, u);CHKERRQ(ierr); 369 ierr = VecViewFromOptions(u, NULL, "-sol_vec_view");CHKERRQ(ierr); 370 371 ierr = VecDestroy(&u);CHKERRQ(ierr); 372 ierr = VecDestroy(&r);CHKERRQ(ierr); 373 ierr = DMDestroy(&dm);CHKERRQ(ierr); 374 ierr = SNESDestroy(&snes);CHKERRQ(ierr); 375 ierr = PetscBagDestroy(&user.bag);CHKERRQ(ierr); 376 ierr = PetscFinalize(); 377 return ierr; 378 } 379 380 /*TEST 381 382 # Convergence 383 test: 384 suffix: 2d_quad_q1_p0_conv 385 requires: !single 386 args: -dm_plex_simplex 0 -dm_plex_separate_marker -dm_refine 1 \ 387 -vel_petscspace_degree 1 -pres_petscspace_degree 0 \ 388 -snes_convergence_estimate -convest_num_refine 2 -snes_error_if_not_converged \ 389 -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \ 390 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \ 391 -fieldsplit_velocity_pc_type lu \ 392 -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi 393 test: 394 suffix: 2d_quad_q1_p0_conv_u0 395 requires: !single 396 args: -dm_plex_simplex 0 -dm_plex_separate_marker -dm_refine 1 -u_0 0.125 \ 397 -vel_petscspace_degree 1 -pres_petscspace_degree 0 \ 398 -snes_convergence_estimate -convest_num_refine 2 -snes_error_if_not_converged \ 399 -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \ 400 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \ 401 -fieldsplit_velocity_pc_type lu \ 402 -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi 403 test: 404 suffix: 2d_quad_q1_p0_conv_u0_alpha 405 requires: !single 406 args: -dm_plex_simplex 0 -dm_plex_separate_marker -dm_refine 1 -u_0 0.125 -alpha 0.3927 \ 407 -vel_petscspace_degree 1 -pres_petscspace_degree 0 \ 408 -snes_convergence_estimate -convest_num_refine 2 -snes_error_if_not_converged \ 409 -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \ 410 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \ 411 -fieldsplit_velocity_pc_type lu \ 412 -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi 413 test: 414 suffix: 2d_quad_q1_p0_conv_gmg_vanka 415 requires: !single long_runtime 416 args: -dm_plex_simplex 0 -dm_plex_separate_marker -dm_plex_box_faces 2,2 -dm_refine_hierarchy 1 \ 417 -vel_petscspace_degree 1 -pres_petscspace_degree 0 \ 418 -snes_convergence_estimate -convest_num_refine 1 -snes_error_if_not_converged \ 419 -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \ 420 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \ 421 -fieldsplit_velocity_pc_type mg \ 422 -fieldsplit_velocity_mg_levels_pc_type patch -fieldsplit_velocity_mg_levels_pc_patch_exclude_subspaces 1 \ 423 -fieldsplit_velocity_mg_levels_pc_patch_construct_codim 0 -fieldsplit_velocity_mg_levels_pc_patch_construct_type vanka \ 424 -fieldsplit_pressure_ksp_rtol 1e-5 -fieldsplit_pressure_pc_type jacobi 425 test: 426 suffix: 2d_tri_p2_p1_conv 427 requires: triangle !single 428 args: -dm_plex_separate_marker -dm_refine 1 \ 429 -vel_petscspace_degree 2 -pres_petscspace_degree 1 \ 430 -dmsnes_check .001 -snes_error_if_not_converged \ 431 -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \ 432 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \ 433 -fieldsplit_velocity_pc_type lu \ 434 -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi 435 test: 436 suffix: 2d_tri_p2_p1_conv_u0_alpha 437 requires: triangle !single 438 args: -dm_plex_separate_marker -dm_refine 0 -u_0 0.125 -alpha 0.3927 \ 439 -vel_petscspace_degree 2 -pres_petscspace_degree 1 \ 440 -dmsnes_check .001 -snes_error_if_not_converged \ 441 -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \ 442 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \ 443 -fieldsplit_velocity_pc_type lu \ 444 -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi 445 test: 446 suffix: 2d_tri_p2_p1_conv_gmg_vcycle 447 requires: triangle !single 448 args: -dm_plex_separate_marker -dm_plex_box_faces 2,2 -dm_refine_hierarchy 1 \ 449 -vel_petscspace_degree 2 -pres_petscspace_degree 1 \ 450 -dmsnes_check .001 -snes_error_if_not_converged \ 451 -ksp_type fgmres -ksp_gmres_restart 10 -ksp_rtol 1.0e-9 -ksp_error_if_not_converged \ 452 -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_factorization_type full \ 453 -fieldsplit_velocity_pc_type mg \ 454 -fieldsplit_pressure_ksp_rtol 1e-10 -fieldsplit_pressure_pc_type jacobi 455 TEST*/ 456