Lines Matching +full:- +full:j
3 …mall stepsize and stored into a binary file, e.g. -ts_type ssp -ts_time_step 1e-5 -ts_max_steps 30…
4 Errors can be computed in the following runs with -simulation -f reference.bin
7 -ts_rk_dtratio is the ratio between larger time step size and small time step size
8 -ts_rk_multirate_type has three choices: none (for single step RK)
13 static const char help[] = "1D periodic Finite Volume solver in slope-limiter form with semidiscret…
14 " advection - Variable coefficient scalar advection\n"
20 " you should type -simulation -f file.bin.\n"
21 " you can choose the number of grids by -da_grid_x.\n"
22 … " you can choose the value of a by -physics_advect_a1 and -physics_advect_a2.\n";
36 PetscReal range = xmax - xmin; in RangeMod()
40 /* --------------------------------- Advection ----------------------------------- */
52 speed = ctx->a; in PhysicsRiemann_Advect()
68 if (x < 0) speeds[0] = ctx->a[0]; /* x is discontinuous point of a */ in PhysicsCharacteristic_Advect()
69 else speeds[0] = ctx->a[1]; in PhysicsCharacteristic_Advect()
76 PetscReal *a = ctx->a, x0; in PhysicsSample_Advect()
82 x0 = x - a[0] * t; in PhysicsSample_Advect()
85 x0 = RangeMod(x - a[0] * t, xmin, xmax); in PhysicsSample_Advect()
92 u[0] = (x0 < 0) ? 1 : -1; in PhysicsSample_Advect()
95 u[0] = (x0 < 0) ? -1 : 1; in PhysicsSample_Advect()
110 u[0] = (x0 < 0) ? 0 : ((x0 < 1) ? x0 : ((x0 < 2) ? 2 - x0 : 0)); in PhysicsSample_Advect()
121 x0 = x - a[1] * t; in PhysicsSample_Advect()
124 x0 = RangeMod(x - a[1] * t, xmin, xmax); in PhysicsSample_Advect()
131 u[0] = (x0 < 0) ? 1 : -1; in PhysicsSample_Advect()
134 u[0] = (x0 < 0) ? -1 : 1; in PhysicsSample_Advect()
149 u[0] = (x0 < 0) ? 0 : ((x0 < 1) ? x0 : ((x0 < 2) ? 2 - x0 : 0)); in PhysicsSample_Advect()
167 ctx->physics.sample = PhysicsSample_Advect; in PhysicsCreate_Advect()
168 ctx->physics.riemann = PhysicsRiemann_Advect; in PhysicsCreate_Advect()
169 ctx->physics.characteristic = PhysicsCharacteristic_Advect; in PhysicsCreate_Advect()
170 ctx->physics.destroy = PhysicsDestroy_SimpleFree; in PhysicsCreate_Advect()
171 ctx->physics.user = user; in PhysicsCreate_Advect()
172 ctx->physics.dof = 1; in PhysicsCreate_Advect()
173 PetscCall(PetscStrallocpy("u", &ctx->physics.fieldname[0])); in PhysicsCreate_Advect()
174 user->a[0] = 1; in PhysicsCreate_Advect()
175 user->a[1] = 1; in PhysicsCreate_Advect()
176 PetscOptionsBegin(ctx->comm, ctx->prefix, "Options for advection", ""); in PhysicsCreate_Advect()
178 PetscCall(PetscOptionsReal("-physics_advect_a1", "Speed1", "", user->a[0], &user->a[0], NULL)); in PhysicsCreate_Advect()
179 PetscCall(PetscOptionsReal("-physics_advect_a2", "Speed2", "", user->a[1], &user->a[1], NULL)); in PhysicsCreate_Advect()
185 /* --------------------------------- Finite Volume Solver for slow parts --------------------------…
190 PetscInt i, j, k, Mx, dof, xs, xm, len_slow; in FVRHSFunctionslow() local
200 hx = (ctx->xmax - ctx->xmin) / Mx; in FVRHSFunctionslow()
208 PetscCall(ISGetSize(ctx->iss, &len_slow)); in FVRHSFunctionslow()
210 if (ctx->bctype == FVBC_OUTFLOW) { in FVRHSFunctionslow()
211 for (i = xs - 2; i < 0; i++) { in FVRHSFunctionslow()
212 for (j = 0; j < dof; j++) x[i * dof + j] = x[j]; in FVRHSFunctionslow()
215 for (j = 0; j < dof; j++) x[i * dof + j] = x[(xs + xm - 1) * dof + j]; in FVRHSFunctionslow()
218 for (i = xs - 1; i < xs + xm + 1; i++) { in FVRHSFunctionslow()
222 /* Determine the right eigenvectors R, where A = R \Lambda R^{-1} */ in FVRHSFunctionslow()
223 …PetscCall((*ctx->physics.characteristic)(ctx->physics.user, dof, &x[i * dof], ctx->R, ctx->Rinv, c… in FVRHSFunctionslow()
224 /* Evaluate jumps across interfaces (i-1, i) and (i, i+1), put in characteristic basis */ in FVRHSFunctionslow()
225 PetscCall(PetscArrayzero(ctx->cjmpLR, 2 * dof)); in FVRHSFunctionslow()
226 cjmpL = &ctx->cjmpLR[0]; in FVRHSFunctionslow()
227 cjmpR = &ctx->cjmpLR[dof]; in FVRHSFunctionslow()
228 for (j = 0; j < dof; j++) { in FVRHSFunctionslow()
230 jmpL = x[(i + 0) * dof + j] - x[(i - 1) * dof + j]; in FVRHSFunctionslow()
231 jmpR = x[(i + 1) * dof + j] - x[(i + 0) * dof + j]; in FVRHSFunctionslow()
233 cjmpL[k] += ctx->Rinv[k + j * dof] * jmpL; in FVRHSFunctionslow()
234 cjmpR[k] += ctx->Rinv[k + j * dof] * jmpR; in FVRHSFunctionslow()
240 (*ctx->limit)(&info, cjmpL, cjmpR, ctx->cslope); in FVRHSFunctionslow()
241 for (j = 0; j < dof; j++) ctx->cslope[j] /= hx; /* rescale to a slope */ in FVRHSFunctionslow()
242 for (j = 0; j < dof; j++) { in FVRHSFunctionslow()
244 for (k = 0; k < dof; k++) tmp += ctx->R[j + k * dof] * ctx->cslope[k]; in FVRHSFunctionslow()
245 slope[i * dof + j] = tmp; in FVRHSFunctionslow()
254 uL = &ctx->uLR[0]; in FVRHSFunctionslow()
255 uR = &ctx->uLR[dof]; in FVRHSFunctionslow()
256 for (j = 0; j < dof; j++) { in FVRHSFunctionslow()
257 uL[j] = x[(i - 1) * dof + j] + slope[(i - 1) * dof + j] * hx / 2; in FVRHSFunctionslow()
258 uR[j] = x[(i - 0) * dof + j] - slope[(i - 0) * dof + j] * hx / 2; in FVRHSFunctionslow()
260 …PetscCall((*ctx->physics.riemann)(ctx->physics.user, dof, uL, uR, ctx->flux, &maxspeed, ctx->xmin … in FVRHSFunctionslow()
263 for (j = 0; j < dof; j++) f[(i - 1) * dof + j] -= ctx->flux[j] / hx; in FVRHSFunctionslow()
266 for (j = 0; j < dof; j++) f[i * dof + j] += ctx->flux[j] / hx; in FVRHSFunctionslow()
270 uL = &ctx->uLR[0]; in FVRHSFunctionslow()
271 uR = &ctx->uLR[dof]; in FVRHSFunctionslow()
272 for (j = 0; j < dof; j++) { in FVRHSFunctionslow()
273 uL[j] = x[(i - 1) * dof + j] + slope[(i - 1) * dof + j] * hx / 2; in FVRHSFunctionslow()
274 uR[j] = x[(i - 0) * dof + j] - slope[(i - 0) * dof + j] * hx / 2; in FVRHSFunctionslow()
276 …PetscCall((*ctx->physics.riemann)(ctx->physics.user, dof, uL, uR, ctx->flux, &maxspeed, ctx->xmin … in FVRHSFunctionslow()
279 for (j = 0; j < dof; j++) f[(i - 1) * dof + j] -= ctx->flux[j] / hx; in FVRHSFunctionslow()
290 /* --------------------------------- Finite Volume Solver for fast parts -------------------------…
295 PetscInt i, j, k, Mx, dof, xs, xm, len_slow; in FVRHSFunctionfast() local
305 hx = (ctx->xmax - ctx->xmin) / Mx; in FVRHSFunctionfast()
313 PetscCall(ISGetSize(ctx->iss, &len_slow)); in FVRHSFunctionfast()
315 if (ctx->bctype == FVBC_OUTFLOW) { in FVRHSFunctionfast()
316 for (i = xs - 2; i < 0; i++) { in FVRHSFunctionfast()
317 for (j = 0; j < dof; j++) x[i * dof + j] = x[j]; in FVRHSFunctionfast()
320 for (j = 0; j < dof; j++) x[i * dof + j] = x[(xs + xm - 1) * dof + j]; in FVRHSFunctionfast()
323 for (i = xs - 1; i < xs + xm + 1; i++) { in FVRHSFunctionfast()
326 if (i > len_slow - 2) { in FVRHSFunctionfast()
327 …PetscCall((*ctx->physics.characteristic)(ctx->physics.user, dof, &x[i * dof], ctx->R, ctx->Rinv, c… in FVRHSFunctionfast()
328 PetscCall(PetscArrayzero(ctx->cjmpLR, 2 * dof)); in FVRHSFunctionfast()
329 cjmpL = &ctx->cjmpLR[0]; in FVRHSFunctionfast()
330 cjmpR = &ctx->cjmpLR[dof]; in FVRHSFunctionfast()
331 for (j = 0; j < dof; j++) { in FVRHSFunctionfast()
333 jmpL = x[(i + 0) * dof + j] - x[(i - 1) * dof + j]; in FVRHSFunctionfast()
334 jmpR = x[(i + 1) * dof + j] - x[(i + 0) * dof + j]; in FVRHSFunctionfast()
336 cjmpL[k] += ctx->Rinv[k + j * dof] * jmpL; in FVRHSFunctionfast()
337 cjmpR[k] += ctx->Rinv[k + j * dof] * jmpR; in FVRHSFunctionfast()
343 (*ctx->limit)(&info, cjmpL, cjmpR, ctx->cslope); in FVRHSFunctionfast()
344 for (j = 0; j < dof; j++) ctx->cslope[j] /= hx; /* rescale to a slope */ in FVRHSFunctionfast()
345 for (j = 0; j < dof; j++) { in FVRHSFunctionfast()
347 for (k = 0; k < dof; k++) tmp += ctx->R[j + k * dof] * ctx->cslope[k]; in FVRHSFunctionfast()
348 slope[i * dof + j] = tmp; in FVRHSFunctionfast()
357 uL = &ctx->uLR[0]; in FVRHSFunctionfast()
358 uR = &ctx->uLR[dof]; in FVRHSFunctionfast()
359 for (j = 0; j < dof; j++) { in FVRHSFunctionfast()
360 uL[j] = x[(i - 1) * dof + j] + slope[(i - 1) * dof + j] * hx / 2; in FVRHSFunctionfast()
361 uR[j] = x[(i - 0) * dof + j] - slope[(i - 0) * dof + j] * hx / 2; in FVRHSFunctionfast()
363 …PetscCall((*ctx->physics.riemann)(ctx->physics.user, dof, uL, uR, ctx->flux, &maxspeed, ctx->xmin … in FVRHSFunctionfast()
366 for (j = 0; j < dof; j++) f[(i - len_slow - 1) * dof + j] -= ctx->flux[j] / hx; in FVRHSFunctionfast()
369 for (j = 0; j < dof; j++) f[(i - len_slow) * dof + j] += ctx->flux[j] / hx; in FVRHSFunctionfast()
373 uL = &ctx->uLR[0]; in FVRHSFunctionfast()
374 uR = &ctx->uLR[dof]; in FVRHSFunctionfast()
375 for (j = 0; j < dof; j++) { in FVRHSFunctionfast()
376 uL[j] = x[(i - 1) * dof + j] + slope[(i - 1) * dof + j] * hx / 2; in FVRHSFunctionfast()
377 uR[j] = x[(i - 0) * dof + j] - slope[(i - 0) * dof + j] * hx / 2; in FVRHSFunctionfast()
379 …PetscCall((*ctx->physics.riemann)(ctx->physics.user, dof, uL, uR, ctx->flux, &maxspeed, ctx->xmin … in FVRHSFunctionfast()
382 for (j = 0; j < dof; j++) f[(i - len_slow) * dof + j] += ctx->flux[j] / hx; in FVRHSFunctionfast()
396 PetscInt i, j, k, Mx, dof, xs, xm, len_slow1, len_slow2; in FVRHSFunctionslow2() local
406 hx = (ctx->xmax - ctx->xmin) / Mx; in FVRHSFunctionslow2()
414 PetscCall(ISGetSize(ctx->iss, &len_slow1)); in FVRHSFunctionslow2()
415 PetscCall(ISGetSize(ctx->iss2, &len_slow2)); in FVRHSFunctionslow2()
416 if (ctx->bctype == FVBC_OUTFLOW) { in FVRHSFunctionslow2()
417 for (i = xs - 2; i < 0; i++) { in FVRHSFunctionslow2()
418 for (j = 0; j < dof; j++) x[i * dof + j] = x[j]; in FVRHSFunctionslow2()
421 for (j = 0; j < dof; j++) x[i * dof + j] = x[(xs + xm - 1) * dof + j]; in FVRHSFunctionslow2()
424 for (i = xs - 1; i < xs + xm + 1; i++) { in FVRHSFunctionslow2()
427 if (i < len_slow1 + len_slow2 + 1 && i > len_slow1 - 2) { in FVRHSFunctionslow2()
428 …PetscCall((*ctx->physics.characteristic)(ctx->physics.user, dof, &x[i * dof], ctx->R, ctx->Rinv, c… in FVRHSFunctionslow2()
429 PetscCall(PetscArrayzero(ctx->cjmpLR, 2 * dof)); in FVRHSFunctionslow2()
430 cjmpL = &ctx->cjmpLR[0]; in FVRHSFunctionslow2()
431 cjmpR = &ctx->cjmpLR[dof]; in FVRHSFunctionslow2()
432 for (j = 0; j < dof; j++) { in FVRHSFunctionslow2()
434 jmpL = x[(i + 0) * dof + j] - x[(i - 1) * dof + j]; in FVRHSFunctionslow2()
435 jmpR = x[(i + 1) * dof + j] - x[(i + 0) * dof + j]; in FVRHSFunctionslow2()
437 cjmpL[k] += ctx->Rinv[k + j * dof] * jmpL; in FVRHSFunctionslow2()
438 cjmpR[k] += ctx->Rinv[k + j * dof] * jmpR; in FVRHSFunctionslow2()
444 (*ctx->limit)(&info, cjmpL, cjmpR, ctx->cslope); in FVRHSFunctionslow2()
445 for (j = 0; j < dof; j++) ctx->cslope[j] /= hx; /* rescale to a slope */ in FVRHSFunctionslow2()
446 for (j = 0; j < dof; j++) { in FVRHSFunctionslow2()
448 for (k = 0; k < dof; k++) tmp += ctx->R[j + k * dof] * ctx->cslope[k]; in FVRHSFunctionslow2()
449 slope[i * dof + j] = tmp; in FVRHSFunctionslow2()
458 uL = &ctx->uLR[0]; in FVRHSFunctionslow2()
459 uR = &ctx->uLR[dof]; in FVRHSFunctionslow2()
460 for (j = 0; j < dof; j++) { in FVRHSFunctionslow2()
461 uL[j] = x[(i - 1) * dof + j] + slope[(i - 1) * dof + j] * hx / 2; in FVRHSFunctionslow2()
462 uR[j] = x[(i - 0) * dof + j] - slope[(i - 0) * dof + j] * hx / 2; in FVRHSFunctionslow2()
464 …PetscCall((*ctx->physics.riemann)(ctx->physics.user, dof, uL, uR, ctx->flux, &maxspeed, ctx->xmin … in FVRHSFunctionslow2()
467 for (j = 0; j < dof; j++) f[(i - len_slow1 - 1) * dof + j] -= ctx->flux[j] / hx; in FVRHSFunctionslow2()
470 for (j = 0; j < dof; j++) f[(i - len_slow1) * dof + j] += ctx->flux[j] / hx; in FVRHSFunctionslow2()
474 uL = &ctx->uLR[0]; in FVRHSFunctionslow2()
475 uR = &ctx->uLR[dof]; in FVRHSFunctionslow2()
476 for (j = 0; j < dof; j++) { in FVRHSFunctionslow2()
477 uL[j] = x[(i - 1) * dof + j] + slope[(i - 1) * dof + j] * hx / 2; in FVRHSFunctionslow2()
478 uR[j] = x[(i - 0) * dof + j] - slope[(i - 0) * dof + j] * hx / 2; in FVRHSFunctionslow2()
480 …PetscCall((*ctx->physics.riemann)(ctx->physics.user, dof, uL, uR, ctx->flux, &maxspeed, ctx->xmin … in FVRHSFunctionslow2()
483 for (j = 0; j < dof; j++) f[(i - len_slow1 - 1) * dof + j] -= ctx->flux[j] / hx; in FVRHSFunctionslow2()
487 uL = &ctx->uLR[0]; in FVRHSFunctionslow2()
488 uR = &ctx->uLR[dof]; in FVRHSFunctionslow2()
489 for (j = 0; j < dof; j++) { in FVRHSFunctionslow2()
490 uL[j] = x[(i - 1) * dof + j] + slope[(i - 1) * dof + j] * hx / 2; in FVRHSFunctionslow2()
491 uR[j] = x[(i - 0) * dof + j] - slope[(i - 0) * dof + j] * hx / 2; in FVRHSFunctionslow2()
493 …PetscCall((*ctx->physics.riemann)(ctx->physics.user, dof, uL, uR, ctx->flux, &maxspeed, ctx->xmin … in FVRHSFunctionslow2()
496 for (j = 0; j < dof; j++) f[(i - len_slow1) * dof + j] += ctx->flux[j] / hx; in FVRHSFunctionslow2()
511 PetscInt i, j, k, Mx, dof, xs, xm, len_slow1, len_slow2; in FVRHSFunctionfast2() local
521 hx = (ctx->xmax - ctx->xmin) / Mx; in FVRHSFunctionfast2()
529 PetscCall(ISGetSize(ctx->iss, &len_slow1)); in FVRHSFunctionfast2()
530 PetscCall(ISGetSize(ctx->iss2, &len_slow2)); in FVRHSFunctionfast2()
532 if (ctx->bctype == FVBC_OUTFLOW) { in FVRHSFunctionfast2()
533 for (i = xs - 2; i < 0; i++) { in FVRHSFunctionfast2()
534 for (j = 0; j < dof; j++) x[i * dof + j] = x[j]; in FVRHSFunctionfast2()
537 for (j = 0; j < dof; j++) x[i * dof + j] = x[(xs + xm - 1) * dof + j]; in FVRHSFunctionfast2()
540 for (i = xs - 1; i < xs + xm + 1; i++) { in FVRHSFunctionfast2()
543 if (i > len_slow1 + len_slow2 - 2) { in FVRHSFunctionfast2()
544 …PetscCall((*ctx->physics.characteristic)(ctx->physics.user, dof, &x[i * dof], ctx->R, ctx->Rinv, c… in FVRHSFunctionfast2()
545 PetscCall(PetscArrayzero(ctx->cjmpLR, 2 * dof)); in FVRHSFunctionfast2()
546 cjmpL = &ctx->cjmpLR[0]; in FVRHSFunctionfast2()
547 cjmpR = &ctx->cjmpLR[dof]; in FVRHSFunctionfast2()
548 for (j = 0; j < dof; j++) { in FVRHSFunctionfast2()
550 jmpL = x[(i + 0) * dof + j] - x[(i - 1) * dof + j]; in FVRHSFunctionfast2()
551 jmpR = x[(i + 1) * dof + j] - x[(i + 0) * dof + j]; in FVRHSFunctionfast2()
553 cjmpL[k] += ctx->Rinv[k + j * dof] * jmpL; in FVRHSFunctionfast2()
554 cjmpR[k] += ctx->Rinv[k + j * dof] * jmpR; in FVRHSFunctionfast2()
560 (*ctx->limit)(&info, cjmpL, cjmpR, ctx->cslope); in FVRHSFunctionfast2()
561 for (j = 0; j < dof; j++) ctx->cslope[j] /= hx; /* rescale to a slope */ in FVRHSFunctionfast2()
562 for (j = 0; j < dof; j++) { in FVRHSFunctionfast2()
564 for (k = 0; k < dof; k++) tmp += ctx->R[j + k * dof] * ctx->cslope[k]; in FVRHSFunctionfast2()
565 slope[i * dof + j] = tmp; in FVRHSFunctionfast2()
574 uL = &ctx->uLR[0]; in FVRHSFunctionfast2()
575 uR = &ctx->uLR[dof]; in FVRHSFunctionfast2()
576 for (j = 0; j < dof; j++) { in FVRHSFunctionfast2()
577 uL[j] = x[(i - 1) * dof + j] + slope[(i - 1) * dof + j] * hx / 2; in FVRHSFunctionfast2()
578 uR[j] = x[(i - 0) * dof + j] - slope[(i - 0) * dof + j] * hx / 2; in FVRHSFunctionfast2()
580 …PetscCall((*ctx->physics.riemann)(ctx->physics.user, dof, uL, uR, ctx->flux, &maxspeed, ctx->xmin … in FVRHSFunctionfast2()
583 for (j = 0; j < dof; j++) f[(i - len_slow1 - len_slow2 - 1) * dof + j] -= ctx->flux[j] / hx; in FVRHSFunctionfast2()
586 for (j = 0; j < dof; j++) f[(i - len_slow1 - len_slow2) * dof + j] += ctx->flux[j] / hx; in FVRHSFunctionfast2()
590 uL = &ctx->uLR[0]; in FVRHSFunctionfast2()
591 uR = &ctx->uLR[dof]; in FVRHSFunctionfast2()
592 for (j = 0; j < dof; j++) { in FVRHSFunctionfast2()
593 uL[j] = x[(i - 1) * dof + j] + slope[(i - 1) * dof + j] * hx / 2; in FVRHSFunctionfast2()
594 uR[j] = x[(i - 0) * dof + j] - slope[(i - 0) * dof + j] * hx / 2; in FVRHSFunctionfast2()
596 …PetscCall((*ctx->physics.riemann)(ctx->physics.user, dof, uL, uR, ctx->flux, &maxspeed, ctx->xmin … in FVRHSFunctionfast2()
599 for (j = 0; j < dof; j++) f[(i - len_slow1 - len_slow2) * dof + j] += ctx->flux[j] / hx; in FVRHSFunctionfast2()
630 PetscCall(PetscFunctionListAdd(&limiters, "lax-wendroff", Limit_LaxWendroff)); in main()
631 PetscCall(PetscFunctionListAdd(&limiters, "beam-warming", Limit_BeamWarming)); in main()
642 PetscCall(PetscFunctionListAdd(&limiters, "cada-torrilhon2", Limit_CadaTorrilhon2)); in main()
643 PetscCall(PetscFunctionListAdd(&limiters, "cada-torrilhon3-r0p1", Limit_CadaTorrilhon3R0p1)); in main()
644 PetscCall(PetscFunctionListAdd(&limiters, "cada-torrilhon3-r1", Limit_CadaTorrilhon3R1)); in main()
645 PetscCall(PetscFunctionListAdd(&limiters, "cada-torrilhon3-r10", Limit_CadaTorrilhon3R10)); in main()
646 PetscCall(PetscFunctionListAdd(&limiters, "cada-torrilhon3-r100", Limit_CadaTorrilhon3R100)); in main()
654 ctx.xmin = -1.0; in main()
657 PetscCall(PetscOptionsReal("-xmin", "X min", "", ctx.xmin, &ctx.xmin, NULL)); in main()
658 PetscCall(PetscOptionsReal("-xmax", "X max", "", ctx.xmax, &ctx.xmax, NULL)); in main()
659 …PetscCall(PetscOptionsFList("-limit", "Name of flux limiter to use", "", limiters, lname, lname, s… in main()
660 …PetscCall(PetscOptionsInt("-draw", "Draw solution vector, bitwise OR of (1=initial,2=final,4=final… in main()
661 …PetscCall(PetscOptionsString("-view_final", "Write final solution in ASCII MATLAB format to given … in main()
662 …PetscCall(PetscOptionsInt("-initial", "Initial condition (depends on the physics)", "", ctx.initia… in main()
663 …PetscCall(PetscOptionsBool("-simulation", "Compare errors with reference solution", "", ctx.simula… in main()
664 PetscCall(PetscOptionsReal("-cfl", "CFL number to time step at", "", ctx.cfl, &ctx.cfl, NULL)); in main()
665 …PetscCall(PetscOptionsEnum("-bc_type", "Boundary condition", "", FVBCTypes, (PetscEnum)ctx.bctype,… in main()
666 …PetscCall(PetscOptionsBool("-recursive_split", "Split the domain recursively", "", ctx.recursive, … in main()
687 /* The names will be shown in the title bars when run with -ts_monitor_draw_solution */ in main()
693 …SetUniformCoordinates(da, ctx.xmin + 0.5 * (ctx.xmax - ctx.xmin) / Mx, ctx.xmax + 0.5 * (ctx.xmax … in main()
704 …/*-------------------------------- create index for slow parts and fast parts --------------------… in main()
708 …if (ctx.xmin + i * (ctx.xmax - ctx.xmin) / (PetscReal)Mx + 0.5 * (ctx.xmax - ctx.xmin) / (PetscRea… in main()
716 /* Create a time-stepping object */ in main()
730 …PetscReal coord = ctx.xmin + i * (ctx.xmax - ctx.xmin) / (PetscReal)Mx + 0.5 * (ctx.xmax - ctx.xmi… in main()
731 if (coord >= 0 && coord < ctx.xmin + (ctx.xmax - ctx.xmin) * 3 / 4.) in main()
733 if (coord >= ctx.xmin + (ctx.xmax - ctx.xmin) * 3 / 4.) in main()
771 mass_difference = (ctx.xmax - ctx.xmin) / (PetscScalar)Mx * (mass_final - mass_initial); in main()
780 PetscCall(PetscOptionsGetString(NULL, NULL, "-f", filename, sizeof(filename), &flg)); in main()
781 … PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_USER, "Must indicate binary file with the -f option"); in main()
786 PetscCall(VecAYPX(XR, -1, X)); in main()
790 …PetscCall(PetscPrintf(comm, "Error ||x-x_e||_1 %g ||x-x_e||_sup %g\n", (double)nrm1, (double)nrms… in main()
802 PetscCall(VecAYPX(Y, -1, X)); in main()
844 …-da_grid_x 60 -initial 1 -xmin -1 -xmax 1 -limit mc -physics_advect_a1 1 -physics_advect_a2 2 -ts_…
848 …-da_grid_x 60 -initial 1 -xmin -1 -xmax 1 -limit mc -physics_advect_a1 1 -physics_advect_a2 2 -ts_…
853 …-da_grid_x 60 -initial 1 -xmin -1 -xmax 1 -limit mc -physics_advect_a1 1 -physics_advect_a2 2 -ts_…
857 …-da_grid_x 60 -initial 1 -xmin -1 -xmax 1 -limit mc -physics_advect_a1 1 -physics_advect_a2 2 -ts_…
862 …-da_grid_x 60 -initial 1 -xmin -1 -xmax 1 -limit mc -physics_advect_a1 1 -physics_advect_a2 2 -ts_…
866 …-da_grid_x 60 -initial 1 -xmin -1 -xmax 1 -limit mc -physics_advect_a1 1 -physics_advect_a2 2 -ts_…