33#ifndef __H__UG__LIB_DISC__SPATIAL_DISC__DISC_HELPER__FINITE_VOLUME_HIGHER_ORDER_GEOMETRY__
34#define __H__UG__LIB_DISC__SPATIAL_DISC__DISC_HELPER__FINITE_VOLUME_HIGHER_ORDER_GEOMETRY__
69template <
int TOrder,
typename TElem,
int TWorldDim,
70 int TQuadOrder = TOrder + 1>
75 static const int p = TOrder;
85 static const int dim = ref_elem_type::dim;
92 static const int order = TOrder;
143 static const size_t nsh = local_shape_fct_set_type::nsh;
157 static const size_t nip = scvf_quad_rule_type::nip;
160 static const size_t numCo = traits::NumCornersOfSCVF;
169 inline size_t to()
const {
return To;}
239 friend class FVGeometry<TOrder, TElem, TWorldDim, TQuadOrder>;
275 static const size_t nip = scv_quad_rule_type::nip;
278 static const size_t numCo = traits::NumCornersOfSCV;
351 friend class FVGeometry<TOrder, TElem, TWorldDim, TQuadOrder>;
380 static const size_t nip = scvf_quad_rule_type::nip;
383 static const size_t numCo = traits::NumCornersOfSCVF;
462 friend class FVGeometry<TOrder, TElem, TWorldDim, TQuadOrder>;
501 size_t orderQuad = TQuadOrder);
566 size_t vDoFID[ref_elem_type::numCorners];
599 typename std::map<int, std::vector<BF> >::const_iterator it;
602 num += (*it).second.size();
609 typename std::map<int, std::vector<BF> >::const_iterator it;
612 else return (*it).second.size();
616 inline const BF&
bf(
int si,
size_t i)
const
618 typename std::map<int, std::vector<BF> >::const_iterator it;
621 return (*it).second[i];
616 inline const BF&
bf(
int si,
size_t i)
const {
…}
625 inline const std::vector<BF>&
bf(
int si)
const
627 typename std::map<int, std::vector<BF> >::const_iterator it;
625 inline const std::vector<BF>&
bf(
int si)
const {
…}
677template <
int TWorldDim,
int TDim = TWorldDim>
685 static const int dim = TDim;
709 static const size_t numCo = traits::NumCornersOfSCVF;
718 inline size_t to()
const {
return To;}
816 std::vector<MathMatrix<worldDim,dim> >
vJtInv;
826 static const size_t numCo = traits::NumCornersOfSCV;
920 std::vector<MathMatrix<worldDim,dim> >
vJtInv;
930 static const size_t numCo = traits::NumCornersOfSCVF;
1035 std::vector<MathMatrix<worldDim,dim> >
vJtInv;
1060 const LFEID& lfeID,
size_t orderQuad,
1116 static const int maxMid = traits::maxNumSCVF +1;
1155 typename std::map<int, std::vector<BF> >::const_iterator it;
1158 num += (*it).second.size();
1165 typename std::map<int, std::vector<BF> >::const_iterator it;
1168 else return (*it).second.size();
1172 inline const BF&
bf(
int si,
size_t i)
const
1175 typename std::map<int, std::vector<BF> >::const_iterator it;
1178 return (*it).second[i];
1172 inline const BF&
bf(
int si,
size_t i)
const {
…}
1182 inline const std::vector<BF>&
bf(
int si)
const
1184 typename std::map<int, std::vector<BF> >::const_iterator it;
1187 return (*it).second;
1182 inline const std::vector<BF>&
bf(
int si)
const {
…}
boundary face
Definition fvho_geom.h:927
const MathVector< dim > & local_corner(size_t co) const
return local corner number i
Definition fvho_geom.h:1000
const number * shape_vector(size_t ip) const
vector of shape functions in ip point
Definition fvho_geom.h:975
const MathVector< dim > & local_grad(size_t ip, size_t sh) const
value of local gradient of shape function i in integration point
Definition fvho_geom.h:979
std::vector< number > vDetJ
Definition fvho_geom.h:1036
size_t nodeId
Definition fvho_geom.h:1012
BF()
Definition fvho_geom.h:933
size_t num_corners() const
number of corners, that bound the scvf
Definition fvho_geom.h:997
number shape(size_t ip, size_t sh) const
value of shape function i in integration point
Definition fvho_geom.h:971
std::vector< std::vector< MathVector< dim > > > vvLocalGrad
Definition fvho_geom.h:1033
size_t num_sh() const
number of shape functions
Definition fvho_geom.h:968
const MathVector< dim > & local_ip(size_t ip) const
local integration point of scvf
Definition fvho_geom.h:952
size_t node_id() const
index of SubControlVolume of the bf
Definition fvho_geom.h:936
const MathMatrix< worldDim, dim > & JTInv(size_t ip) const
Transposed Inverse of Jacobian in integration point.
Definition fvho_geom.h:960
size_t nip
Definition fvho_geom.h:1023
const MathVector< worldDim > & global_ip(size_t ip) const
global integration point of scvf
Definition fvho_geom.h:956
std::vector< MathMatrix< worldDim, dim > > vJtInv
Definition fvho_geom.h:1035
std::vector< MathVector< dim > > vLocalIP
Definition fvho_geom.h:1024
std::vector< MathVector< worldDim > > vGlobalIP
Definition fvho_geom.h:1025
const MathVector< worldDim > & global_grad(size_t ip, size_t sh) const
value of global gradient of shape function i in integration point
Definition fvho_geom.h:988
number Vol
Definition fvho_geom.h:1028
static const size_t numCo
Number of corners of bf.
Definition fvho_geom.h:930
MidID vMidID[numCo]
Definition fvho_geom.h:1020
std::vector< std::vector< MathVector< worldDim > > > vvGlobalGrad
Definition fvho_geom.h:1034
const MathVector< worldDim > * global_grad_vector(size_t ip) const
vector of global gradients in ip point
Definition fvho_geom.h:993
number weight(size_t ip) const
integration weight
Definition fvho_geom.h:948
number detJ(size_t ip) const
Determinant of Jacobian in integration point.
Definition fvho_geom.h:964
const number * vWeight
Definition fvho_geom.h:1026
size_t num_ip() const
number of integration points on scvf
Definition fvho_geom.h:945
number volume() const
volume of bf
Definition fvho_geom.h:942
MathVector< dim > vLocPos[numCo]
Definition fvho_geom.h:1018
std::vector< std::vector< number > > vvShape
Definition fvho_geom.h:1032
MathVector< worldDim > Normal
Definition fvho_geom.h:1027
MathVector< worldDim > vGloPos[numCo]
Definition fvho_geom.h:1019
const MathVector< worldDim > & global_corner(size_t co) const
return glbal corner number i
Definition fvho_geom.h:1004
const MathVector< worldDim > & normal() const
outer normal on bf. Norm is equal to area
Definition fvho_geom.h:939
const MathVector< dim > * local_grad_vector(size_t ip) const
vector of local gradients in ip point
Definition fvho_geom.h:984
size_t nsh
Definition fvho_geom.h:1031
Sub-Control Volume Face structure.
Definition fvho_geom.h:706
size_t From
Definition fvho_geom.h:792
const MathVector< dim > * local_grad_vector(size_t ip) const
vector of local gradients in ip point
Definition fvho_geom.h:763
std::vector< number > vDetJ
Definition fvho_geom.h:817
const MathVector< worldDim > & global_grad(size_t ip, size_t sh) const
value of global gradient of shape function i in integration point
Definition fvho_geom.h:767
const MathVector< dim > & local_corner(size_t co) const
return local corner number i
Definition fvho_geom.h:779
MathVector< worldDim > vGloPos[numCo]
Definition fvho_geom.h:802
size_t to() const
index of SubControlVolume on one side of the scvf
Definition fvho_geom.h:718
size_t num_sh() const
number of shape functions
Definition fvho_geom.h:747
const MathVector< worldDim > & global_ip(size_t ip) const
global integration point of scvf
Definition fvho_geom.h:732
static const size_t numCo
Number of corners of scvf.
Definition fvho_geom.h:709
MidID vMidID[numCo]
Definition fvho_geom.h:803
size_t num_ip() const
number of integration points on scvf
Definition fvho_geom.h:721
number shape(size_t ip, size_t sh) const
value of shape function i in integration point
Definition fvho_geom.h:750
const MathVector< dim > & local_ip(size_t ip) const
local integration point of scvf
Definition fvho_geom.h:728
const MathVector< dim > & local_grad(size_t ip, size_t sh) const
value of local gradient of shape function i in integration point
Definition fvho_geom.h:758
std::vector< MathVector< dim > > vLocalIP
Definition fvho_geom.h:807
size_t To
Definition fvho_geom.h:792
std::vector< std::vector< MathVector< worldDim > > > vvGlobalGrad
Definition fvho_geom.h:815
std::vector< std::vector< MathVector< dim > > > vvLocalGrad
Definition fvho_geom.h:814
size_t from() const
index of SubControlVolume on one side of the scvf
Definition fvho_geom.h:715
size_t nip
Definition fvho_geom.h:806
number weight(size_t ip) const
integration weight
Definition fvho_geom.h:724
MathVector< worldDim > Normal
Definition fvho_geom.h:795
const MathVector< worldDim > & normal() const
normal on scvf (points direction "from"->"to"), normalized
Definition fvho_geom.h:736
const MathVector< worldDim > & global_corner(size_t co) const
return glbal corner number i
Definition fvho_geom.h:783
const MathVector< worldDim > * global_grad_vector(size_t ip) const
vector of gloabl gradients in ip point
Definition fvho_geom.h:772
std::vector< MathMatrix< worldDim, dim > > vJtInv
Definition fvho_geom.h:816
SCVF()
Definition fvho_geom.h:712
const number * shape_vector(size_t ip) const
vector of shape functions in ip point
Definition fvho_geom.h:754
std::vector< MathVector< worldDim > > vGlobalIP
Definition fvho_geom.h:808
size_t num_corners() const
number of corners, that bound the scvf
Definition fvho_geom.h:776
MathVector< dim > vLocPos[numCo]
Definition fvho_geom.h:801
const number * vWeight
Definition fvho_geom.h:809
size_t nsh
Definition fvho_geom.h:812
number detJ(size_t ip) const
Determinant of Jacobian in integration point.
Definition fvho_geom.h:743
const MathMatrix< worldDim, dim > & JTInv(size_t ip) const
Transposed Inverse of Jacobian in integration point.
Definition fvho_geom.h:739
std::vector< std::vector< number > > vvShape
Definition fvho_geom.h:813
std::vector< number > vDetJMap
Definition fvho_geom.h:818
sub control volume structure
Definition fvho_geom.h:823
std::vector< number > vDetJMap
Definition fvho_geom.h:922
const MathMatrix< worldDim, dim > & JTInv(size_t ip) const
Transposed Inverse of Jacobian in integration point.
Definition fvho_geom.h:861
MathVector< dim > vLocPos[numCo]
Definition fvho_geom.h:905
SCV()
Definition fvho_geom.h:829
size_t nip
Definition fvho_geom.h:910
const MathVector< dim > & local_ip(size_t ip) const
local integration point of scv
Definition fvho_geom.h:853
std::vector< std::vector< number > > vvShape
Definition fvho_geom.h:917
const MathVector< worldDim > & global_corner(size_t co) const
return glbal corner number i
Definition fvho_geom.h:839
number shape(size_t ip, size_t sh) const
value of shape function i in integration point
Definition fvho_geom.h:872
const MathVector< dim > * local_grad_vector(size_t ip) const
vector of local gradients in ip point
Definition fvho_geom.h:885
const MathVector< worldDim > & global_grad(size_t ip, size_t sh) const
value of global gradient of shape function i in integration point
Definition fvho_geom.h:889
size_t num_sh() const
number of shape functions
Definition fvho_geom.h:869
const MathVector< dim > & local_corner(size_t co) const
return local corner number i
Definition fvho_geom.h:835
const MathVector< dim > & local_grad(size_t ip, size_t sh) const
value of local gradient of shape function i in integration point
Definition fvho_geom.h:880
const number * vWeight
Definition fvho_geom.h:913
std::vector< MathMatrix< worldDim, dim > > vJtInv
Definition fvho_geom.h:920
size_t num_corners() const
number of corners, that bound the scvf
Definition fvho_geom.h:832
MidID vMidID[numCo]
Definition fvho_geom.h:907
static const size_t numCo
Number of corners of scvf.
Definition fvho_geom.h:826
const number * shape_vector(size_t ip) const
vector of shape functions in ip point
Definition fvho_geom.h:876
std::vector< MathVector< worldDim > > vGlobalIP
Definition fvho_geom.h:912
const MathVector< worldDim > * global_grad_vector(size_t ip) const
vector of global gradients in ip point
Definition fvho_geom.h:894
std::vector< number > vDetJ
Definition fvho_geom.h:921
MathVector< worldDim > vGloPos[numCo]
Definition fvho_geom.h:906
std::vector< MathVector< dim > > vLocalIP
Definition fvho_geom.h:911
const MathVector< worldDim > & global_ip(size_t ip) const
global integration point
Definition fvho_geom.h:857
number weight(size_t ip) const
weigth of integration point
Definition fvho_geom.h:849
std::vector< std::vector< MathVector< dim > > > vvLocalGrad
Definition fvho_geom.h:918
size_t num_ip() const
number of integration points
Definition fvho_geom.h:846
size_t nodeId
Definition fvho_geom.h:902
number detJ(size_t ip) const
Determinant of Jacobian in integration point.
Definition fvho_geom.h:865
size_t nsh
Definition fvho_geom.h:916
size_t node_id() const
node id that this scv is associated to
Definition fvho_geom.h:843
std::vector< std::vector< MathVector< worldDim > > > vvGlobalGrad
Definition fvho_geom.h:919
Geometry and shape functions for any order Vertex-Centered Finite Volume.
Definition fvho_geom.h:679
DimFVGeometry()
construct object and initialize local values and sizes
Definition fvho_geom.cpp:1396
const MathVector< dim > * scv_local_ips() const
returns all ips of scv as they appear in scv loop
Definition fvho_geom.h:1099
void reset_curr_elem()
Definition fvho_geom.h:1190
std::vector< MathVector< dim > > m_vLocSCV_IP
Definition fvho_geom.h:1112
const BF & bf(int si, size_t i) const
returns the boundary face i for subsetIndex
Definition fvho_geom.h:1172
traits::scvf_type scvf_type
type of SubControlVolumeFace
Definition fvho_geom.h:692
void clear_boundary_subsets()
reset all boundary subsets
Definition fvho_geom.h:1147
void update_local(ReferenceObjectID roid, const LFEID &lfeID, size_t orderQuad)
update local data
Definition fvho_geom.cpp:1402
size_t num_scvf_ips() const
returns number of all scvf ips
Definition fvho_geom.h:1087
size_t num_scv_ips() const
returns number of all scv ips
Definition fvho_geom.h:1096
std::vector< SCVF > m_vSCVF
SubControlVolumeFaces (size: numSCVF)
Definition fvho_geom.h:1222
size_t m_numSCVFIP
number of scvf ip
Definition fvho_geom.h:1228
size_t m_numSCV
current number of scv
Definition fvho_geom.h:1232
size_t m_numSubElem
current number of subelements
Definition fvho_geom.h:1210
size_t num_scvf() const
number of SubControlVolumeFaces
Definition fvho_geom.h:1069
traits::scv_type scv_type
type of SubControlVolume
Definition fvho_geom.h:695
void add_boundary_subset(int subsetIndex)
add subset that is interpreted as boundary subset.
Definition fvho_geom.h:1141
static const bool usesHangingNodes
Hanging node flag: this Geometry does not support hanging nodes.
Definition fvho_geom.h:698
size_t num_scv() const
number of SubControlVolumes
Definition fvho_geom.h:1076
const SCVF & scvf(size_t i) const
const access to SubControlVolumeFace number i
Definition fvho_geom.h:1072
std::vector< BF > m_vEmptyVectorBF
Definition fvho_geom.h:1194
void update_boundary_faces(GridObject *pElem, const MathVector< worldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
update boundary data for given element
Definition fvho_geom.cpp:1731
void remove_boundary_subset(int subsetIndex)
removes subset that is interpreted as boundary subset.
Definition fvho_geom.h:1144
size_t num_bf() const
number of all boundary faces
Definition fvho_geom.h:1153
LFEID m_lfeID
current trial space
Definition fvho_geom.h:1207
static const int worldDim
dimension of world
Definition fvho_geom.h:688
std::vector< MathVector< worldDim > > m_vGlobSCVF_IP
Definition fvho_geom.h:1107
const MathVector< worldDim > * scv_global_ips() const
returns all ips of scv as they appear in scv loop
Definition fvho_geom.h:1102
const MathVector< worldDim > * scvf_global_ips() const
returns all ips of scvf as they appear in scv loop
Definition fvho_geom.h:1093
int m_quadOrderSCVF
quadrature order
Definition fvho_geom.h:1225
std::vector< SubElement > m_vSubElem
subelements (size: numSubElem)
Definition fvho_geom.h:1137
const std::vector< BF > & bf(int si) const
returns reference to vector of boundary faces for subsetIndex
Definition fvho_geom.h:1182
int m_orderShape
current order
Definition fvho_geom.h:1204
static const bool staticLocalData
flag indicating if local data may change
Definition fvho_geom.h:701
void update(GridObject *pElem, const MathVector< worldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
update data for given element
Definition fvho_geom.h:1052
std::vector< MathVector< dim > > m_vLocSCVF_IP
Definition fvho_geom.h:1108
const SCV & scv(size_t i) const
const access to SubControlVolume number i
Definition fvho_geom.h:1079
size_t m_nsh
number of shape functions
Definition fvho_geom.h:1213
size_t num_bf(int si) const
number of boundary faces on subset 'subsetIndex'
Definition fvho_geom.h:1163
size_t num_sh() const
number of shape functions
Definition fvho_geom.h:1083
std::vector< MathVector< worldDim > > m_vGlobSCV_IP
Definition fvho_geom.h:1111
std::map< int, std::vector< BF > > m_mapVectorBF
Definition fvho_geom.h:1193
size_t m_numSCVPerSubElem
number of SCV per SubElement
Definition fvho_geom.h:1235
size_t m_numSCVIP
number of scv ip
Definition fvho_geom.h:1244
size_t num_boundary_subsets()
number of registered boundary subsets
Definition fvho_geom.h:1150
size_t m_numSCVFPerSubElem
current number of SCVF per SubElement
Definition fvho_geom.h:1219
void update_local(ReferenceObjectID roid, const LFEID &lfeID)
Definition fvho_geom.h:1046
size_t m_numSCVF
current number of scvf
Definition fvho_geom.h:1216
static const int dim
dimension of reference element
Definition fvho_geom.h:685
std::vector< SCV > m_vSCV
SubControlVolumes (size: numSCV)
Definition fvho_geom.h:1238
ReferenceObjectID m_roid
current reference object id
Definition fvho_geom.h:1201
int m_quadOrderSCV
current quadrature order scv
Definition fvho_geom.h:1241
fv1_dim_traits< TDim, TWorldDim > traits
traits used
Definition fvho_geom.h:682
GridObject * m_pElem
pointer to current element
Definition fvho_geom.h:1198
static const int maxMid
Definition fvho_geom.h:1116
const MathVector< dim > * scvf_local_ips() const
returns all ips of scvf as they appear in scv loop
Definition fvho_geom.h:1090
boundary face
Definition fvho_geom.h:377
const MathVector< dim > & local_ip(size_t ip) const
local integration point of scvf
Definition fvho_geom.h:405
MathMatrix< worldDim, dim > vJtInv[nip]
Definition fvho_geom.h:486
const MathVector< dim > & local_grad(size_t ip, size_t sh) const
value of local gradient of shape function i in integration point
Definition fvho_geom.h:432
MathVector< worldDim > Normal
Definition fvho_geom.h:479
MathVector< dim > vvLocalGrad[nip][nsh]
Definition fvho_geom.h:484
const MathMatrix< worldDim, dim > & JTInv(size_t ip) const
Transposed Inverse of Jacobian in integration point.
Definition fvho_geom.h:413
number vDetJ[nip]
Definition fvho_geom.h:487
MathVector< worldDim > vGlobalIP[nip]
Definition fvho_geom.h:477
number weight(size_t ip) const
integration weight
Definition fvho_geom.h:401
size_t nodeId
Definition fvho_geom.h:465
size_t num_corners() const
number of corners, that bound the scvf
Definition fvho_geom.h:450
const MathVector< dim > & local_corner(size_t co) const
return local corner number i
Definition fvho_geom.h:453
static const size_t nip
number of integration points
Definition fvho_geom.h:380
const MathVector< worldDim > * global_grad_vector(size_t ip) const
vector of global gradients in ip point
Definition fvho_geom.h:446
size_t num_sh() const
number of shape functions
Definition fvho_geom.h:421
MathVector< dim > vLocPos[numCo]
Definition fvho_geom.h:471
MathVector< dim > vLocalIP[nip]
Definition fvho_geom.h:476
size_t node_id() const
index of SubControlVolume of the bf
Definition fvho_geom.h:389
const number * vWeight
Definition fvho_geom.h:478
MathVector< worldDim > vvGlobalGrad[nip][nsh]
Definition fvho_geom.h:485
const number * shape_vector(size_t ip) const
vector of shape functions in ip point
Definition fvho_geom.h:428
MidID vMidID[numCo]
Definition fvho_geom.h:473
const MathVector< dim > * local_grad_vector(size_t ip) const
vector of local gradients in ip point
Definition fvho_geom.h:437
const MathVector< worldDim > & global_corner(size_t co) const
return glbal corner number i
Definition fvho_geom.h:457
const MathVector< worldDim > & global_ip(size_t ip) const
global integration point of scvf
Definition fvho_geom.h:409
number vvShape[nip][nsh]
Definition fvho_geom.h:483
BF()
Definition fvho_geom.h:386
number detJ(size_t ip) const
Determinant of Jacobian in integration point.
Definition fvho_geom.h:417
number volume() const
volume of bf
Definition fvho_geom.h:395
number shape(size_t ip, size_t sh) const
value of shape function i in integration point
Definition fvho_geom.h:424
const MathVector< worldDim > & normal() const
outer normal on bf. Norm is equal to area
Definition fvho_geom.h:392
const MathVector< worldDim > & global_grad(size_t ip, size_t sh) const
value of global gradient of shape function i in integration point
Definition fvho_geom.h:441
size_t num_ip() const
number of integration points on scvf
Definition fvho_geom.h:398
static const size_t numCo
Number of corners of bf.
Definition fvho_geom.h:383
MathVector< worldDim > vGloPos[numCo]
Definition fvho_geom.h:472
number Vol
Definition fvho_geom.h:480
Sub-Control Volume Face structure.
Definition fvho_geom.h:154
const MathVector< dim > * local_grad_vector(size_t ip) const
vector of local gradients in ip point
Definition fvho_geom.h:214
number detJ(size_t ip) const
Determinant of Jacobian in integration point.
Definition fvho_geom.h:194
size_t num_corners() const
number of corners, that bound the scvf
Definition fvho_geom.h:227
MathVector< worldDim > vGlobalIP[nip]
Definition fvho_geom.h:258
size_t To
Definition fvho_geom.h:243
MidID vMidID[numCo]
Definition fvho_geom.h:254
const MathMatrix< worldDim, dim > & JTInv(size_t ip) const
Transposed Inverse of Jacobian in integration point.
Definition fvho_geom.h:190
const number * vWeight
Definition fvho_geom.h:259
size_t to() const
index of SubControlVolume on one side of the scvf
Definition fvho_geom.h:169
size_t num_ip() const
number of integration points on scvf
Definition fvho_geom.h:172
const MathVector< dim > & local_corner(size_t co) const
return local corner number i
Definition fvho_geom.h:230
const MathVector< dim > & local_grad(size_t ip, size_t sh) const
value of local gradient of shape function i in integration point
Definition fvho_geom.h:209
number vDetJ[nip]
Definition fvho_geom.h:266
size_t from() const
index of SubControlVolume on one side of the scvf
Definition fvho_geom.h:166
const MathVector< worldDim > & global_ip(size_t ip) const
global integration point of scvf
Definition fvho_geom.h:183
number weight(size_t ip) const
integration weight
Definition fvho_geom.h:175
static const size_t numCo
Number of corners of scvf.
Definition fvho_geom.h:160
size_t From
Definition fvho_geom.h:243
MathVector< dim > vLocPos[numCo]
Definition fvho_geom.h:252
const number * shape_vector(size_t ip) const
vector of shape functions in ip point
Definition fvho_geom.h:205
MathMatrix< worldDim, dim > vJtInv[nip]
Definition fvho_geom.h:265
const MathVector< worldDim > & global_corner(size_t co) const
return global corner number i
Definition fvho_geom.h:234
MathVector< worldDim > vGloPos[numCo]
Definition fvho_geom.h:253
MathVector< worldDim > vvGlobalGrad[nip][nsh]
Definition fvho_geom.h:264
const MathVector< worldDim > & global_grad(size_t ip, size_t sh) const
value of global gradient of shape function i in integration point
Definition fvho_geom.h:218
const MathVector< worldDim > * global_grad_vector(size_t ip) const
vector of global gradients in ip point
Definition fvho_geom.h:223
number vvShape[nip][nsh]
Definition fvho_geom.h:262
number shape(size_t ip, size_t sh) const
value of shape function i in integration point
Definition fvho_geom.h:201
number vDetJMap[nip]
Definition fvho_geom.h:267
MathVector< worldDim > Normal
Definition fvho_geom.h:246
SCVF()
Definition fvho_geom.h:163
static const size_t nip
Number of integration points.
Definition fvho_geom.h:157
const MathVector< worldDim > & normal() const
normal on scvf (points direction "from"->"to"), normalized
Definition fvho_geom.h:187
MathVector< dim > vLocalIP[nip]
Definition fvho_geom.h:257
MathVector< dim > vvLocalGrad[nip][nsh]
Definition fvho_geom.h:263
const MathVector< dim > & local_ip(size_t ip) const
local integration point of scvf
Definition fvho_geom.h:179
size_t num_sh() const
number of shape functions
Definition fvho_geom.h:198
sub control volume structure
Definition fvho_geom.h:272
const MathVector< worldDim > & global_corner(size_t co) const
return glbal corner number i
Definition fvho_geom.h:291
static const size_t numCo
Number of corners of scvf.
Definition fvho_geom.h:278
number vDetJMap[nip]
Definition fvho_geom.h:372
MathVector< worldDim > vGloPos[numCo]
Definition fvho_geom.h:363
const MathVector< dim > & local_ip(size_t ip) const
local integration point of scv
Definition fvho_geom.h:305
MathVector< dim > vLocalIP[nip]
Definition fvho_geom.h:357
MathVector< dim > vvLocalGrad[nip][nsh]
Definition fvho_geom.h:368
const MathVector< worldDim > & global_grad(size_t ip, size_t sh) const
value of global gradient of shape function i in integration point
Definition fvho_geom.h:341
static const size_t nip
Number of integration points.
Definition fvho_geom.h:275
const number * vWeight
Definition fvho_geom.h:359
number detJ(size_t ip) const
Determinant of Jacobian in integration point.
Definition fvho_geom.h:317
size_t node_id() const
node id that this scv is associated to
Definition fvho_geom.h:295
MathMatrix< worldDim, dim > vJtInv[nip]
Definition fvho_geom.h:370
const MathVector< dim > * local_grad_vector(size_t ip) const
vector of local gradients in ip point
Definition fvho_geom.h:337
size_t num_corners() const
number of corners, that bound the scvf
Definition fvho_geom.h:284
MathVector< worldDim > vvGlobalGrad[nip][nsh]
Definition fvho_geom.h:369
const MathVector< worldDim > * global_grad_vector(size_t ip) const
vector of global gradients in ip point
Definition fvho_geom.h:346
MathVector< dim > vLocPos[numCo]
Definition fvho_geom.h:362
const MathMatrix< worldDim, dim > & JTInv(size_t ip) const
Transposed Inverse of Jacobian in integration point.
Definition fvho_geom.h:313
const MathVector< worldDim > & global_ip(size_t ip) const
global integration point
Definition fvho_geom.h:309
size_t num_sh() const
number of shape functions
Definition fvho_geom.h:321
const MathVector< dim > & local_grad(size_t ip, size_t sh) const
value of local gradient of shape function i in integration point
Definition fvho_geom.h:332
number shape(size_t ip, size_t sh) const
value of shape function i in integration point
Definition fvho_geom.h:324
const number * shape_vector(size_t ip) const
vector of shape functions in ip point
Definition fvho_geom.h:328
number vvShape[nip][nsh]
Definition fvho_geom.h:367
number vDetJ[nip]
Definition fvho_geom.h:371
MidID vMidID[numCo]
Definition fvho_geom.h:364
SCV()
Definition fvho_geom.h:281
MathVector< worldDim > vGlobalIP[nip]
Definition fvho_geom.h:358
number weight(size_t ip) const
weigth of integration point
Definition fvho_geom.h:301
const MathVector< dim > & local_corner(size_t co) const
return local corner number i
Definition fvho_geom.h:287
size_t num_ip() const
number of integration points
Definition fvho_geom.h:298
size_t nodeId
Definition fvho_geom.h:354
base class for all FVGeometries
Definition fv_geom_base.h:76
Geometry and shape functions for any order Vertex-Centered Finite Volume.
Definition fvho_geom.h:72
void update_local_data()
update local data
Definition fvho_geom.cpp:997
TElem * m_pElem
pointer to current element
Definition fvho_geom.h:641
const scvf_quad_rule_type & m_rSCVFQuadRule
Quad Rule scvf.
Definition fvho_geom.h:662
size_t num_scv() const
number of SubControlVolumes
Definition fvho_geom.h:520
SCV m_vSCV[numSCV]
SubControlVolumes.
Definition fvho_geom.h:650
FVGeometry()
construct object and initialize local values and sizes
Definition fvho_geom.cpp:957
static const bool staticLocalData
flag indicating if local data may change
Definition fvho_geom.h:149
MathVector< worldDim > m_vGlobSCV_IP[numSCVIP]
Definition fvho_geom.h:555
void update_boundary_faces(GridObject *elem, const MathVector< worldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
update boundary data for given element
Definition fvho_geom.cpp:1251
static const size_t numSCVFPerSubElem
number of SCVF per SubElement
Definition fvho_geom.h:105
static const int worldDim
dimension of world
Definition fvho_geom.h:88
size_t num_boundary_subsets()
number of registered boundary subsets
Definition fvho_geom.h:594
MathVector< dim > m_vLocSCV_IP[numSCVIP]
Definition fvho_geom.h:556
SCVF m_vSCVF[numSCVF]
SubControlVolumeFaces.
Definition fvho_geom.h:647
MathVector< worldDim > m_vGlobSCVF_IP[numSCVFIP]
Definition fvho_geom.h:551
void update(GridObject *elem, const MathVector< worldDim > *vCornerCoords, const ISubsetHandler *ish=NULL)
update data for given element
Definition fvho_geom.cpp:1134
static const int p
small abbreviation for order
Definition fvho_geom.h:75
std::vector< BF > m_vEmptyVectorBF
Definition fvho_geom.h:637
const MathVector< dim > * scvf_local_ips() const
returns all ips of scvf as they appear in scv loop
Definition fvho_geom.h:534
std::map< int, std::vector< BF > > m_mapVectorBF
Definition fvho_geom.h:636
size_t num_scvf() const
number of SubControlVolumeFaces
Definition fvho_geom.h:513
traits::scv_type scv_type
type of SubControlVolume
Definition fvho_geom.h:121
void add_boundary_subset(int subsetIndex)
add subset that is interpreted as boundary subset.
Definition fvho_geom.h:585
static const size_t nsh
number of shape functions
Definition fvho_geom.h:143
ReferenceMapping< ref_elem_type, worldDim > m_rMapping
Reference Mapping.
Definition fvho_geom.h:653
static const size_t numSCVFIP
number of scvf ip
Definition fvho_geom.h:117
void clear_boundary_subsets()
reset all boundary subsets
Definition fvho_geom.h:591
const SCVF & scvf(size_t i) const
const access to SubControlVolumeFace number i
Definition fvho_geom.h:516
const local_shape_fct_set_type & m_rTrialSpace
Shape function set.
Definition fvho_geom.h:659
GaussQuadrature< scvf_type, quadOrderSCVF > scvf_quad_rule_type
type of quadrature rule
Definition fvho_geom.h:114
static const int order
order
Definition fvho_geom.h:92
static const size_t numSCV
number of SubControlVolumes
Definition fvho_geom.h:127
static const int quadOrderSCV
quadrature order
Definition fvho_geom.h:130
GaussQuadrature< scv_type, quadOrderSCV > scv_quad_rule_type
type of quadrature rule
Definition fvho_geom.h:133
const MathVector< dim > * scv_local_ips() const
returns all ips of scv as they appear in scv loop
Definition fvho_geom.h:543
static const size_t numSCVPerSubElem
number of SCV per SubElement
Definition fvho_geom.h:124
void remove_boundary_subset(int subsetIndex)
removes subset that is interpreted as boundary subset.
Definition fvho_geom.h:588
SubElement m_vSubElem[numSubElem]
subelements
Definition fvho_geom.h:581
size_t num_bf() const
number of all boundary faces
Definition fvho_geom.h:597
MathVector< dim > m_vLocSCVF_IP[numSCVFIP]
Definition fvho_geom.h:552
static const int dim
dimension of reference element
Definition fvho_geom.h:85
fv1_traits< ref_elem_type, worldDim > traits
traits used
Definition fvho_geom.h:99
LagrangeLSFS< ref_elem_type, p > local_shape_fct_set_type
type of Shape function used
Definition fvho_geom.h:140
static const size_t numSCVIP
number of scv ip
Definition fvho_geom.h:136
size_t num_sh() const
number of shape functions
Definition fvho_geom.h:527
const scv_quad_rule_type & m_rSCVQuadRule
Quad Rule scv.
Definition fvho_geom.h:665
size_t num_scvf_ips() const
returns number of all scvf ips
Definition fvho_geom.h:531
static const size_t numSCVF
number of SubControlVolumeFaces
Definition fvho_geom.h:108
void reset_curr_elem()
Definition fvho_geom.h:633
static const size_t numSubElem
number of subelements
Definition fvho_geom.h:95
const ref_elem_type & m_rRefElem
Reference Element.
Definition fvho_geom.h:656
void update_local(ReferenceObjectID roid, const LFEID &lfeID=LFEID(LFEID::LAGRANGE, worldDim, 1), size_t orderQuad=TQuadOrder)
update Geometry for roid
Definition fvho_geom.cpp:968
MathVector< dim > m_vLocCorner[ref_elem_type::numCorners]
corners of reference element
Definition fvho_geom.h:644
reference_element_traits< TElem >::reference_element_type ref_elem_type
type of reference element
Definition fvho_geom.h:82
size_t num_bf(int si) const
number of boundary faces on subset 'subsetIndex'
Definition fvho_geom.h:607
traits::scvf_type scvf_type
type of SubControlVolumeFace
Definition fvho_geom.h:102
TElem elem_type
type of element
Definition fvho_geom.h:79
const MathVector< worldDim > * scvf_global_ips() const
returns all ips of scvf as they appear in scv loop
Definition fvho_geom.h:537
static const int quadOrderSCVF
quadrature order
Definition fvho_geom.h:111
const MathVector< worldDim > * scv_global_ips() const
returns all ips of scv as they appear in scv loop
Definition fvho_geom.h:546
const std::vector< BF > & bf(int si) const
returns reference to vector of boundary faces for subsetIndex
Definition fvho_geom.h:625
size_t num_scv_ips() const
returns number of all scv ips
Definition fvho_geom.h:540
static const int maxMid
Definition fvho_geom.h:560
const SCV & scv(size_t i) const
const access to SubControlVolume number i
Definition fvho_geom.h:523
const BF & bf(int si, size_t i) const
returns the boundary face i for subsetIndex
Definition fvho_geom.h:616
static const bool usesHangingNodes
Hanging node flag: this Geometry does not support hanging nodes.
Definition fvho_geom.h:146
fixed order gauss quadrature
Definition gauss_quad.h:44
The base class for all geometric objects, such as vertices, edges, faces, volumes,...
Definition grid_base_objects.h:157
Definition subset_handler_interface.h:223
Identifier for Local Finite Elements.
Definition local_finite_element_id.h:98
int order() const
returns the order of the local finite element
Definition local_finite_element_id.h:126
@ LAGRANGE
Definition local_finite_element_id.h:104
Lagrange Shape Function Set without virtual functions and fixed order.
Definition lagrange.h:48
A class for fixed size, dense matrices.
Definition math_matrix.h:63
a mathematical Vector with N entries.
Definition math_vector.h:97
Definition reference_mapping.h:65
#define UG_ASSERT(expr, msg)
Definition assert.h:70
#define UG_THROW(msg)
Definition error.h:57
double number
Definition types.h:124
ReferenceObjectID
these ids are used to identify the shape of a geometric object.
Definition grid_base_objects.h:74
Definition fvho_geom.h:1120
MathVector< worldDim > vvGloMid[dim+1][maxMid]
Definition fvho_geom.h:1127
size_t vDoFID[traits::maxNumSCV]
Definition fvho_geom.h:1122
std::vector< int > vElemBndSide
Definition fvho_geom.h:1133
bool isBndElem
Definition fvho_geom.h:1130
MathVector< dim > vvLocMid[dim+1][maxMid]
Definition fvho_geom.h:1126
Definition fvho_geom.h:564
MathVector< worldDim > vvGloMid[dim+1][maxMid]
Definition fvho_geom.h:571
bool isBndElem
Definition fvho_geom.h:574
std::vector< int > vElemBndSide
Definition fvho_geom.h:577
size_t vDoFID[ref_elem_type::numCorners]
Definition fvho_geom.h:566
MathVector< dim > vvLocMid[dim+1][maxMid]
Definition fvho_geom.h:570
helper class to store dimension and id of a midpoint of a sub-element
Definition fv_geom_base.h:41
Definition metaprogramming_util.h:165
Traits for Finite Volumes for a generic element of the fixed dimensionalities.
Definition fv_util.h:1325
Traits for Finite Volumes (dummy implementation, s. the instantiations below)
Definition fv_util.h:178
traits for reference elements
Definition reference_element_traits.h:48