Plugins
|
Finite Volume Element Discretization for the compressible Navier-Stokes Equation. More...
#include <compressible_navier_stokes_fv1.h>
Public Member Functions | |
template<typename TElem , typename TFVGeom > | |
void | add_def_A_elem (LocalVector &d, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
adds the stiffness part to the local defect More... | |
template<typename TElem , typename TFVGeom > | |
void | add_def_M_elem (LocalVector &d, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
adds the mass part to the local defect More... | |
template<typename TElem , typename TFVGeom > | |
void | add_jac_A_elem (LocalMatrix &J, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
adds the stiffness part to the local jacobian More... | |
template<typename TElem , typename TFVGeom > | |
void | add_jac_M_elem (LocalMatrix &J, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
adds the mass part to the local jacobian More... | |
template<typename TElem , typename TFVGeom > | |
void | add_rhs_elem (LocalVector &d, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
adds the source part to the local defect More... | |
SmartPtr< CplUserData< number, dim > > | adiabatic_index () |
returns adiabatic index More... | |
virtual std::string | disc_type () const |
returns string identifying disc type More... | |
template<typename TElem , typename TFVGeom > | |
void | fsh_elem_loop () |
finishes the element loop More... | |
SmartPtr< CplUserData< number, dim > > | kinematic_viscosity () |
returns kinematic viscosity More... | |
template<typename TFVGeom > | |
number | mach_number_blending (MathVector< dim > &UpwindVel, const TFVGeom &geo, size_t ip, const MathVector< dim > &StdVel, number kinVisco) |
computes the Mach-number blended Upwind velocity More... | |
template<typename TElem , typename TFVGeom > | |
void | prep_elem (const LocalVector &u, GridObject *elem, const ReferenceObjectID roid, const MathVector< dim > vCornerCoords[]) |
prepares the element for evaluation More... | |
template<typename TElem , typename TFVGeom > | |
void | prep_elem_loop (const ReferenceObjectID roid, const int si) |
prepares the element loop More... | |
template<typename TElem , typename TFVGeom > | |
void | prep_timestep_elem (const number time, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
prepares the element loop More... | |
void | prepare_setting (const std::vector< LFEID > &vLfeID, bool bNonRegularGrid) |
type of trial space for each function used More... | |
void | set_adiabatic_index (SmartPtr< CplUserData< number, dim > > user) |
sets the adiabatic index More... | |
void | set_kinematic_viscosity (SmartPtr< CplUserData< number, dim > > user) |
sets the kinematic viscosity More... | |
void | set_source (SmartPtr< CplUserData< MathVector< dim >, dim > > user) |
sets the source function More... | |
void | set_upwind (const std::string &name) |
sets the upwind based on a string identifier More... | |
void | set_upwind (SmartPtr< INavierStokesUpwind< dim > > spUpwind) |
sets an upwinding for the convective term of momentum equation More... | |
DataImport< MathVector< dim >, dim > | sourceSCVF () |
returns scvf source More... | |
CompressibleNavierStokesFV1 (const char *functions, const char *subsets) | |
CompressibleNavierStokesFV1 (const std::vector< std::string > &vFct, const std::vector< std::string > &vSubset) | |
Public Member Functions inherited from ug::NavierStokes::CompressibleNavierStokesBase< TDomain > | |
CompressibleNavierStokesBase (const char *functions, const char *subsets) | |
CompressibleNavierStokesBase (const std::vector< std::string > &vFct, const std::vector< std::string > &vSubset) | |
void | set_adiabatic_index (number val) |
flag if using Mach-number Blending More... | |
void | set_mach_number_blend (bool machNrBlend) |
sets if Mach-number blending is used in momentum equation More... | |
virtual bool | requests_local_time_series () |
returns if local time series is needed More... | |
Public Member Functions inherited from ug::NavierStokes::NavierStokesBase< TDomain > | |
void | set_exact_jacobian (bool bExactJacobian) |
sets if the exact jacobian is computed (fixpoint approximation else) More... | |
void | set_exact_jacobian (number fullNewtonFactor) |
NavierStokesBase (const char *functions, const char *subsets) | |
NavierStokesBase (const std::vector< std::string > &vFct, const std::vector< std::string > &vSubset) | |
void | set_kinematic_viscosity (number val) |
sets the kinematic viscosity More... | |
void | set_source (const std::vector< number > &vSource) |
sets the source function More... | |
Public Member Functions inherited from ug::IElemDisc< TDomain > | |
void | add_elem_modifier (SmartPtr< IElemDiscModifier< TDomain > > elemModifier) |
std::vector< SmartPtr< IElemDiscModifier< TDomain > > > & | get_elem_modifier () |
IElemDisc (const char *functions, const char *subsets) | |
IElemDisc (const std::vector< std::string > &vFct, const std::vector< std::string > &vSubset) | |
Public Member Functions inherited from IElemAssembleFuncs< IElemDisc< TDomain >, TDomain > | |
virtual void | add_def_A_elem (LocalVector &d, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
virtual void | add_def_A_expl_elem (LocalVector &d, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
virtual void | add_def_M_elem (LocalVector &d, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
virtual void | add_jac_A_elem (LocalMatrix &J, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
virtual void | add_jac_M_elem (LocalMatrix &J, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
virtual void | add_rhs_elem (LocalVector &rhs, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
IElemDisc< TDomain > & | asLeaf () |
void | check_roid (ReferenceObjectID roid, int discType) |
void | do_add_def_A_elem (LocalVector &d, LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
void | do_add_def_A_expl_elem (LocalVector &d, LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
void | do_add_def_M_elem (LocalVector &d, LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
void | do_add_jac_A_elem (LocalMatrix &J, LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
void | do_add_jac_M_elem (LocalMatrix &J, LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
void | do_add_rhs_elem (LocalVector &rhs, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
void | do_fsh_elem_loop () |
void | do_fsh_timestep (const number time, VectorProxyBase *u, size_t algebra_id) |
void | do_fsh_timestep_elem (const number time, LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
void | do_prep_elem (LocalVector &u, GridObject *elem, const ReferenceObjectID roid, const MathVector< dim > vCornerCoords[]) |
void | do_prep_elem_loop (const ReferenceObjectID roid, const int si) |
void | do_prep_timestep (number future_time, const number time, VectorProxyBase *u, size_t algebra_id) |
void | do_prep_timestep_elem (const number time, LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
virtual void | fsh_timestep (number time, VectorProxyBase *u) |
virtual void | fsh_timestep_elem (const number time, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
IElemAssembleFuncs () | |
virtual void | post_assemble_loop () |
virtual void | prep_assemble_loop () |
virtual void | prep_elem (const LocalVector &u, GridObject *elem, const ReferenceObjectID roid, const MathVector< dim > vCornerCoords[]) |
virtual void | prep_elem_loop (const ReferenceObjectID roid, const int si) |
virtual void | prep_timestep (number future_time, number time, VectorProxyBase *u) |
virtual void | prep_timestep_elem (const number time, const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
void | set_roid (ReferenceObjectID id, int discType) |
virtual | ~IElemAssembleFuncs () |
Public Member Functions inherited from ug::IElemError< class > | |
IElemError (const char *functions, const char *subsets) | |
IElemError (const std::vector< std::string > &vFct, const std::vector< std::string > &vSubset) | |
Public Member Functions inherited from ug::IElemDiscBase< class > | |
SmartPtr< ApproximationSpace< TDomain > > | approx_space () |
ConstSmartPtr< ApproximationSpace< TDomain > > | approx_space () const |
void | check_setup (bool bNonRegularGrid) |
void | clear_imports () |
SmartPtr< TDomain > | domain () |
ConstSmartPtr< TDomain > | domain () const |
const FunctionGroup & | function_group () const |
ConstSmartPtr< FunctionPattern > | function_pattern () const |
IDataImport< dim > & | get_import (size_t i) |
IElemDiscBase (const char *functions="", const char *subsets="") | |
IElemDiscBase (const std::vector< std::string > &vFct, const std::vector< std::string > &vSubset) | |
bool | is_time_dependent () const |
bool | local_time_series_needed () |
const LocalVectorTimeSeries * | local_time_solutions () const |
const FunctionIndexMapping & | map () const |
number | mass_scale () const |
number | mass_scale () const |
number | mass_scale (const size_t timePoint) const |
number | mass_scale (const size_t timePoint) const |
const std::vector< number > & | mass_scales () const |
const std::vector< number > & | mass_scales () const |
size_t | num_fct () const |
size_t | num_imports () const |
size_t | num_subsets () const |
void | register_import (IDataImport< dim > &Imp) |
void | set_approximation_space (SmartPtr< ApproximationSpace< TDomain > > approxSpace) |
void | set_functions (const std::string &functions) |
void | set_functions (const std::vector< std::string > &functions) |
void | set_stationary () |
void | set_stationary (bool bStationaryForced=true) |
void | set_subsets (const std::string &subsets) |
void | set_subsets (const std::vector< std::string > &subsets) |
void | set_time_dependent (LocalVectorTimeSeries &locTimeSeries, const std::vector< number > &vScaleMass, const std::vector< number > &vScaleStiff) |
void | set_time_independent () |
void | set_time_point (const size_t timePoint) |
number | stiff_scale () const |
number | stiff_scale () const |
number | stiff_scale (const size_t timePoint) const |
number | stiff_scale (const size_t timePoint) const |
const std::vector< number > & | stiff_scales () const |
const std::vector< number > & | stiff_scales () const |
TDomain::subset_handler_type & | subset_handler () |
const TDomain::subset_handler_type & | subset_handler () const |
const std::vector< std::string > & | symb_fcts () const |
const std::vector< std::string > & | symb_subsets () const |
number | time () const |
size_t | time_point () const |
virtual int | type () const |
virtual bool | use_hanging () const |
virtual | ~IElemDiscBase () |
Public Member Functions inherited from IElemEstimatorFuncs< IElemDisc< TDomain >, TDomain > | |
IElemDisc< TDomain > & | asLeaf () |
void | check_roid (ReferenceObjectID roid, int discType) |
virtual void | compute_err_est_A_elem (const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[], const number &scale) |
virtual void | compute_err_est_M_elem (const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[], const number &scale) |
virtual void | compute_err_est_rhs_elem (GridObject *elem, const MathVector< dim > vCornerCoords[], const number &scale) |
void | do_compute_err_est_A_elem (LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[], const number &scale) |
void | do_compute_err_est_M_elem (LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[], const number &scale) |
void | do_compute_err_est_rhs_elem (GridObject *elem, const MathVector< dim > vCornerCoords[], const number &scale) |
void | do_fsh_err_est_elem_loop () |
void | do_prep_err_est_elem (LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
void | do_prep_err_est_elem_loop (const ReferenceObjectID roid, const int si) |
virtual SmartPtr< IErrEstData< TDomain > > | err_est_data () |
bool | err_est_enabled () const |
virtual void | fsh_err_est_elem_loop () |
IElemEstimatorFuncs () | |
virtual void | prep_err_est_elem (const LocalVector &u, GridObject *elem, const MathVector< dim > vCornerCoords[]) |
virtual void | prep_err_est_elem_loop (const ReferenceObjectID roid, const int si) |
void | set_error_estimator (SmartPtr< IErrEstData< TDomain > > ee) |
void | set_roid (ReferenceObjectID id, int discType) |
virtual | ~IElemEstimatorFuncs () |
Static Public Attributes | |
static const int | dim = base_type::dim |
World dimension. More... | |
Static Public Attributes inherited from ug::NavierStokes::CompressibleNavierStokesBase< TDomain > | |
static const int | dim = base_type::dim |
World dimension. More... | |
Static Public Attributes inherited from ug::NavierStokes::NavierStokesBase< TDomain > | |
static const int | dim = base_type::dim |
World dimension. More... | |
Static Public Attributes inherited from ug::IElemDisc< TDomain > | |
static const int | dim |
Static Public Attributes inherited from IElemAssembleFuncs< IElemDisc< TDomain >, TDomain > | |
static const int | dim |
Static Public Attributes inherited from ug::IElemError< class > | |
static const int | dim |
Static Public Attributes inherited from ug::IElemDiscBase< class > | |
static const int | dim |
Static Public Attributes inherited from IElemEstimatorFuncs< IElemDisc< TDomain >, TDomain > | |
static const int | dim |
Protected Types | |
typedef CompressibleNavierStokesBase< TDomain > | base_type |
Base class type. More... | |
typedef CompressibleNavierStokesFV1< TDomain > | this_type |
own type More... | |
Protected Types inherited from ug::NavierStokes::CompressibleNavierStokesBase< TDomain > | |
typedef NavierStokesBase< TDomain > | base_type |
Base class type. More... | |
typedef CompressibleNavierStokesBase< TDomain > | this_type |
own type More... | |
Protected Types inherited from ug::NavierStokes::NavierStokesBase< TDomain > | |
typedef IElemDisc< TDomain > | base_type |
Base class type. More... | |
typedef NavierStokesBase< TDomain > | this_type |
own type More... | |
Static Protected Attributes | |
static const size_t | _P_ = dim |
abbreviation for pressure More... | |
static const size_t | _Rho_ = dim + 1 |
abbreviation for density More... | |
Private Member Functions | |
void | register_all_funcs (bool bHang) |
template<typename TElem , typename TFVGeom > | |
void | register_func () |
Additional Inherited Members | |
Public Types inherited from ug::IElemDisc< TDomain > | |
typedef IElemAssembleFuncs< IElemDisc< TDomain >, TDomain > | assemble_base_type |
typedef IElemError< TDomain > | base_type |
typedef TDomain | domain_type |
typedef IElemEstimatorFuncs< IElemDisc< TDomain >, TDomain > | estimator_base_type |
Public Types inherited from IElemAssembleFuncs< IElemDisc< TDomain >, TDomain > | |
typedef TDomain | domain_type |
typedef IElemDisc< TDomain > | leaf_type |
Public Types inherited from ug::IElemError< class > | |
typedef TDomain | domain_type |
typedef IElemEstimatorFuncs< IElemDisc< TDomain >, TDomain > | estimator_base_type |
Public Types inherited from ug::IElemDiscBase< class > | |
typedef TDomain | domain_type |
typedef TDomain::position_type | position_type |
Public Types inherited from IElemEstimatorFuncs< IElemDisc< TDomain >, TDomain > | |
typedef TDomain | domain_type |
typedef IElemDisc< TDomain > | leaf_type |
Finite Volume Element Discretization for the compressible Navier-Stokes Equation.
This class implements the IElemDisc interface to provide element local assemblings for the compressible Navier-Stokes equation for an thermical and caloric ideal gas.
The unknowns of the equation are:
The equation takes the form
\begin{align*} \frac{\partial \rho \vec{u}}{\partial t} - \nabla \left( \rho \nu (\nabla \vec{u} + (\nabla \vec{u})^T) - \frac{2}{3} Id * \nabla \vec{u} \right) + \nabla \cdot \left( \rho \vec{u} \vec{u}^T \right) + \nabla p &= \vec{f}_m\\ \frac{\partial \rho}{\partial t} + \nabla \cdot (\rho \vec{u}) &= f_c \\ \frac{1}{\gamma - 1}\frac{\partial \rho}{\partial t} + \frac{1}{2} \frac{\partial \rho \vec{u}^2}{\partial t} + \frac{\gamma}{\gamma - 1} \nabla \cdot (\rho \vec{u}^T) + \frac{1}{2} \nabla \cdot (\rho \vec{u}^2 \vec{u}^T) - \nabla \left( \vec{u}^T \cdot \tau \right) &= f_e \\ \end{align*}
with
The first equation models the conservation of momentum and is therefore referred to as the Momentum equation. The second equation models the conservation of mass and is known as the Mass continuity equation or simply Continuity equation. The third equation models the conservation of energy and is therefore referred to as the Energy equation.
These three equations are normalized by reference quantities (c.f. A. Gordner 'Numerische Simulation nichtlinearer Aeroakustik bei kleinen Machzahlen'),
with
TDomain | Domain |
TAlgebra | Algebra |
|
protected |
Base class type.
|
protected |
own type
ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::CompressibleNavierStokesFV1 | ( | const char * | functions, |
const char * | subsets | ||
) |
Constructor (setting default values)
References ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::init().
ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::CompressibleNavierStokesFV1 | ( | const std::vector< std::string > & | vFct, |
const std::vector< std::string > & | vSubset | ||
) |
Constructor (setting default values)
References ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::init().
void ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::add_def_A_elem | ( | LocalVector & | d, |
const LocalVector & | u, | ||
GridObject * | elem, | ||
const MathVector< dim > | vCornerCoords[] | ||
) |
adds the stiffness part to the local defect
This function adds the contribution of one element of the stiffness part to the local defect.
The local solution \( \mathbf{u} \) passes the p1 conform lagrange dofs of the element, i.e. the unknowns in the corners of the element.
Define the functional matrix of the velocity at some integration point by
\begin{align*} \nabla \vec{u}|_{ip} := \left(\left.\frac{\partial u_i}{\partial x_j} \right|_{ip}\right)_{i,j = 1,\dots,dim} \end{align*}
and the diffusive flux to
\begin{align*} \vec{f}^{\text{diff}}|_{ip} := - \nu \left( \nabla \vec{u}|_{ip} + (\nabla \vec{u}|_{ip})^T) \right) \cdot \vec{n}_{ip} \end{align*}
Suppose that a procedure producing some upwind velocity \( \vec{u}^{up}\) is available for each ip. Then, we define the convective flux to
\begin{align*} \vec{f}^{\text{conv}}|_{ip} := \left( \vec{u}^{up}|_{ip} \vec{u}^{up}|_{ip}^T \right) \cdot \vec{n}_{ip} \end{align*}
Define the pressure at the ip to be interpolated as
\begin{align*} p|_{ip} := \sum_{sh=0}^{n_{sh}-1} \mathbf{u}(\text{\_P\_}, sh) \phi_{sh} \end{align*}
The contribution added for the momentum equation is:
\begin{align*} \forall d=0,\dots, \text{dim}-1, \quad \forall sh=0,\dots, n_{sh}-1: \qquad \mathbf{d}(d, sh) &\text{ +=} \int_{\partial B_{sh} \cap T} \vec{f}^{\text{diff}}_d + \vec{f}^{\text{conv}}_d + p \vec{n}_d \;dS \approx \sum_{i=1}^{n_{\partial B_{sh}}} |\partial B^i_{sh} \cap T| \left( \left.\vec{f}^{\text{diff}}_d\right|_{ip} + \left.\vec{f}^{\text{conv}}_d\right|_{ip} + p|_{ip} \vec{n}|_{ip} \right)\\ \end{align*}
Assume, that some stabilized velocity \( \vec{u}^{stab}\) is available for each ip. Then, the contribution added for the continuity equation is:
\begin{align*} \forall sh=0,\dots, n_{sh}-1: \qquad \mathbf{d}(\text{\_P\_}, sh) &\text{ +=} \int_{\partial B_{sh} \cap T} \vec{u}^{stab} \cdot \vec{n} \;dS \approx \sum_{i=1}^{n_{\partial B_{sh}}} |\partial B^i_{sh} \cap T| \; \vec{u}^{stab}|_{ip} \cdot \vec{n}|_{ip}\\ \end{align*}
[in,out] | d | local defect with values added on output |
[in] | u | local vector of current solution |
[in] | time | current time step |
m_bLaplace)
References dim, ug::GeomProvider< class >::get(), ug::MatVecMult(), ug::TransposedMatVecMultAdd(), UG_ASSERT, UG_THROW, ug::NavierStokes::INavierStokesUpwind< dim >::upwind_vel(), ug::VecProd(), ug::VecScaleAdd(), and ug::VecSet().
void ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::add_def_M_elem | ( | LocalVector & | d, |
const LocalVector & | u, | ||
GridObject * | elem, | ||
const MathVector< dim > | vCornerCoords[] | ||
) |
adds the mass part to the local defect
This function adds the contribution of one element of the mass part to the local defect.
The local solution \( \mathbf{u} \) passes the p1 conform lagrange dofs of the element, i.e. the unknowns in the corners of the element.
The contribution added is:
\begin{align*} \forall d=0,\dots, \text{dim}-1, \quad \forall sh=0,\dots, n_{sh}-1: \qquad \mathbf{d}(d, sh) &\text{ +=} \int_{B_{sh} \cap T} \vec{u}_d dV \approx |B_{sh} \cap T| \mathbf{u}(d, sh)\\ \mathbf{d}(\text{\_P\_, sh}) &\text{ +=} 0 \end{align*}
[in,out] | d | local defect with values added on output |
[in] | u | local vector of current solution |
[in] | time | current time step |
References dim, ug::GeomProvider< class >::get(), UG_ASSERT, and ug::VecProd().
void ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::add_jac_A_elem | ( | LocalMatrix & | J, |
const LocalVector & | u, | ||
GridObject * | elem, | ||
const MathVector< dim > | vCornerCoords[] | ||
) |
adds the stiffness part to the local jacobian
This function adds the local contributions of the stiffness part to the local jacobian.
For the definition of \( \vec{f}^{\text{diff}}|_{ip}\) and \( \vec{f}^{\text{conv}}|_{ip}\) see add_def_A_elem.
The derivative of the diffusive flux is given by
\begin{align*} \frac{\partial \vec{f}^{\text{diff}}_{d_1}|_{ip}}{\partial \vec{u}_{d_2}^{sh}} &= - \nu \left( \delta_{d_1, d_2} \sum_{k=1}^{\text{dim}} \left. \frac{\partial \phi^{sh}}{\partial x_k}\right|_{ip} \vec{n_k}|_{ip} + \left. \frac{\partial \phi^{sh}}{\partial x_{d_1}}\right|_{ip} \vec{n_{d_2}}|_{ip} \right)\\ \frac{\partial \vec{f}^{\text{diff}}_{d_1}|_{ip}}{\partial p^{sh}} &= 0 \end{align*}
For the derivative of the convective term, we can use a fixpoint linearization if we only differentiate the first factor of \( \vec{f}^{\text{conv}}_{d_1} = (\vec{u}^{\text{blend}} (\vec{u}^{\text{blend}})^T)_{d_1 k} \vec{n}_k = \vec{u}^{\text{blend}}_{d_1} \sum_{k=1}^{\text{dim}} \vec{u}^{\text{blend}}_k \vec{n}_k = \vec{u}^{\text{blend}}_{d_1} \langle \vec{u}^{\text{blend}}, \vec{n} \rangle \)
This gives
\begin{align*} \frac{\partial \vec{f}^{\text{conv}}_{d_1}}{\partial \vec{u}_{d_2}^{sh}} &= \langle \vec{u}^{\text{blend}}, \vec{n} \rangle \frac{\partial \vec{u}^{\text{blend}}_{d_1}}{\partial \vec{u}_{d_2}^{sh}} &= \langle \vec{u}^{\text{blend}}, \vec{n} \rangle \left( (1-\omega) \delta_{d_1 d_2} \phi^{sh} + \omega \frac{\partial \vec{u}^{\text{up}}_{d_1}}{\partial \vec{u}_{d_2}^{sh}} \right) \end{align*}
The derivative of the pressure term is given by
\begin{align*} \frac{\partial p \vec{n}|_{ip}}{\partial \vec{u}_{d_2}^{sh}} &= 0\\ \frac{\partial p \vec{n}|_{ip}}{\partial p^{sh}} &= \phi^{sh}|_{ip} \vec{n} \end{align*}
The derivative of the continuity term is given by
\begin{align*} \frac{\partial \vec{u}^{\text{stab}}|_{ip} \vec{n}|_{ip}}{\partial \vec{u}_{d_1}^{sh}} &= \sum_{d_2=1}^{dim} \left. \frac{\partial u_{d_2}^{\text{stab}}} {\partial \vec{u}_{d_1}^{sh}} \right|_{ip} \vec{n}_{d_2}\\ \frac{\partial \vec{u}^{\text{stab}}|_{ip} \vec{n}|_{ip}}{\partial p^{sh}} &=\sum_{d_2=1}^{dim} \left. \frac{\partial u_{d_2}^{\text{stab}}} {\partial p^{sh}} \right|_{ip} \vec{n}_{d_2}\\ \end{align*}
[in,out] | J | local jacobian with values added on output |
[in] | u | local vector of current solution |
[in] | time | current time step |
void ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::add_jac_M_elem | ( | LocalMatrix & | J, |
const LocalVector & | u, | ||
GridObject * | elem, | ||
const MathVector< dim > | vCornerCoords[] | ||
) |
adds the mass part to the local jacobian
This function adds the contribution of one element of the mass part to the local jacobian.
The local solution \( \mathbf{u} \) passes the p1 conform lagrange dofs of the element, i.e. the unknowns in the corners of the element.
The contribution added is:
\begin{align*} \forall d=0,\dots, \text{dim}-1, \quad \forall sh=0,\dots, n_{sh}-1: \qquad \mathcal{J}(d, sh, d, sh) &\text{ +=} \int_{B_{sh} \cap T} dV = |B_{sh} \cap T|\\ \end{align*}
[in,out] | J | local jacobian with values added on output |
[in] | u | local vector of current solution |
[in] | time | current time step |
void ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::add_rhs_elem | ( | LocalVector & | d, |
GridObject * | elem, | ||
const MathVector< dim > | vCornerCoords[] | ||
) |
adds the source part to the local defect
This function adds the contribution of one element of the source part to the local defect, if a source function \( \vec{f} \) has been specified. Otherwise a zero source function is assumed.
The contribution added is:
\begin{align*} \forall d=0,\dots, \text{dim}-1, \quad \forall sh=0,\dots, n_{sh}-1: \qquad \mathbf{d}(d, sh) &\text{ +=} \int_{B_{sh} \cap T} \vec{f}_d \; dV \approx |B_{sh} \cap T| \; \vec{f}_d|_{ip} \end{align*}
[in,out] | d | local defect with values added on output |
[in] | time | current time step |
References dim, ug::GeomProvider< class >::get(), and UG_ASSERT.
|
inlinevirtual |
returns adiabatic index
Implements ug::NavierStokes::CompressibleNavierStokesBase< TDomain >.
References ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::m_imAdiabaticIndex, and DataImport< number, dim >::user_data().
|
inlinevirtual |
returns string identifying disc type
Implements ug::NavierStokes::CompressibleNavierStokesBase< TDomain >.
|
virtual |
finishes the element loop
Reimplemented from IElemAssembleFuncs< IElemDisc< TDomain >, TDomain >.
|
protected |
|
inlinevirtual |
returns kinematic viscosity
Implements ug::NavierStokes::CompressibleNavierStokesBase< TDomain >.
References ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::m_imKinViscosity, and DataImport< number, dim >::user_data().
|
inline |
computes the Mach-number blended Upwind velocity
This function compute the peclet blended velocity. First, the standard interpolated (no upwind) velocity is computed at the integration point of the scvf:
\begin{align*} \vec{u}^{\text{std}}_d|_{ip} := \sum_{sh=0}^{n_{sh}-1} \mathbf{u}(d, sh) \; \phi_{sh} \end{align*}
Defining the Peclet Number
\begin{align*} Pe := \frac{\vec{u}^{\text{std}}|_{ip} \cdot \vec{n}|_{ip} \; L}{\nu} \end{align*}
with the Viscosity \( \nu \) and the diffusion length \( L \), a weighting factor is computed:
\begin{align*} \omega := \frac{Pe^2}{5 + Pe^2} \end{align*}
and the blended velocity is computed to
\begin{align*} \vec{u}^{\text{blend}}|_{ip} := \omega \vec{u}^{\text{up}} + (1-\omega) \vec{u}^{\text{std}} \end{align*}
[in,out] | UpwindVel | upwind vel on entry, blended vel on exit |
[in] | scvf | SubControlVolumeFace of the Velocity |
[in] | StdVel | standard interpolation vel |
[in] | kinVisco | kinematic Viscosity at scvf |
References ug::VecDistance(), ug::VecProd(), ug::VecScaleAdd(), and ug::VecTwoNormSq().
void ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::prep_elem | ( | const LocalVector & | u, |
GridObject * | elem, | ||
const ReferenceObjectID | roid, | ||
const MathVector< dim > | vCornerCoords[] | ||
) |
prepares the element for evaluation
This function prepare a specific element for assembling. This function is always called before any of the assemble_... functions is performed. Here, the Finite Volume Geometry is prepared and global positions of the evaluation points of the DataImports are requested (and local positions in case of hanging node assembling)
[in] | elem | element to prepare the data for |
[in] | u | current local solution vector |
[in] | glob_ind | global indices of the local vector components |
References dim, ug::GeomProvider< class >::get(), subset_handler(), and UG_CATCH_THROW.
void ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::prep_elem_loop | ( | const ReferenceObjectID | roid, |
const int | si | ||
) |
prepares the element loop
This function is used to prepare the element loop. It gets called, before the loop over all elements starts. Here, general checks are performed:
In addition, local coordinates of evaluation points are requested by the DataImports in case of element-fixed points.
References dim, ug::GeomProvider< class >::get(), and UG_THROW.
void ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::prep_timestep_elem | ( | const number | time, |
const LocalVector & | u, | ||
GridObject * | elem, | ||
const MathVector< dim > | vCornerCoords[] | ||
) |
prepares the element loop
This function is used to set all reference quantities (like e.g. the reynolds number) which are necessary the formulate the dimensionless Navier-Stokes-equations
References dim, num_vertices(), and ug::VecTwoNorm().
|
virtual |
type of trial space for each function used
Implements ug::IElemDiscBase< class >.
References dim, ug::LFEID::LAGRANGE, type(), and UG_THROW.
|
private |
register utils
|
private |
register utils
|
virtual |
sets the adiabatic index
Implements ug::NavierStokes::CompressibleNavierStokesBase< TDomain >.
|
virtual |
sets the kinematic viscosity
Implements ug::NavierStokes::CompressibleNavierStokesBase< TDomain >.
|
virtual |
sets the source function
Implements ug::NavierStokes::CompressibleNavierStokesBase< TDomain >.
|
inline |
sets the upwind based on a string identifier
References ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::m_spConvUpwind.
|
inline |
sets an upwinding for the convective term of momentum equation
References ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::m_spConvUpwind.
|
inline |
returns scvf source
References ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::m_imSourceSCVF.
|
staticprotected |
abbreviation for pressure
|
staticprotected |
abbreviation for density
|
static |
World dimension.
|
protected |
factor for exact jacobian, (1 for exact jacobian, 0 for fix point)
|
protected |
Data import for adiabatic index.
Referenced by ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::adiabatic_index().
|
protected |
Data import for kinematic viscosity.
Referenced by ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::kinematic_viscosity().
|
protected |
Data import for source.
|
protected |
|
protected |
|
protected |
|
protected |
reference quantities for getting dimensionless equations
|
protected |
|
protected |
|
protected |
|
protected |
Upwinding for velocity in convective term of momentum equation.
Referenced by ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >::set_upwind().