ug4
ug::FV1IBGeometry< TElem, TWorldDim > Class Template Reference

Geometry and shape functions for 1st order Vertex-Centered Finite Volume. More...

#include <fv1ib_geom.h>

+ Inheritance diagram for ug::FV1IBGeometry< TElem, TWorldDim >:

Classes

class  BF
 boundary face More...
 
class  SCV
 sub control volume structure More...
 
class  SCVF
 Sub-Control Volume Face structure. More...
 

Public Types

typedef TElem elem_type
 type of element More...
 
typedef LagrangeP1< ref_elem_typelocal_shape_fct_set_type
 type of Shape function used More...
 
typedef reference_element_traits< TElem >::reference_element_type ref_elem_type
 type of reference element More...
 
typedef traits::scv_type scv_type
 type of SubControlVolume More...
 
typedef fv1_traits< ref_elem_type, TWorldDim > traits
 used traits More...
 

Public Member Functions

void adapt (GridObject *elem, const MathVector< TWorldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
 adapt integration points for elements cut by the inner boundary More...
 
void adapt_integration_points (GridObject *elem, const MathVector< TWorldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
 adapt integration points for elements cut by the inner boundary More...
 
void adapt_normals (GridObject *elem, const MathVector< TWorldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
 adapt normals for elements cut by the inner boundary More...
 
void add_boundary_subset (int subsetIndex)
 add subset that is interpreted as boundary subset. More...
 
const std::vector< BF > & bf (int si) const
 returns reference to vector of boundary faces for subsetIndex More...
 
const BFbf (int si, size_t i) const
 returns the boundary face i for subsetIndex More...
 
void clear_boundary_subsets ()
 reset all boundary subsets More...
 
const MathVector< worldDim > * corners () const
 get vector of the global coordinates of corners for current element More...
 
 FV1IBGeometry ()
 construct object and initialize local values and sizes More...
 
size_t num_bf () const
 number of all boundary faces More...
 
size_t num_bf (int si) const
 number of boundary faces on subset 'subsetIndex' More...
 
size_t num_boundary_subsets ()
 number of registered boundary subsets More...
 
size_t num_scv () const
 number of SubControlVolumes More...
 
size_t num_scv_ips () const
 returns number of all scv ips More...
 
size_t num_scvf () const
 number of SubControlVolumeFaces More...
 
size_t num_scvf_ips () const
 returns number of all scvf ips More...
 
size_t num_sh () const
 number of shape functions More...
 
void remove_boundary_subset (int subsetIndex)
 removes subset that is interpreted as boundary subset. More...
 
const SCVscv (size_t i) const
 const access to SubControlVolume number i More...
 
const MathVector< worldDim > * scv_global_ips () const
 returns all ips of scv as they appear in scv loop More...
 
const MathVector< dim > * scv_local_ips () const
 returns all ips of scv as they appear in scv loop More...
 
const SCVFscvf (size_t i) const
 const access to SubControlVolumeFace number i More...
 
const MathVector< worldDim > * scvf_global_ips () const
 returns all ips of scvf as they appear in scv loop More...
 
const MathVector< dim > * scvf_local_ips () const
 returns all ips of scvf as they appear in scv loop More...
 
void update (GridObject *elem, const MathVector< worldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
 update data for given element More...
 
void update_boundary_faces (GridObject *elem, const MathVector< worldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
 update boundary data for given element More...
 
void update_local_data ()
 update local data More...
 

Public Attributes

MathVector< worldDimm_MidPoint
 
MathVector< worldDimvSCVF_PosOnEdge [numSCVF]
 new members for adaptions More...
 

Static Public Attributes

static const int dim = ref_elem_type::dim
 dimension of reference element More...
 
static const size_t nip = 1
 number of integration points More...
 
static const size_t nsh = local_shape_fct_set_type::nsh
 number of shape functions More...
 
static const size_t numSCV
 number of SubControlVolumes More...
 
static const size_t numSCVF
 number of SubControlVolumeFaces More...
 
static const int order = 1
 order More...
 
static const bool staticLocalData = true
 flag indicating if local data may change More...
 
static const bool usesHangingNodes = false
 Hanging node flag: this Geometry does not support hanging nodes. More...
 
static const int worldDim = TWorldDim
 dimension of world More...
 

Protected Attributes

ReferenceMapping< ref_elem_type, worldDimm_mapping
 Reference Mapping. More...
 
std::map< int, std::vector< BF > > m_mapVectorBF
 
TElem * m_pElem
 pointer to current element More...
 
const ref_elem_typem_rRefElem
 Reference Element. More...
 
const local_shape_fct_set_typem_rTrialSpace
 Shape function set. More...
 
std::vector< BFm_vEmptyVectorBF
 
MathVector< worldDimm_vGlobSCVF_IP [numSCVF]
 
MathVector< dimm_vLocSCVF_IP [numSCVF]
 
SCV m_vSCV [numSCV]
 SubControlVolumes. More...
 
SCVF m_vSCVF [numSCVF]
 SubControlVolumeFaces. More...
 
MathVector< worldDimm_vvGloMid [dim+1][maxMid]
 
MathVector< dimm_vvLocMid [dim+1][maxMid]
 local and global geom object midpoints for each dimension More...
 

Static Protected Attributes

static const int maxMid = numSCVF + 1
 max number of geom objects in all dimensions More...
 

Detailed Description

template<typename TElem, int TWorldDim>
class ug::FV1IBGeometry< TElem, TWorldDim >

Geometry and shape functions for 1st order Vertex-Centered Finite Volume.

Template Parameters
TElemElement type
TWorldDim(physical) world dimension

Member Typedef Documentation

◆ elem_type

template<typename TElem , int TWorldDim>
typedef TElem ug::FV1IBGeometry< TElem, TWorldDim >::elem_type

type of element

◆ local_shape_fct_set_type

template<typename TElem , int TWorldDim>
typedef LagrangeP1<ref_elem_type> ug::FV1IBGeometry< TElem, TWorldDim >::local_shape_fct_set_type

type of Shape function used

◆ ref_elem_type

template<typename TElem , int TWorldDim>
typedef reference_element_traits<TElem>::reference_element_type ug::FV1IBGeometry< TElem, TWorldDim >::ref_elem_type

type of reference element

◆ scv_type

template<typename TElem , int TWorldDim>
typedef traits::scv_type ug::FV1IBGeometry< TElem, TWorldDim >::scv_type

type of SubControlVolume

◆ traits

template<typename TElem , int TWorldDim>
typedef fv1_traits<ref_elem_type, TWorldDim> ug::FV1IBGeometry< TElem, TWorldDim >::traits

used traits

Constructor & Destructor Documentation

◆ FV1IBGeometry()

template<typename TElem , int TWorldDim>
FV1IBGeometry::FV1IBGeometry ( )

construct object and initialize local values and sizes

References ug::FV1IBGeometry< TElem, TWorldDim >::update_local_data().

Member Function Documentation

◆ adapt()

template<typename TElem , int TWorldDim>
void FV1IBGeometry::adapt ( GridObject elem,
const MathVector< TWorldDim > *  vCornerCoords,
const ISubsetHandler ish = NULL 
)

adapt integration points for elements cut by the inner boundary

References UG_LOG, ug::VecDot(), and ug::VecSubtract().

◆ adapt_integration_points()

template<typename TElem , int TWorldDim>
void FV1IBGeometry::adapt_integration_points ( GridObject elem,
const MathVector< TWorldDim > *  vCornerCoords,
const ISubsetHandler ish = NULL 
)

adapt integration points for elements cut by the inner boundary

References UG_THROW.

◆ adapt_normals()

template<typename TElem , int TWorldDim>
void FV1IBGeometry::adapt_normals ( GridObject elem,
const MathVector< TWorldDim > *  vCornerCoords,
const ISubsetHandler ish = NULL 
)

adapt normals for elements cut by the inner boundary

BEACHTE: normal on scvf (points direction "from"->"to")

References UG_LOG, and UG_THROW.

◆ add_boundary_subset()

template<typename TElem , int TWorldDim>
void ug::FV1IBGeometry< TElem, TWorldDim >::add_boundary_subset ( int  subsetIndex)
inline

add subset that is interpreted as boundary subset.

References ug::FV1IBGeometry< TElem, TWorldDim >::m_mapVectorBF.

◆ bf() [1/2]

template<typename TElem , int TWorldDim>
const std::vector<BF>& ug::FV1IBGeometry< TElem, TWorldDim >::bf ( int  si) const
inline

returns reference to vector of boundary faces for subsetIndex

References ug::FV1IBGeometry< TElem, TWorldDim >::m_mapVectorBF, and ug::FV1IBGeometry< TElem, TWorldDim >::m_vEmptyVectorBF.

◆ bf() [2/2]

template<typename TElem , int TWorldDim>
const BF& ug::FV1IBGeometry< TElem, TWorldDim >::bf ( int  si,
size_t  i 
) const
inline

◆ clear_boundary_subsets()

template<typename TElem , int TWorldDim>
void ug::FV1IBGeometry< TElem, TWorldDim >::clear_boundary_subsets ( )
inline

reset all boundary subsets

References ug::FV1IBGeometry< TElem, TWorldDim >::m_mapVectorBF.

◆ corners()

template<typename TElem , int TWorldDim>
const MathVector<worldDim>* ug::FV1IBGeometry< TElem, TWorldDim >::corners ( ) const
inline

get vector of the global coordinates of corners for current element

References ug::FV1IBGeometry< TElem, TWorldDim >::m_vvGloMid.

◆ num_bf() [1/2]

template<typename TElem , int TWorldDim>
size_t ug::FV1IBGeometry< TElem, TWorldDim >::num_bf ( ) const
inline

◆ num_bf() [2/2]

template<typename TElem , int TWorldDim>
size_t ug::FV1IBGeometry< TElem, TWorldDim >::num_bf ( int  si) const
inline

number of boundary faces on subset 'subsetIndex'

References ug::FV1IBGeometry< TElem, TWorldDim >::m_mapVectorBF.

◆ num_boundary_subsets()

template<typename TElem , int TWorldDim>
size_t ug::FV1IBGeometry< TElem, TWorldDim >::num_boundary_subsets ( )
inline

number of registered boundary subsets

References ug::FV1IBGeometry< TElem, TWorldDim >::m_mapVectorBF.

◆ num_scv()

template<typename TElem , int TWorldDim>
size_t ug::FV1IBGeometry< TElem, TWorldDim >::num_scv ( ) const
inline

◆ num_scv_ips()

template<typename TElem , int TWorldDim>
size_t ug::FV1IBGeometry< TElem, TWorldDim >::num_scv_ips ( ) const
inline

returns number of all scv ips

References ug::FV1IBGeometry< TElem, TWorldDim >::numSCV.

◆ num_scvf()

template<typename TElem , int TWorldDim>
size_t ug::FV1IBGeometry< TElem, TWorldDim >::num_scvf ( ) const
inline

number of SubControlVolumeFaces

References ug::FV1IBGeometry< TElem, TWorldDim >::numSCVF.

Referenced by ug::FV1IBGeometry< TElem, TWorldDim >::scvf().

◆ num_scvf_ips()

template<typename TElem , int TWorldDim>
size_t ug::FV1IBGeometry< TElem, TWorldDim >::num_scvf_ips ( ) const
inline

returns number of all scvf ips

References ug::FV1IBGeometry< TElem, TWorldDim >::numSCVF.

◆ num_sh()

template<typename TElem , int TWorldDim>
size_t ug::FV1IBGeometry< TElem, TWorldDim >::num_sh ( ) const
inline

number of shape functions

References ug::FV1IBGeometry< TElem, TWorldDim >::nsh.

◆ remove_boundary_subset()

template<typename TElem , int TWorldDim>
void ug::FV1IBGeometry< TElem, TWorldDim >::remove_boundary_subset ( int  subsetIndex)
inline

removes subset that is interpreted as boundary subset.

References ug::FV1IBGeometry< TElem, TWorldDim >::m_mapVectorBF.

◆ scv()

template<typename TElem , int TWorldDim>
const SCV& ug::FV1IBGeometry< TElem, TWorldDim >::scv ( size_t  i) const
inline

◆ scv_global_ips()

template<typename TElem , int TWorldDim>
const MathVector<worldDim>* ug::FV1IBGeometry< TElem, TWorldDim >::scv_global_ips ( ) const
inline

returns all ips of scv as they appear in scv loop

References ug::FV1IBGeometry< TElem, TWorldDim >::m_vvGloMid.

◆ scv_local_ips()

template<typename TElem , int TWorldDim>
const MathVector<dim>* ug::FV1IBGeometry< TElem, TWorldDim >::scv_local_ips ( ) const
inline

returns all ips of scv as they appear in scv loop

References ug::FV1IBGeometry< TElem, TWorldDim >::m_vvLocMid.

◆ scvf()

template<typename TElem , int TWorldDim>
const SCVF& ug::FV1IBGeometry< TElem, TWorldDim >::scvf ( size_t  i) const
inline

◆ scvf_global_ips()

template<typename TElem , int TWorldDim>
const MathVector<worldDim>* ug::FV1IBGeometry< TElem, TWorldDim >::scvf_global_ips ( ) const
inline

returns all ips of scvf as they appear in scv loop

References ug::FV1IBGeometry< TElem, TWorldDim >::m_vGlobSCVF_IP.

◆ scvf_local_ips()

template<typename TElem , int TWorldDim>
const MathVector<dim>* ug::FV1IBGeometry< TElem, TWorldDim >::scvf_local_ips ( ) const
inline

returns all ips of scvf as they appear in scv loop

References ug::FV1IBGeometry< TElem, TWorldDim >::m_vLocSCVF_IP.

◆ update()

template<typename TElem , int TWorldDim>
void FV1IBGeometry::update ( GridObject elem,
const MathVector< worldDim > *  vCornerCoords,
const ISubsetHandler ish = NULL 
)

update data for given element

References ug::AveragePositions(), ug::MatVecMult(), and UG_ASSERT.

◆ update_boundary_faces()

◆ update_local_data()

template<typename TElem , int TWorldDim>
void FV1IBGeometry::update_local_data ( )

Member Data Documentation

◆ dim

template<typename TElem , int TWorldDim>
const int ug::FV1IBGeometry< TElem, TWorldDim >::dim = ref_elem_type::dim
static

dimension of reference element

◆ m_mapping

template<typename TElem , int TWorldDim>
ReferenceMapping<ref_elem_type, worldDim> ug::FV1IBGeometry< TElem, TWorldDim >::m_mapping
protected

Reference Mapping.

◆ m_mapVectorBF

◆ m_MidPoint

template<typename TElem , int TWorldDim>
MathVector<worldDim> ug::FV1IBGeometry< TElem, TWorldDim >::m_MidPoint

◆ m_pElem

template<typename TElem , int TWorldDim>
TElem* ug::FV1IBGeometry< TElem, TWorldDim >::m_pElem
protected

pointer to current element

◆ m_rRefElem

template<typename TElem , int TWorldDim>
const ref_elem_type& ug::FV1IBGeometry< TElem, TWorldDim >::m_rRefElem
protected

Reference Element.

◆ m_rTrialSpace

template<typename TElem , int TWorldDim>
const local_shape_fct_set_type& ug::FV1IBGeometry< TElem, TWorldDim >::m_rTrialSpace
protected

Shape function set.

◆ m_vEmptyVectorBF

template<typename TElem , int TWorldDim>
std::vector<BF> ug::FV1IBGeometry< TElem, TWorldDim >::m_vEmptyVectorBF
protected

◆ m_vGlobSCVF_IP

template<typename TElem , int TWorldDim>
MathVector<worldDim> ug::FV1IBGeometry< TElem, TWorldDim >::m_vGlobSCVF_IP[numSCVF]
protected

◆ m_vLocSCVF_IP

template<typename TElem , int TWorldDim>
MathVector<dim> ug::FV1IBGeometry< TElem, TWorldDim >::m_vLocSCVF_IP[numSCVF]
protected

◆ m_vSCV

template<typename TElem , int TWorldDim>
SCV ug::FV1IBGeometry< TElem, TWorldDim >::m_vSCV[numSCV]
protected

SubControlVolumes.

Referenced by ug::FV1IBGeometry< TElem, TWorldDim >::scv().

◆ m_vSCVF

template<typename TElem , int TWorldDim>
SCVF ug::FV1IBGeometry< TElem, TWorldDim >::m_vSCVF[numSCVF]
protected

SubControlVolumeFaces.

Referenced by ug::FV1IBGeometry< TElem, TWorldDim >::scvf().

◆ m_vvGloMid

template<typename TElem , int TWorldDim>
MathVector<worldDim> ug::FV1IBGeometry< TElem, TWorldDim >::m_vvGloMid[dim+1][maxMid]
protected

◆ m_vvLocMid

template<typename TElem , int TWorldDim>
MathVector<dim> ug::FV1IBGeometry< TElem, TWorldDim >::m_vvLocMid[dim+1][maxMid]
protected

local and global geom object midpoints for each dimension

Referenced by ug::FV1IBGeometry< TElem, TWorldDim >::scv_local_ips().

◆ maxMid

template<typename TElem , int TWorldDim>
const int ug::FV1IBGeometry< TElem, TWorldDim >::maxMid = numSCVF + 1
staticprotected

max number of geom objects in all dimensions

◆ nip

template<typename TElem , int TWorldDim>
const size_t ug::FV1IBGeometry< TElem, TWorldDim >::nip = 1
static

◆ nsh

template<typename TElem , int TWorldDim>
const size_t ug::FV1IBGeometry< TElem, TWorldDim >::nsh = local_shape_fct_set_type::nsh
static

◆ numSCV

template<typename TElem , int TWorldDim>
const size_t ug::FV1IBGeometry< TElem, TWorldDim >::numSCV
static
Initial value:
= (ref_elem_type::REFERENCE_OBJECT_ID != ROID_PYRAMID)
? ref_elem_type::numCorners : 8
@ ROID_PYRAMID
Definition: grid_base_objects.h:83

number of SubControlVolumes

Referenced by ug::FV1IBGeometry< TElem, TWorldDim >::num_scv(), and ug::FV1IBGeometry< TElem, TWorldDim >::num_scv_ips().

◆ numSCVF

template<typename TElem , int TWorldDim>
const size_t ug::FV1IBGeometry< TElem, TWorldDim >::numSCVF
static
Initial value:
= (ref_elem_type::REFERENCE_OBJECT_ID != ROID_PYRAMID)
? ref_elem_type::numEdges : 12

number of SubControlVolumeFaces

Referenced by ug::FV1IBGeometry< TElem, TWorldDim >::num_scvf(), and ug::FV1IBGeometry< TElem, TWorldDim >::num_scvf_ips().

◆ order

template<typename TElem , int TWorldDim>
const int ug::FV1IBGeometry< TElem, TWorldDim >::order = 1
static

order

◆ staticLocalData

template<typename TElem , int TWorldDim>
const bool ug::FV1IBGeometry< TElem, TWorldDim >::staticLocalData = true
static

flag indicating if local data may change

◆ usesHangingNodes

template<typename TElem , int TWorldDim>
const bool ug::FV1IBGeometry< TElem, TWorldDim >::usesHangingNodes = false
static

Hanging node flag: this Geometry does not support hanging nodes.

◆ vSCVF_PosOnEdge

template<typename TElem , int TWorldDim>
MathVector<worldDim> ug::FV1IBGeometry< TElem, TWorldDim >::vSCVF_PosOnEdge[numSCVF]

new members for adaptions

◆ worldDim

template<typename TElem , int TWorldDim>
const int ug::FV1IBGeometry< TElem, TWorldDim >::worldDim = TWorldDim
static

dimension of world


The documentation for this class was generated from the following files: