Lines Matching +full:docs +full:- +full:review
3 …-dependent Navier-Stokes equations of compressible gas dynamics in a static Eulerian three-dimensi…
4 Moreover, the Navier-Stokes example has been developed using PETSc, so that the pointwise physics (…
6 ## The Navier-Stokes Equations
9 The compressible Navier-Stokes equations in conservative form are
14 …abla \cdot \left( \frac{\bm{U} \otimes \bm{U}}{\rho} + P \bm{I}_3 -\bm\sigma \right) - \rho \bm{b}…
15 …t} + \nabla \cdot \left( \frac{(E + P)\bm{U}}{\rho} -\bm{u} \cdot \bm{\sigma} - k \nabla T \right)…
17 $$ (eq-ns)
19 …stress tensor, with $\mu$ the dynamic viscosity coefficient, and $\lambda = - 2/3$ the Stokes hypo…
20 In equations {eq}`eq-ns`, $\rho$ represents the volume mass density, $U$ the momentum density (defi…
23 P = \left( {c_p}/{c_v} -1\right) \left( E - {\bm{U}\cdot\bm{U}}/{(2 \rho)} \right) \, ,
24 $$ (eq-state)
28 The system {eq}`eq-ns` can be rewritten in vector form
31 \frac{\partial \bm{q}}{\partial t} + \nabla \cdot \bm{F}(\bm{q}) -S(\bm{q}) = 0 \, ,
32 $$ (eq-vector-ns)
34 for the state variables 5-dimensional vector
62 - \bm{\sigma} \\
63 - \bm{u} \cdot \bm{\sigma} - k \nabla T
72 $$ (eq-ns-flux)
83 We use tensor-product bases $\psi_{kji} = h_i(X_0)h_j(X_1)h_k(X_2)$.
85 To obtain a finite element discretization, we first multiply the strong form {eq}`eq-vector-ns` by …
88 … \left(\frac{\partial \bm{q}_N}{\partial t} + \nabla \cdot \bm{F}(\bm{q}_N) - \bm{S}(\bm{q}_N) \ri…
97 \int_{\Omega} \bm v \cdot \left( \frac{\partial \bm{q}_N}{\partial t} - \bm{S}(\bm{q}_N) \right) \…
98 - \int_{\Omega} \nabla \bm v \!:\! \bm{F}(\bm{q}_N)\,dV & \\
102 $$ (eq-weak-vector-ns)
111 <!-- TODO: This should be reframed in terms of PETSc TS's F(t, u, \dot u) = G(t, u) rather than spe…
115 #### Explicit Time-Stepping Method
116 <!-- TODO: This should talk about the mass operator and the options associated with it (i.e. lumpe…
117 …formulation is solved with the adaptive Runge-Kutta-Fehlberg (RKF4-5) method by default (any expli…
138 f(t^n, \bm{q}_N^n) = - [\nabla \cdot \bm{F}(\bm{q}_N)]^n + [S(\bm{q}_N)]^n \, .
141 #### Implicit Time-Stepping Method
143 …d using the option `-implicit` is solved with Backward Differentiation Formula (BDF) method by def…
148 $$ (eq-ts-implicit-ns)
156 …e integration scheme (backward difference formulas, generalized alpha, implicit Runge-Kutta, etc.).
157 Each nonlinear system {eq}`eq-ts-implicit-ns` will correspond to a weak form, as explained below.
158 …ing how difficult a given problem is to solve, we consider the Jacobian of {eq}`eq-ts-implicit-ns`,
164 …ted by the second term, which is a sort of "mass matrix", and typically well-conditioned independe…
166 Both terms are significant for time-accurate simulation and the setup costs of strong preconditione…
168 …me stepping solvers can be found in the [TS User Guide](https://petsc.org/release/docs/manual/ts/).
171 We solve {eq}`eq-weak-vector-ns` using a Galerkin discretization (default) or a stabilized method, …
173 …-dominated problems (any time the cell Péclet number is larger than 1), and those tend to blow up …
174 …llows {cite}`hughesetal2010`, which offers a comprehensive review of stabilization and shock-captu…
176 - **SUPG** (streamline-upwind/Petrov-Galerkin)
178 …ighted residual of the strong form {eq}`eq-vector-ns` is added to the Galerkin formulation {eq}`eq…
183 …\int_{\Omega} \bm v \cdot \left( \frac{\partial \bm{q}_N}{\partial t} - \bm{S}(\bm{q}_N) \right) …
184 - \int_{\Omega} \nabla \bm v \!:\! \bm{F}(\bm{q}_N)\,dV & \\
187 \nabla \cdot \bm{F} \, (\bm{q}_N) - \bm{S}(\bm{q}_N) \right) \,dV &= 0
190 $$ (eq-weak-vector-ns-supg)
192 This stabilization technique can be selected using the option `-stab supg`.
194 - **SU** (streamline-upwind)
196 …This method is a simplified version of *SUPG* {eq}`eq-weak-vector-ns-supg` which is developed for …
200 …\int_{\Omega} \bm v \cdot \left( \frac{\partial \bm{q}_N}{\partial t} - \bm{S}(\bm{q}_N) \right) …
201 - \int_{\Omega} \nabla \bm v \!:\! \bm{F}(\bm{q}_N)\,dV & \\
206 $$ (eq-weak-vector-ns-su)
208 This stabilization technique can be selected using the option `-stab su`.
210 In both {eq}`eq-weak-vector-ns-su` and {eq}`eq-weak-vector-ns-supg`, $\bm\tau \in \mathbb R^{5\time…
211 … be explained via an ansatz for subgrid state fluctuations $\tilde{\bm q} = -\bm\tau \bm r$ where …
212 …riational form can be readily expressed by differentiating $\bm F_{\text{adv}}$ of {eq}`eq-ns-flux`
219 (\diff\bm U \otimes \bm U + \bm U \otimes \diff\bm U)/\rho - (\bm U \otimes \bm U)/\rho^2 \diff\rho…
220 (E + P)\diff\bm U/\rho + (\diff E + \diff P)\bm U/\rho - (E + P) \bm U/\rho^2 \diff\rho
225 where $\diff P$ is defined by differentiating {eq}`eq-state`.
228 …m X} = \nabla_{\bm x}\bm X \cdot \bm u$, with units of reference length (non-dimensional) per seco…
241 $$ (eq-peclet)
243 For scalar advection-diffusion, the stabilization is a scalar
247 $$ (eq-tau-advdiff)
249 where $\xi(\mathrm{Pe}) = \coth \mathrm{Pe} - 1/\mathrm{Pe}$ approaches 1 at large local Péclet num…
250 Note that $\tau$ has units of time and, in the transport-dominated limit, is proportional to elemen…
251 For advection-diffusion, $\bm F(q) = \bm u q$, and thus the SU stabilization term is
255 $$ (eq-su-stabilize-advdiff)
257 where the term in parentheses is a rank-1 diffusivity tensor that has been pulled back to the refer…
258 See {cite}`hughesetal2010` equations 15-17 and 34-36 for further discussion of this formulation.
260 #### Navier-Stokes $\tau$ definition
262 For the Navier-Stokes and Euler equations, {cite}`whiting2003hierarchical` defines a $5\times 5$ di…
267 The Navier-Stokes code in this example uses the following formulation for $\tau_c$, $\tau_m$, $\tau…
287 #### Advection-Diffusion $\tau$ definition
289 For Advection-Diffusion, we first examine a 1D definition given by:
295 …oefficients, $h$ is the element length, and $\textrm{minreg}_n \{x_j\} = (\sum_j x_j^{-n})^{-1/n}$.
296 …mpatible with higher dimensional domains, we use a similar system from the Navier-Stokes equations.
302 …^2 \bm u \cdot (\bm u \cdot \bm g) + \left(C_d \kappa\right)^2 \Vert \bm g \Vert_F^2\right]^{-1/2}
307 The scaling coefficients are set via `-Ctau_t`, `-Ctau_a`, and `-Ctau_d`, respectively.
314 $$ (eq-tau-conservative)
321 \Lambda_i = [u_i - a, u_i, u_i, u_i, u_i+a],
322 $$ (eq-eigval-advdiff)
330 $$ (eq-wavespeed)
336 …three types of problems/physical models that can be selected at run time via the option `-problem`.
337 …-advection`, the problem of the transport of energy in a uniform vector velocity field, {ref}`prob…
341 The strong residual in the SUPG operator in {eq}`eq-weak-vector-ns-supg` and {eq}`eq-weak-vector-ns…
344 Additionally, libCEED does not currently support calculating double-derivatives.
355 For compressible Navier-Stokes, this requires projecting 12 scalars-per-node: 4 conserved scalars (…
357 This method can be selected with `-div_diff_flux_projection_method indirect`.
369 …m F_{\text{diff}}(\bm{q}_N)$ at quadrature points, so we apply integration-by-parts to achieve a c…
373 - \int_{\Omega} \nabla \bm v \!:\! \bm{F}_{\text{diff}}(\bm{q}_N)\,dV
378 For compressible Navier-Stokes, this means only projecting only 4 scalars-per-node.
380 The projection can be enabled using `-div_diff_flux_projection_method direct`.
384 The linear solve for the projection can be controlled via `-div_diff_flux_projection_ksp*` flags.
388 …-resolved (the smallest length scale resolved by the grid is much larger than the smallest physica…
389 This is known as large-eddy simulation (LES), as only the "large" scales of turbulence are resolved.
390 This filtering operation results in an extra stress-like term, $\bm{\tau}^r$, representing the effe…
394 \frac{\partial \bm{\overline q}}{\partial t} + \nabla \cdot \bm{\overline F}(\bm{\overline q}) -S(\…
395 $$ (eq-vector-les)
407 $$ (eq-les-flux)
414 :::{list-table} SGS Model Options
415 :header-rows: 1
417 * - Option
418 - Description
419 - Default value
420 - Unit
422 * - `-sgs_model_type`
423 - Type of subgrid stress model to use. Currently only `data_driven` is available
424 - `none`
425 - string
428 (sgs-dd-model)=
429 #### Data-Driven SGS Model
431 The data-driven SGS model implemented here uses a small neural network to compute the SGS term.
436 The slope parameter for the Leaky ReLU function is set via `-sgs_model_dd_leakyrelu_alpha`.
437 …he network are assumed to be normalized on a min-max scale, so they must be rescaled by the origin…
438 Parameters for the neural network are put into files in a directory found in `-sgs_model_dd_paramet…
441 Note that the weight coefficients are assumed to be in column-major order.
445 The data-driven model parameters in the examples directory are not accurate and are for regression …
448 ##### Data-Driven Model Using External Libraries
450 There are two different modes for using the data-driven model: fused and sequential.
454 To use the fused mode, set `-sgs_model_dd_implementation fused`.
462 To specify the path to the PyTorch model file, use `-sgs_model_dd_torch_model_path`.
463 … automatically from the libCEED backend chosen, but can be overridden with `-sgs_model_dd_torch_mo…
464 …ce on host while using a GPU libCEED backend (e.g. `/gpu/cuda`), then host-to-device transfers (an…
466 The sequential mode is available using a libCEED based inference evaluation via `-sgs_model_dd_impl…
468 :::{list-table} Data-driven SGS Model Options
469 :header-rows: 1
471 * - Option
472 - Description
473 - Default value
474 - Unit
476 * - `-sgs_model_dd_leakyrelu_alpha`
477 - Slope parameter for Leaky ReLU activation function. `0` corresponds to normal ReLU
478 - 0
479 -
481 * - `-sgs_model_dd_parameter_dir`
482 - Path to directory with data-driven model parameters (weights, biases, etc.)
483 - `./dd_sgs_parameters`
484 - string
486 * - `-sgs_model_dd_model_implementation`
487 …- Which computational implementation to use for SGS DD model (`fused`, `sequential_ceed`, `sequent…
488 - `fused`
489 - string
491 * - `-sgs_model_dd_torch_model_path`
492 - Path to the PyTorch `*.pt` file containing the DD inference model
493 -
494 - string
496 * - `-sgs_model_dd_torch_model_device`
497 - What hardware to perform the model inference on (`cpu`, `cuda`, `hip`, `xpu`)
498 - Default matches the libCEED backend
499 - string
506 (bc-stg)=
510 Below follows a re-description of the formulation to match the present notation, and then a descrip…
521 \bm{\hat{x}}^n &= \left[(x - U_0 t)\max(2\kappa_{\min}/\kappa^n, 0.1) , y, z \right]^T
537 \kappa^n = \kappa_{\min} (1 + \alpha)^{n-1} \ , \quad \forall n=1, 2, ... , N
543 …appa^n}{\sum^N_{n=1} E(\kappa^n)\Delta \kappa^n} \ ,\quad \Delta \kappa^n = \kappa^n - \kappa^{n-1}
549 f_\eta = \exp \left[-(12\kappa /\kappa_\eta)^2 \right], \quad
550 f_\mathrm{cut} = \exp \left( - \left [ \frac{4\max(\kappa-0.9\kappa_\mathrm{cut}, 0)}{\kappa_\mathr…
553 …turbulent dissipation frequency, and is given as $2\pi (\nu^3/\varepsilon)^{-1/4}$ with $\nu$ the …
560 …elocity, with the option of weakly enforcing either density or temperature using the `-weakT` flag.
587 lt --Calc-->ke --Calc-->kn
588 y --Calc-->ke
594 u0 --Copy--> u0C[U0]
600 ubar --Copy--> ubarC;
601 y --Copy--> yC;
602 lt --Copy--> ltC;
603 eps --Copy--> epsC;
605 rand --Copy--> randC;
606 rand --> N --Calc--> kn;
607 Rij --Calc--> Cij[C_ij]
611 The spatially-varying terms are then evaluated at each quadrature point on-the-fly, either by inter…
632 | ----------------- | -------- | -------------- | --------- |
647 …e STG boundary condition, the `-bc_inflow` option should be set to the boundary faces that need th…
648 The `-stg_use` flag is then used to enable/disable applying STG to those faces.
650 :::{list-table} STG Runtime Options
651 :header-rows: 1
653 * - Option
654 - Description
655 - Default value
656 - Unit
658 * - `-stg_use`
659 - Enable STG for `bc_inflow` faces
660 - `false`
661 -
663 * - `-stg_inflow_path`
664 - Path to the STGInflow file
665 - `./STGInflow.dat`
666 -
668 * - `-stg_rand_path`
669 - Path to the STGRand file
670 - `./STGRand.dat`
671 -
673 * - `-stg_alpha`
674 - Growth rate of the wavemodes
675 - `1.01`
676 -
678 * - `-stg_u0`
679 - Convective velocity, $U_0$
680 - `0.0`
681 - `m/s`
683 * - `-stg_mean_only`
684 - Only impose the mean velocity (no fluctutations)
685 - `false`
686 -
688 * - `-stg_strong`
689 - Strongly enforce the STG inflow boundary condition
690 - `false`
691 -
693 * - `-stg_fluctuating_IC`
694 - "Extrude" the fluctuations through the domain as an initial condition
695 - `false`
696 -
698 * - `-stg_dx`
699 …- Set the element size in the x direction. Default is calculated for box meshes, assuming equispac…
700 -
701 - `m`
703 * - `-stg_h_scale_factor`
704 - Scale element size for cutoff frequency calculation
705 - $1/p$
706 -
719 S(\bm{q}) = -\sigma(\bm{x})\left.\frac{\partial \bm{q}}{\partial \bm{Y}}\right\rvert_{\bm{q}} \bm{Y…
722 …- P_\mathrm{ref}, \bm{0}, 0]^T$, and $\sigma(\bm{x})$ is a linear ramp starting at `-idl_start` wi…
723 The damping is defined in terms of a pressure-primitive anomaly $\bm Y'$ converted to conservative …
724 $P_\mathrm{ref}$ has a default value equal to `-reference_pressure` flag, with an optional flag `-i…
726 :::{list-table} IDL Runtime Options
727 :header-rows: 1
729 * - Option
730 - Description
731 - Default value
732 - Unit
734 * - `-idl_decay_time`
735 - Characteristic timescale of the pressure deviance decay. The timestep is good starting point
736 - `-1` (disabled)
737 - `s`
739 * - `-idl_start`
740 - Start of IDL in the x direction
741 - `0`
742 - `m`
744 * - `-idl_length`
745 - Length of IDL in the positive x direction
746 - `0`
747 - `m`
749 * - `-idl_pressure`
750 - Pressure used for IDL reference pressure
751 - `-reference_pressure`
752 - `Pa`
756 <!-- TODO: Move the Riemann/freestream description here-->