|
Plugins
|
Finite Volume Element Discretization for the incompressible Navier-Stokes Equation. More...
#include <navier_stokes_base.h>
Inheritance diagram for ug::NavierStokes::NavierStokesBase< TDomain >:Public Member Functions | |
| virtual std::string | disc_type () const =0 |
| returns string identifying disc type | |
| virtual SmartPtr< CplUserData< number, dim > > | kinematic_viscosity ()=0 |
| returns kinematic viscosity | |
| virtual bool | requests_local_time_series () |
| returns if local time series is needed | |
| void | set_exact_jacobian (bool bExactJacobian) |
| sets if the exact jacobian is computed (fixpoint approximation else) | |
| 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) | |
| virtual void | set_kinematic_viscosity (SmartPtr< CplUserData< number, dim > > user)=0 |
| sets the kinematic viscosity | |
| void | set_kinematic_viscosity (number val) |
| sets the kinematic viscosity | |
| virtual void | set_source (SmartPtr< CplUserData< MathVector< dim >, dim > > user)=0 |
| sets the source function | |
| void | set_source (const std::vector< number > &vSource) |
| sets the source function | |
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 ug::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[]) |
| TLeaf & | 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_elem_loop () |
| 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< typename TDomain > | |
| 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< typename TDomain > | |
| 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 |
| virtual void | prepare_setting (const std::vector< LFEID > &vLfeID, bool bNonRegularGrid)=0 |
| 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 ug::IElemEstimatorFuncs< IElemDisc< TDomain >, TDomain > | |
| TLeaf & | 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. | |
Static Public Attributes inherited from ug::IElemDisc< TDomain > | |
| static const int | dim |
Static Public Attributes inherited from ug::IElemAssembleFuncs< IElemDisc< TDomain >, TDomain > | |
| static const int | dim |
Static Public Attributes inherited from ug::IElemError< typename TDomain > | |
| static const int | dim |
Static Public Attributes inherited from ug::IElemDiscBase< typename TDomain > | |
| static const int | dim |
Static Public Attributes inherited from ug::IElemEstimatorFuncs< IElemDisc< TDomain >, TDomain > | |
| static const int | dim |
Protected Types | |
| typedef IElemDisc< TDomain > | base_type |
| Base class type. | |
| typedef NavierStokesBase< TDomain > | this_type |
| own type | |
Finite Volume Element Discretization for the incompressible Navier-Stokes Equation.
This class implements the IElemDisc interface to provide element local assemblings for the incompressible Navier-Stokes equation.
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) \right) + \nabla \cdot \left( \rho \vec{u} \vec{u}^T \right) + \nabla p &= \vec{f}\\ \nabla \cdot (\rho \vec{u}) &= 0 \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.
In component-wise notation, the equation reads
\begin{align*} \frac{\partial \rho u_{d_1}}{\partial t} - \sum_{d_2 = 1}^{\text{dim}} \frac{\partial}{\partial x_{d_2}} \left( \rho \nu \left( \frac{\partial u_{d_1}}{\partial x_{d_2}} + \frac{\partial u_{d_2}}{\partial x_{d_1}} \right)\right) + \sum_{d_2 = 1}^{\text{dim}} \frac{\partial}{\partial x_{d_2}} \left( \rho u_{d_2} u_{d_1} \right) + \frac{\partial p}{\partial x_{d_1}} &= f_{d_1} \qquad \forall \quad d_1 = 1, \dots, \text{dim}\\ \sum_{d = 1}^{\text{dim}} \frac{\partial \rho u_d}{\partial x_{d}} &= 0 \end{align*}
Since the density is assumed to be constant, it set to \(\rho \equiv 1 \).
The finite volume element discretization uses a dual mesh consisting of Control Volumes \(\mathcal{B}_h\), that cover the domain. For each control volume \(B \in \mathcal{B}_h\) the following equation is solved
\begin{align*} \frac{\partial}{\partial t} \int_{B} \begin{pmatrix} \vec{u} \\ 0 \end{pmatrix} dV + \int_{\partial B} \begin{pmatrix} - \rho \nu \left(\nabla \vec{u} + (\nabla \vec{u})^T \right) \vec{n} + \vec{u} \vec{u}^T \vec{n} + p \vec{n} \\ \vec{u} \vec{n} \end{pmatrix} dS = \int_B \begin{pmatrix} \vec{f} \\ 0 \end{pmatrix} dV \end{align*}
where \( \vec{n}\) is the outer normal on the boundary of the control volume and \(\int_B \cdot \; dV \) and \(\int_{\partial B} \cdot \; dS \) are the integration over the control volume and the integration over the boundary of the control volume. By $f B_{sh} $f is denoted the control volume associated to a shape function (i.e. vertices, since we use P1 Lagrange ansatz functions), T usually denotes the Element.
In order to number the local unknowns, we use the following notation.
The access to local unknowns of local vectors and matrices is now given by
\begin{align*} &\mathbf{d}(d_1, sh), \mathbf{d}(\text{\_P\_}, sh),\\ &\mathcal{J}(d_1, sh_1, d_2, sh_2), \mathcal{J}(d_1, sh_1, \text{\_P\_}, sh_2), \end{align*}
| TDomain | Domain |
| TAlgebra | Algebra |
|
protected |
Base class type.
|
protected |
own type
| ug::NavierStokes::NavierStokesBase< TDomain >::NavierStokesBase | ( | const char * | functions, |
| const char * | subsets | ||
| ) |
Constructor (setting default values)
| ug::NavierStokes::NavierStokesBase< TDomain >::NavierStokesBase | ( | const std::vector< std::string > & | vFct, |
| const std::vector< std::string > & | vSubset | ||
| ) |
Constructor (setting default values)
|
pure virtual |
returns string identifying disc type
Implemented in ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >, ug::NavierStokes::NavierStokesFE< TDomain >, ug::NavierStokes::NavierStokesFV< TDomain >, ug::NavierStokes::NavierStokesFV1< TDomain >, ug::NavierStokes::NavierStokesFVCR< TDomain >, ug::NavierStokes::CompressibleNavierStokesBase< TDomain >, and ug::NavierStokes::IncompressibleNavierStokesBase< TDomain >.
|
pure virtual |
returns kinematic viscosity
Implemented in ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >, ug::NavierStokes::NavierStokesFE< TDomain >, ug::NavierStokes::NavierStokesFV< TDomain >, ug::NavierStokes::NavierStokesFV1< TDomain >, ug::NavierStokes::NavierStokesFVCR< TDomain >, ug::NavierStokes::CompressibleNavierStokesBase< TDomain >, and ug::NavierStokes::IncompressibleNavierStokesBase< TDomain >.
|
inlinevirtual |
returns if local time series is needed
Reimplemented from ug::IElemDiscBase< typename TDomain >.
Reimplemented in ug::NavierStokes::CompressibleNavierStokesBase< TDomain >, ug::NavierStokes::NavierStokesFV< TDomain >, and ug::NavierStokes::IncompressibleNavierStokesBase< TDomain >.
|
inline |
sets if the exact jacobian is computed (fixpoint approximation else)
References ug::NavierStokes::NavierStokesBase< TDomain >::m_bFullNewtonFactor.
|
inline |
| void ug::NavierStokes::NavierStokesBase< TDomain >::set_kinematic_viscosity | ( | number | val | ) |
sets the kinematic viscosity
This method sets the kinematic viscosity value.
| [in] | data | kinematic Viscosity |
References make_sp().
|
pure virtual |
sets the kinematic viscosity
This method sets the kinematic viscosity value.
| [in] | data | kinematic Viscosity |
Implemented in ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >, ug::NavierStokes::NavierStokesFE< TDomain >, ug::NavierStokes::NavierStokesFV< TDomain >, ug::NavierStokes::NavierStokesFV1< TDomain >, ug::NavierStokes::NavierStokesFVCR< TDomain >, ug::NavierStokes::CompressibleNavierStokesBase< TDomain >, and ug::NavierStokes::IncompressibleNavierStokesBase< TDomain >.
| void ug::NavierStokes::NavierStokesBase< TDomain >::set_source | ( | const std::vector< number > & | vSource | ) |
sets the source function
This method sets the source value. A zero value is assumed as default.
| [in] | data | source data |
|
pure virtual |
sets the source function
This method sets the source value. A zero value is assumed as default.
| [in] | data | source data |
Implemented in ug::NavierStokes::CompressibleNavierStokesFV1< TDomain >, ug::NavierStokes::NavierStokesFE< TDomain >, ug::NavierStokes::NavierStokesFV< TDomain >, ug::NavierStokes::NavierStokesFV1< TDomain >, ug::NavierStokes::NavierStokesFVCR< TDomain >, ug::NavierStokes::CompressibleNavierStokesBase< TDomain >, and ug::NavierStokes::IncompressibleNavierStokesBase< TDomain >.
|
static |
World dimension.
|
protected |
factor for exact jacobian, (1 for exact jacobian, 0 for fix point)
Referenced by ug::NavierStokes::NavierStokesBase< TDomain >::set_exact_jacobian(), and ug::NavierStokes::NavierStokesBase< TDomain >::set_exact_jacobian().