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