Lines Matching +full:- +full:j
3 The 2D driven cavity problem is solved in a velocity-vorticity formulation.\n\
5 -lidvelocity <lid>, where <lid> = dimensionless velocity of lid\n\
6 -grashof <gr>, where <gr> = dimensionless temperature gradent\n\
7 -prandtl <pr>, where <pr> = dimensionless thermal/momentum diffusity ratio\n\
8 -contours : draw contour plots of solution\n\n";
10 See src/snes/tutorials/ex19.c for the steady-state version.
17 Continuation and Differential-Algebraic Equations, 2003.
27 nonlinear solve (e.g., -snes_type newtonls). The DAE versus ODE variants
28 are controlled using the -parabolic option.
37 a relative tolerance of 1.e-3. On the example problem, setting
38 -snes_type ksponly has only minor impact on number of steps, but
41 See also https://lists.mcs.anl.gov/pipermail/petsc-dev/2010-March/002362.html
44 /* ------------------------------------------------------------------------
51 - Lap(U) - Grad_y(Omega) = 0
52 - Lap(V) + Grad_x(Omega) = 0
53 Omega_t - Lap(Omega) + Div([U*Omega,V*Omega]) - GR*Grad_x(T) = 0
54 T_t - Lap(T) + PR*Div([U*T,V*T]) = 0
59 No-slip, rigid-wall Dirichlet conditions are used for [U,V].
61 vorticity: Omega = - Grad_y(U) + Grad_x(V), where along each
67 A finite difference approximation with the usual 5-point stencil
74 * applied matrix-free via the option -snes_mf
76 without a preconditioner due to ill-conditioning).
78 ------------------------------------------------------------------------- */
84 petscsys.h - base PETSc routines petscvec.h - vectors
85 petscmat.h - matrices
86 petscis.h - index sets petscksp.h - Krylov subspace methods
87 petscviewer.h - viewers petscpc.h - preconditioners
88 petscksp.h - linear solvers petscsnes.h - nonlinear solvers
95 User-defined routines and data structures
113 AppCtx user; /* user-defined work context */ in main()
140 …PetscCall(PetscOptionsReal("-lidvelocity", "Lid velocity, related to Reynolds number", "", user.li… in main()
141 …PetscCall(PetscOptionsReal("-prandtl", "Ratio of viscous to thermal diffusivity", "", user.prandtl… in main()
142 …PetscCall(PetscOptionsReal("-grashof", "Ratio of buoyant to viscous forces", "", user.grashof, &us… in main()
143 …PetscCall(PetscOptionsBool("-parabolic", "Relax incompressibility to make the system parabolic ins… in main()
144 …PetscCall(PetscOptionsReal("-cfl_initial", "Advective CFL for the first time step", "", user.cfl_i… in main()
147 PetscCall(DMDASetFieldName(da, 0, "x-velocity")); in main()
148 PetscCall(DMDASetFieldName(da, 1, "y-velocity")); in main()
152 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - in main()
155 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ in main()
157 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - in main()
159 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ in main()
170 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - in main()
172 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ in main()
184 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - in main()
187 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */ in main()
196 /* ------------------------------------------------------------------- */
199 FormInitialSolution - Forms initial approximation.
202 user - user-defined application context
203 X - vector
206 X - vector
211 PetscInt i, j, mx, xs, ys, xm, ym; in FormInitialSolution() local
216 grashof = user->grashof; in FormInitialSolution()
219 dx = 1.0 / (mx - 1); in FormInitialSolution()
222 Get local grid boundaries (for 2-dimensional DMDA): in FormInitialSolution()
223 xs, ys - starting grid indices (no ghost points) in FormInitialSolution()
224 xm, ym - widths of local grid (no ghost points) in FormInitialSolution()
230 - For default PETSc vectors, VecGetArray() returns a pointer to in FormInitialSolution()
232 - You MUST call VecRestoreArray() when you no longer need access to in FormInitialSolution()
241 for (j = ys; j < ys + ym; j++) { in FormInitialSolution()
243 x[j][i].u = 0.0; in FormInitialSolution()
244 x[j][i].v = 0.0; in FormInitialSolution()
245 x[j][i].omega = 0.0; in FormInitialSolution()
246 x[j][i].temp = (grashof > 0) * i * dx; in FormInitialSolution()
260 PetscInt xints, xinte, yints, yinte, i, j; in FormIFunctionLocal() local
266 grashof = user->grashof; in FormIFunctionLocal()
267 prandtl = user->prandtl; in FormIFunctionLocal()
268 lid = user->lidvelocity; in FormIFunctionLocal()
277 dhx = (PetscReal)(info->mx - 1); in FormIFunctionLocal()
278 dhy = (PetscReal)(info->my - 1); in FormIFunctionLocal()
284 xints = info->xs; in FormIFunctionLocal()
285 xinte = info->xs + info->xm; in FormIFunctionLocal()
286 yints = info->ys; in FormIFunctionLocal()
287 yinte = info->ys + info->ym; in FormIFunctionLocal()
291 j = 0; in FormIFunctionLocal()
294 for (i = info->xs; i < info->xs + info->xm; i++) { in FormIFunctionLocal()
295 f[j][i].u = x[j][i].u; in FormIFunctionLocal()
296 f[j][i].v = x[j][i].v; in FormIFunctionLocal()
297 f[j][i].omega = x[j][i].omega + (x[j + 1][i].u - x[j][i].u) * dhy; in FormIFunctionLocal()
298 f[j][i].temp = x[j][i].temp - x[j + 1][i].temp; in FormIFunctionLocal()
303 if (yinte == info->my) { in FormIFunctionLocal()
304 j = info->my - 1; in FormIFunctionLocal()
305 yinte = yinte - 1; in FormIFunctionLocal()
307 for (i = info->xs; i < info->xs + info->xm; i++) { in FormIFunctionLocal()
308 f[j][i].u = x[j][i].u - lid; in FormIFunctionLocal()
309 f[j][i].v = x[j][i].v; in FormIFunctionLocal()
310 f[j][i].omega = x[j][i].omega + (x[j][i].u - x[j - 1][i].u) * dhy; in FormIFunctionLocal()
311 f[j][i].temp = x[j][i].temp - x[j - 1][i].temp; in FormIFunctionLocal()
320 for (j = info->ys; j < info->ys + info->ym; j++) { in FormIFunctionLocal()
321 f[j][i].u = x[j][i].u; in FormIFunctionLocal()
322 f[j][i].v = x[j][i].v; in FormIFunctionLocal()
323 f[j][i].omega = x[j][i].omega - (x[j][i + 1].v - x[j][i].v) * dhx; in FormIFunctionLocal()
324 f[j][i].temp = x[j][i].temp; in FormIFunctionLocal()
329 if (xinte == info->mx) { in FormIFunctionLocal()
330 i = info->mx - 1; in FormIFunctionLocal()
331 xinte = xinte - 1; in FormIFunctionLocal()
333 for (j = info->ys; j < info->ys + info->ym; j++) { in FormIFunctionLocal()
334 f[j][i].u = x[j][i].u; in FormIFunctionLocal()
335 f[j][i].v = x[j][i].v; in FormIFunctionLocal()
336 f[j][i].omega = x[j][i].omega - (x[j][i].v - x[j][i - 1].v) * dhx; in FormIFunctionLocal()
337 f[j][i].temp = x[j][i].temp - (PetscReal)(grashof > 0); in FormIFunctionLocal()
342 for (j = yints; j < yinte; j++) { in FormIFunctionLocal()
347 vx = x[j][i].u; in FormIFunctionLocal()
350 vxm = .5 * (vx - avx); in FormIFunctionLocal()
351 vy = x[j][i].v; in FormIFunctionLocal()
354 vym = .5 * (vy - avy); in FormIFunctionLocal()
357 u = x[j][i].u; in FormIFunctionLocal()
358 udot = user->parabolic ? xdot[j][i].u : 0.; in FormIFunctionLocal()
359 uxx = (2.0 * u - x[j][i - 1].u - x[j][i + 1].u) * hydhx; in FormIFunctionLocal()
360 uyy = (2.0 * u - x[j - 1][i].u - x[j + 1][i].u) * hxdhy; in FormIFunctionLocal()
361 f[j][i].u = udot + uxx + uyy - .5 * (x[j + 1][i].omega - x[j - 1][i].omega) * hx; in FormIFunctionLocal()
364 u = x[j][i].v; in FormIFunctionLocal()
365 udot = user->parabolic ? xdot[j][i].v : 0.; in FormIFunctionLocal()
366 uxx = (2.0 * u - x[j][i - 1].v - x[j][i + 1].v) * hydhx; in FormIFunctionLocal()
367 uyy = (2.0 * u - x[j - 1][i].v - x[j + 1][i].v) * hxdhy; in FormIFunctionLocal()
368 f[j][i].v = udot + uxx + uyy + .5 * (x[j][i + 1].omega - x[j][i - 1].omega) * hy; in FormIFunctionLocal()
371 u = x[j][i].omega; in FormIFunctionLocal()
372 uxx = (2.0 * u - x[j][i - 1].omega - x[j][i + 1].omega) * hydhx; in FormIFunctionLocal()
373 uyy = (2.0 * u - x[j - 1][i].omega - x[j + 1][i].omega) * hxdhy; in FormIFunctionLocal()
374 …j][i].omega = (xdot[j][i].omega + uxx + uyy + (vxp * (u - x[j][i - 1].omega) + vxm * (x[j][i + 1].… in FormIFunctionLocal()
377 u = x[j][i].temp; in FormIFunctionLocal()
378 uxx = (2.0 * u - x[j][i - 1].temp - x[j][i + 1].temp) * hydhx; in FormIFunctionLocal()
379 uyy = (2.0 * u - x[j - 1][i].temp - x[j + 1][i].temp) * hxdhy; in FormIFunctionLocal()
380 …j][i].temp = (xdot[j][i].temp + uxx + uyy + prandtl * ((vxp * (u - x[j][i - 1].temp) + vxm * (x[j]… in FormIFunctionLocal()
385 Flop count (multiply-adds are counted as 2 operations) in FormIFunctionLocal()
387 PetscCall(PetscLogFlops(84.0 * info->ym * info->xm)); in FormIFunctionLocal()
394 …-da_grid_x 20 -da_grid_y 20 -lidvelocity 100 -grashof 1e3 -ts_max_steps 100 -ts_rtol 1e-3 -ts_atol…
399 …-da_grid_x 20 -da_grid_y 20 -lidvelocity 100 -grashof 1e3 -ts_run_steps 100 -ts_rtol 1e-3 -ts_atol…
405 …-da_grid_x 20 -da_grid_y 20 -lidvelocity 100 -grashof 1e3 -ts_max_steps 100 -ts_rtol 1e-3 -ts_atol…
411 …-da_refine 2 -lidvelocity 100 -grashof 1e3 -ts_max_steps 10 -ts_rtol 1e-3 -ts_atol 1e-3 -pc_type n…
417 args: -da_refine 1 -lidvelocity 100 -grashof 1e3 -ts_max_steps 10 -ts_rtol 1e-3 -ts_atol 1e-3
423 args: -da_refine 1 -lidvelocity 100 -grashof 1e3 -ts_max_steps 10 -ts_rtol 1e-3 -ts_atol 1e-3