v0.14.0
Loading...
Searching...
No Matches
Classes | Public Types | Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes | Friends | List of all members
MoFEM::ForcesAndSourcesCore Struct Reference

structure to get information form mofem into EntitiesFieldData More...

#include <src/finite_elements/ForcesAndSourcesCore.hpp>

Inheritance diagram for MoFEM::ForcesAndSourcesCore:
[legend]
Collaboration diagram for MoFEM::ForcesAndSourcesCore:
[legend]

Classes

struct  UserDataOperator
 

Public Types

typedef boost::function< int(int order_row, int order_col, int order_data)> RuleHookFun
 
typedef boost::function< MoFEMErrorCode(ForcesAndSourcesCore *fe_raw_ptr, int order_row, int order_col, int order_data)> GaussHookFun
 
- Public Types inherited from MoFEM::KspMethod
enum  KSPContext { CTX_SETFUNCTION , CTX_OPERATORS , CTX_KSPNONE }
 pass information about context of KSP/DM for with finite element is computed More...
 
- Public Types inherited from MoFEM::PetscData
enum  DataContext {
  CTX_SET_NONE = 0 , CTX_SET_F = 1 << 0 , CTX_SET_A = 1 << 1 , CTX_SET_B = 1 << 2 ,
  CTX_SET_X = 1 << 3 , CTX_SET_X_T = 1 << 4 , CTX_SET_X_TT = 1 << 6 , CTX_SET_TIME = 1 << 7
}
 
using Switches = std::bitset< 8 >
 
- Public Types inherited from MoFEM::SnesMethod
enum  SNESContext { CTX_SNESSETFUNCTION , CTX_SNESSETJACOBIAN , CTX_SNESNONE }
 
- Public Types inherited from MoFEM::TSMethod
enum  TSContext {
  CTX_TSSETRHSFUNCTION , CTX_TSSETRHSJACOBIAN , CTX_TSSETIFUNCTION , CTX_TSSETIJACOBIAN ,
  CTX_TSTSMONITORSET , CTX_TSNONE
}
 

Public Member Functions

 ForcesAndSourcesCore (Interface &m_field)
 
boost::ptr_deque< UserDataOperator > & getOpPtrVector ()
 Use to push back operator for row operator. More...
 
auto & getElementPolynomialBase ()
 Get the Entity Polynomial Base object. More...
 
auto & getUserPolynomialBase ()
 Get the User Polynomial Base object. More...
 
virtual MoFEMErrorCode preProcess ()
 function is run at the beginning of loop More...
 
virtual MoFEMErrorCode operator() ()
 function is run for every finite element More...
 
virtual MoFEMErrorCode postProcess ()
 function is run at the end of loop More...
 
int getMaxDataOrder () const
 Get max order of approximation for data fields. More...
 
int getMaxRowOrder () const
 Get max order of approximation for field in rows. More...
 
int getMaxColOrder () const
 Get max order of approximation for field in columns. More...
 
auto & getEntData (const FieldSpace space, const EntityType type, const int side)
 Get the entity data. More...
 
auto & getDataOnElementBySpaceArray ()
 Get data on entities and space. More...
 
auto & getDerivedDataOnElementBySpaceArray ()
 Get derived data on entities and space. More...
 
- Public Member Functions inherited from MoFEM::FEMethod
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 FEMethod ()=default
 
auto getFEName () const
 get finite element name More...
 
auto getDataDofsPtr () const
 
auto getDataVectorDofsPtr () const
 
const FieldEntity_vector_viewgetDataFieldEnts () const
 
boost::shared_ptr< FieldEntity_vector_view > & getDataFieldEntsPtr () const
 
auto & getRowFieldEnts () const
 
auto & getRowFieldEntsPtr () const
 
auto & getColFieldEnts () const
 
auto & getColFieldEntsPtr () const
 
auto getRowDofsPtr () const
 
auto getColDofsPtr () const
 
auto getNumberOfNodes () const
 
EntityHandle getFEEntityHandle () const
 
MoFEMErrorCode getNodeData (const std::string field_name, VectorDouble &data, const bool reset_dofs=true)
 
- Public Member Functions inherited from MoFEM::BasicMethod
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 BasicMethod ()
 
virtual ~BasicMethod ()=default
 
int getNinTheLoop () const
 get number of evaluated element in the loop More...
 
int getLoopSize () const
 get loop size More...
 
auto getLoHiFERank () const
 Get lo and hi processor rank of iterated entities. More...
 
auto getLoFERank () const
 Get upper rank in loop for iterating elements. More...
 
auto getHiFERank () const
 Get upper rank in loop for iterating elements. More...
 
unsigned int getFieldBitNumber (std::string field_name) const
 
MoFEMErrorCode copyBasicMethod (const BasicMethod &basic)
 Copy data from other base method to this base method. More...
 
virtual MoFEMErrorCode preProcess ()
 function is run at the beginning of loop More...
 
virtual MoFEMErrorCode operator() ()
 function is run for every finite element More...
 
virtual MoFEMErrorCode postProcess ()
 function is run at the end of loop More...
 
boost::weak_ptr< CacheTuplegetCacheWeakPtr () const
 Get the cache weak ptr object. More...
 
- Public Member Functions inherited from MoFEM::KspMethod
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 KspMethod ()
 
virtual ~KspMethod ()=default
 
MoFEMErrorCode copyKsp (const KspMethod &ksp)
 copy data form another method More...
 
- Public Member Functions inherited from MoFEM::PetscData
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 PetscData ()
 
virtual ~PetscData ()=default
 
MoFEMErrorCode copyPetscData (const PetscData &petsc_data)
 
- Public Member Functions inherited from MoFEM::UnknownInterface
virtual MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const =0
 
template<class IFACE >
MoFEMErrorCode registerInterface (bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE *&iface) const
 Get interface refernce to pointer of interface. More...
 
template<class IFACE >
MoFEMErrorCode getInterface (IFACE **const iface) const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_pointer< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get interface pointer to pointer of interface. More...
 
template<class IFACE , typename boost::enable_if< boost::is_reference< IFACE >, int >::type = 0>
IFACE getInterface () const
 Get reference to interface. More...
 
template<class IFACE >
IFACE * getInterface () const
 Function returning pointer to interface. More...
 
virtual ~UnknownInterface ()=default
 
- Public Member Functions inherited from MoFEM::SnesMethod
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 SnesMethod ()
 
virtual ~SnesMethod ()=default
 
MoFEMErrorCode copySnes (const SnesMethod &snes)
 Copy snes data. More...
 
- Public Member Functions inherited from MoFEM::TSMethod
MoFEMErrorCode query_interface (boost::typeindex::type_index type_index, UnknownInterface **iface) const
 
 TSMethod ()
 
virtual ~TSMethod ()=default
 
MoFEMErrorCode copyTs (const TSMethod &ts)
 Copy TS solver data. More...
 

Public Attributes

InterfacemField
 
RuleHookFun getRuleHook
 Hook to get rule. More...
 
GaussHookFun setRuleHook
 Set function to calculate integration rule. More...
 
MatrixDouble gaussPts
 Matrix of integration points. More...
 
- Public Attributes inherited from MoFEM::FEMethod
std::string feName
 Name of finite element. More...
 
boost::shared_ptr< const NumeredEntFiniteElementnumeredEntFiniteElementPtr
 
boost::function< bool(FEMethod *fe_method_ptr)> exeTestHook
 Tet if element to skip element. More...
 
- Public Attributes inherited from MoFEM::BasicMethod
int nInTheLoop
 number currently of processed method More...
 
int loopSize
 local number oe methods to process More...
 
std::pair< int, int > loHiFERank
 Llo and hi processor rank of iterated entities. More...
 
int rAnk
 processor rank More...
 
int sIze
 number of processors in communicator More...
 
const RefEntity_multiIndexrefinedEntitiesPtr
 container of mofem dof entities More...
 
const RefElement_multiIndexrefinedFiniteElementsPtr
 container of mofem finite element entities More...
 
const ProblemproblemPtr
 raw pointer to problem More...
 
const Field_multiIndexfieldsPtr
 raw pointer to fields container More...
 
const FieldEntity_multiIndexentitiesPtr
 raw pointer to container of field entities More...
 
const DofEntity_multiIndexdofsPtr
 raw pointer container of dofs More...
 
const FiniteElement_multiIndexfiniteElementsPtr
 raw pointer to container finite elements More...
 
const EntFiniteElement_multiIndexfiniteElementsEntitiesPtr
 
const FieldEntityEntFiniteElementAdjacencyMap_multiIndexadjacenciesPtr
 
boost::function< MoFEMErrorCode()> preProcessHook
 Hook function for pre-processing. More...
 
boost::function< MoFEMErrorCode()> postProcessHook
 Hook function for post-processing. More...
 
boost::function< MoFEMErrorCode()> operatorHook
 Hook function for operator. More...
 
boost::movelib::unique_ptr< boolvecAssembleSwitch
 
boost::movelib::unique_ptr< boolmatAssembleSwitch
 
boost::weak_ptr< CacheTuplecacheWeakPtr
 
- Public Attributes inherited from MoFEM::KspMethod
KSPContext ksp_ctx
 Context. More...
 
KSP ksp
 KSP solver. More...
 
Vec & ksp_f
 
Mat & ksp_A
 
Mat & ksp_B
 
- Public Attributes inherited from MoFEM::PetscData
Switches data_ctx
 
Vec f
 
Mat A
 
Mat B
 
Vec x
 
Vec x_t
 
Vec x_tt
 
- Public Attributes inherited from MoFEM::SnesMethod
SNESContext snes_ctx
 
SNES snes
 snes solver More...
 
Vec & snes_x
 state vector More...
 
Vec & snes_f
 residual More...
 
Mat & snes_A
 jacobian matrix More...
 
Mat & snes_B
 preconditioner of jacobian matrix More...
 
- Public Attributes inherited from MoFEM::TSMethod
TS ts
 time solver More...
 
TSContext ts_ctx
 
PetscInt ts_step
 time step number More...
 
PetscReal ts_a
 shift for U_t (see PETSc Time Solver) More...
 
PetscReal ts_aa
 shift for U_tt shift for U_tt More...
 
PetscReal ts_t
 time More...
 
PetscReal ts_dt
 time step size More...
 
Vec & ts_u
 state vector More...
 
Vec & ts_u_t
 time derivative of state vector More...
 
Vec & ts_u_tt
 second time derivative of state vector More...
 
Vec & ts_F
 residual vector More...
 
Mat & ts_A
 
Mat & ts_B
 Preconditioner for ts_A. More...
 

Protected Member Functions

MoFEMErrorCode getEntitySense (const EntityType type, boost::ptr_vector< EntitiesFieldData::EntData > &data) const
 get sense (orientation) of entity More...
 
MoFEMErrorCode getEntityDataOrder (const EntityType type, const FieldSpace space, boost::ptr_vector< EntitiesFieldData::EntData > &data) const
 Get the entity data order. More...
 
template<EntityType type>
MoFEMErrorCode getEntitySense (EntitiesFieldData &data) const
 Get the entity sense (orientation) More...
 
template<EntityType type>
MoFEMErrorCode getEntityDataOrder (EntitiesFieldData &data, const FieldSpace space) const
 Get the entity data order for given space. More...
 
MoFEMErrorCode getFaceNodes (EntitiesFieldData &data) const
 Get nodes on faces. More...
 
MoFEMErrorCode getSpacesAndBaseOnEntities (EntitiesFieldData &data) const
 Get field approximation space and base on entities. More...
 
virtual int getRule (int order_row, int order_col, int order_data)
 another variant of getRule More...
 
virtual MoFEMErrorCode setGaussPts (int order_row, int order_col, int order_data)
 set user specific integration rule More...
 
MoFEMErrorCode calHierarchicalBaseFunctionsOnElement (const FieldApproximationBase b)
 Calculate base functions. More...
 
MoFEMErrorCode calHierarchicalBaseFunctionsOnElement ()
 Calculate base functions. More...
 
MoFEMErrorCode calBernsteinBezierBaseFunctionsOnElement ()
 Calculate Bernstein-Bezier base. More...
 
MoFEMErrorCode createDataOnElement (EntityType type)
 Create a entity data on element object. More...
 
MoFEMErrorCode loopOverOperators ()
 Iterate user data operators. More...
 
Indices
template<typename EXTRACTOR >
MoFEMErrorCode getNodesIndices (const int bit_number, FieldEntity_vector_view &ents_field, VectorInt &nodes_indices, VectorInt &local_nodes_indices, EXTRACTOR &&extractor) const
 get node indices More...
 
MoFEMErrorCode getRowNodesIndices (EntitiesFieldData &data, const int bit_number) const
 get row node indices from FENumeredDofEntity_multiIndex More...
 
MoFEMErrorCode getColNodesIndices (EntitiesFieldData &data, const int bit_number) const
 get col node indices from FENumeredDofEntity_multiIndex More...
 
template<typename EXTRACTOR >
MoFEMErrorCode getEntityIndices (EntitiesFieldData &data, const int bit_number, FieldEntity_vector_view &ents_field, const EntityType type_lo, const EntityType type_hi, EXTRACTOR &&extractor) const
 
MoFEMErrorCode getEntityRowIndices (EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
 
MoFEMErrorCode getEntityColIndices (EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
 
MoFEMErrorCode getNoFieldIndices (const int bit_number, boost::shared_ptr< FENumeredDofEntity_multiIndex > dofs, VectorInt &nodes_indices) const
 get NoField indices More...
 
MoFEMErrorCode getNoFieldRowIndices (EntitiesFieldData &data, const int bit_number) const
 get col NoField indices More...
 
MoFEMErrorCode getNoFieldColIndices (EntitiesFieldData &data, const int bit_number) const
 get col NoField indices More...
 
Data
MoFEMErrorCode getBitRefLevelOnData ()
 
MoFEMErrorCode getNoFieldFieldData (const int bit_number, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs, VectorFieldEntities &ent_field) const
 Get field data on nodes. More...
 
MoFEMErrorCode getNoFieldFieldData (EntitiesFieldData &data, const int bit_number) const
 
MoFEMErrorCode getNodesFieldData (EntitiesFieldData &data, const int bit_number) const
 Get data on nodes. More...
 
MoFEMErrorCode getEntityFieldData (EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
 
Data form NumeredDofEntity_multiIndex
MoFEMErrorCode getProblemNodesIndices (const std::string &field_name, const NumeredDofEntity_multiIndex &dofs, VectorInt &nodes_indices) const
 get indices of nodal indices which are declared for problem but not this particular element More...
 
MoFEMErrorCode getProblemTypeIndices (const std::string &field_name, const NumeredDofEntity_multiIndex &dofs, EntityType type, int side_number, VectorInt &indices) const
 get indices by type (generic function) which are declared for problem but not this particular element More...
 
MoFEMErrorCode getProblemNodesRowIndices (const std::string &field_name, VectorInt &nodes_indices) const
 
MoFEMErrorCode getProblemTypeRowIndices (const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
 
MoFEMErrorCode getProblemNodesColIndices (const std::string &field_name, VectorInt &nodes_indices) const
 
MoFEMErrorCode getProblemTypeColIndices (const std::string &field_name, EntityType type, int side_number, VectorInt &indices) const
 
Deprecated (do not use)
virtual int getRule (int order)
 
virtual MoFEMErrorCode setGaussPts (int order)
 

Protected Attributes

const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACEdataOnElement
 Entity data on element entity rows fields. More...
 
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACEderivedDataOnElement
 Entity data on element entity columns fields. More...
 
EntitiesFieldDatadataNoField
 
EntitiesFieldDatadataH1
 
EntitiesFieldDatadataHcurl
 
EntitiesFieldDatadataHdiv
 
EntitiesFieldDatadataL2
 
boost::ptr_deque< UserDataOperatoropPtrVector
 Vector of finite element users data operators. More...
 
EntityType lastEvaluatedElementEntityType
 Last evaluated type of element entity. More...
 
MatrixDouble coordsAtGaussPts
 coordinated at gauss points More...
 
double elementMeasure
 

Private Member Functions

MoFEMErrorCode setSideFEPtr (const ForcesAndSourcesCore *side_fe_ptr)
 Set the pointer to face element on the side. More...
 
MoFEMErrorCode setRefineFEPtr (const ForcesAndSourcesCore *refine_fe_ptr)
 Set the pointer to face element refined. More...
 

Private Attributes

boost::shared_ptr< BaseFunctionelementPolynomialBasePtr
 Pointer to entity polynomial base. More...
 
boost::shared_ptr< BaseFunctionuserPolynomialBasePtr
 Pointer to user polynomail base. More...
 
ForcesAndSourcesCoresidePtrFE
 Element to integrate on the sides. More...
 
ForcesAndSourcesCorerefinePtrFE
 Element to integrate parent or child. More...
 

Friends

class UserDataOperator
 
class VolumeElementForcesAndSourcesCoreOnSide
 
class FaceElementForcesAndSourcesCoreOnSide
 
class FaceElementForcesAndSourcesCoreOnChildParent
 
class EdgeElementForcesAndSourcesCoreOnChildParent
 
class VolumeElementForcesAndSourcesCoreOnContactPrismSide
 

Additional Inherited Members

- Static Public Member Functions inherited from MoFEM::UnknownInterface
static MoFEMErrorCode getLibVersion (Version &version)
 Get library version. More...
 
static MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version)
 Get database major version. More...
 
static MoFEMErrorCode setFileVersion (moab::Interface &moab, Version version=Version(MoFEM_VERSION_MAJOR, MoFEM_VERSION_MINOR, MoFEM_VERSION_BUILD))
 Get database major version. More...
 
static MoFEMErrorCode getInterfaceVersion (Version &version)
 Get database major version. More...
 
- Static Public Attributes inherited from MoFEM::PetscData
static constexpr Switches CtxSetNone = PetscData::Switches(CTX_SET_NONE)
 
static constexpr Switches CtxSetF = PetscData::Switches(CTX_SET_F)
 
static constexpr Switches CtxSetA = PetscData::Switches(CTX_SET_A)
 
static constexpr Switches CtxSetB = PetscData::Switches(CTX_SET_B)
 
static constexpr Switches CtxSetX = PetscData::Switches(CTX_SET_X)
 
static constexpr Switches CtxSetX_T = PetscData::Switches(CTX_SET_X_T)
 
static constexpr Switches CtxSetX_TT = PetscData::Switches(CTX_SET_X_TT)
 
static constexpr Switches CtxSetTime = PetscData::Switches(CTX_SET_TIME)
 

Detailed Description

structure to get information form mofem into EntitiesFieldData

Examples
dynamic_first_order_con_law.cpp, hanging_node_approx.cpp, hello_world.cpp, prism_elements_from_surface.cpp, prism_polynomial_approximation.cpp, and quad_polynomial_approximation.cpp.

Definition at line 22 of file ForcesAndSourcesCore.hpp.

Member Typedef Documentation

◆ GaussHookFun

typedef boost::function<MoFEMErrorCode(ForcesAndSourcesCore *fe_raw_ptr, int order_row, int order_col, int order_data)> MoFEM::ForcesAndSourcesCore::GaussHookFun

Definition at line 33 of file ForcesAndSourcesCore.hpp.

◆ RuleHookFun

typedef boost::function<int(int order_row, int order_col, int order_data)> MoFEM::ForcesAndSourcesCore::RuleHookFun

Definition at line 28 of file ForcesAndSourcesCore.hpp.

Constructor & Destructor Documentation

◆ ForcesAndSourcesCore()

MoFEM::ForcesAndSourcesCore::ForcesAndSourcesCore ( Interface m_field)

Definition at line 40 of file ForcesAndSourcesCore.cpp.

41 :
42
43 mField(m_field), getRuleHook(0), setRuleHook(0),
45
46 boost::make_shared<EntitiesFieldData>(MBENTITYSET), // NOSPACE,
47 boost::make_shared<EntitiesFieldData>(MBENTITYSET), // NOFIELD
48 boost::make_shared<EntitiesFieldData>(MBENTITYSET), // H1
49 boost::make_shared<EntitiesFieldData>(MBENTITYSET), // HCURL
50 boost::make_shared<EntitiesFieldData>(MBENTITYSET), // HDIV
51 boost::make_shared<EntitiesFieldData>(MBENTITYSET) // L2
52
53 },
55
56 nullptr,
57 boost::make_shared<DerivedEntitiesFieldData>(
58 dataOnElement[NOFIELD]), // NOFIELD
59 boost::make_shared<DerivedEntitiesFieldData>(dataOnElement[H1]), // H1
60 boost::make_shared<DerivedEntitiesFieldData>(
61 dataOnElement[HCURL]), // HCURL
62 boost::make_shared<DerivedEntitiesFieldData>(
63 dataOnElement[HDIV]), // HDIV
64 boost::make_shared<DerivedEntitiesFieldData>(dataOnElement[L2]) // L2
65
66 },
70 lastEvaluatedElementEntityType(MBMAXTYPE), sidePtrFE(nullptr),
71 refinePtrFE(nullptr) {
72
73 dataOnElement[NOSPACE]->dataOnEntities[MBENTITYSET].push_back(
74 new EntitiesFieldData::EntData());
75
76 dataOnElement[NOFIELD]->dataOnEntities[MBENTITYSET].push_back(
77 new EntitiesFieldData::EntData());
78 derivedDataOnElement[NOFIELD]->dataOnEntities[MBENTITYSET].push_back(
79 new EntitiesFieldData::EntData());
80}
@ L2
field with C-1 continuity
Definition: definitions.h:88
@ NOFIELD
scalar or vector of scalars describe (no true field)
Definition: definitions.h:84
@ H1
continuous field
Definition: definitions.h:85
@ NOSPACE
Definition: definitions.h:83
@ HCURL
field with continuous tangents
Definition: definitions.h:86
@ HDIV
field with continuous normal traction
Definition: definitions.h:87
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > derivedDataOnElement
Entity data on element entity columns fields.
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
EntityType lastEvaluatedElementEntityType
Last evaluated type of element entity.
GaussHookFun setRuleHook
Set function to calculate integration rule.
ForcesAndSourcesCore * refinePtrFE
Element to integrate parent or child.
RuleHookFun getRuleHook
Hook to get rule.

Member Function Documentation

◆ calBernsteinBezierBaseFunctionsOnElement()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::calBernsteinBezierBaseFunctionsOnElement ( )
protected

Calculate Bernstein-Bezier base.

Returns
MoFEMErrorCode

Definition at line 1106 of file ForcesAndSourcesCore.cpp.

1106 {
1108
1109 const auto ele_type = numeredEntFiniteElementPtr->getEntType();
1110
1111 auto get_nodal_base_data = [&](EntitiesFieldData &data, auto field_ptr) {
1113 auto &space = data.dataOnEntities[MBVERTEX][0].getSpace();
1114 auto &base = data.dataOnEntities[MBVERTEX][0].getBase();
1115 auto &bb_node_order = data.dataOnEntities[MBVERTEX][0].getBBNodeOrder();
1116
1117 auto &field_ents = getDataFieldEnts();
1118 auto bit_number = field_ptr->getBitNumber();
1119 const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
1120 bit_number, get_id_for_min_type<MBVERTEX>());
1121 auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
1122 cmp_uid_lo);
1123 if (lo != field_ents.end()) {
1124 const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
1125 bit_number, get_id_for_max_type<MBVERTEX>());
1126 auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
1127 if (lo != hi) {
1128
1129 for (auto it = lo; it != hi; ++it)
1130 if (auto first_e = it->lock()) {
1131 space = first_e->getSpace();
1132 base = first_e->getApproxBase();
1133 const int num_nodes = getNumberOfNodes();
1134 bb_node_order.resize(num_nodes, false);
1135 bb_node_order.clear();
1136
1137 std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
1138
1139 for (; it != hi; ++it) {
1140 if (auto e = it->lock()) {
1141 const auto &sn = e->getSideNumberPtr();
1142 const int side_number = sn->side_number;
1143 const int brother_side_number = sn->brother_side_number;
1144 if (brother_side_number != -1)
1145 brother_ents_vec.emplace_back(e);
1146 bb_node_order[side_number] = e->getMaxOrder();
1147 }
1148 }
1149
1150 for (auto &it : brother_ents_vec) {
1151 if (const auto e = it.lock()) {
1152 const auto &sn = e->getSideNumberPtr();
1153 const int side_number = sn->side_number;
1154 const int brother_side_number = sn->brother_side_number;
1155 bb_node_order[brother_side_number] = bb_node_order[side_number];
1156 }
1157 }
1158
1159 break;
1160 }
1161 }
1162 }
1163
1165 };
1166
1167 auto get_entity_base_data = [&](EntitiesFieldData &data, auto field_ptr,
1168 const EntityType type_lo,
1169 const EntityType type_hi) {
1171 for (EntityType t = MBEDGE; t != MBPOLYHEDRON; ++t) {
1172 for (auto &dat : data.dataOnEntities[t]) {
1173 dat.getOrder() = 0;
1174 dat.getBase() = NOBASE;
1175 dat.getSpace() = NOSPACE;
1176 dat.getFieldData().resize(0, false);
1177 dat.getFieldDofs().resize(0, false);
1178 }
1179 }
1180
1181 auto &field_ents = getDataFieldEnts();
1182 auto bit_number = field_ptr->getBitNumber();
1183 const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
1184 bit_number, get_id_for_min_type(type_lo));
1185 auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
1186 cmp_uid_lo);
1187 if (lo != field_ents.end()) {
1188 const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
1189 bit_number, get_id_for_max_type(type_hi));
1190 auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
1191
1192 std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
1193 for (; lo != hi; ++lo) {
1194 if (auto e = lo->lock()) {
1195 if (auto cache = e->entityCacheDataDofs.lock()) {
1196 if (cache->loHi[0] != cache->loHi[1]) {
1197 if (const auto side = e->getSideNumberPtr()->side_number;
1198 side >= 0) {
1199 const EntityType type = e->getEntType();
1200 auto &dat = data.dataOnEntities[type][side];
1201 const int brother_side =
1202 e->getSideNumberPtr()->brother_side_number;
1203 if (brother_side != -1)
1204 brother_ents_vec.emplace_back(e);
1205 dat.getBase() = e->getApproxBase();
1206 dat.getSpace() = e->getSpace();
1207 const auto ent_order = e->getMaxOrder();
1208 dat.getOrder() =
1209 dat.getOrder() > ent_order ? dat.getOrder() : ent_order;
1210 }
1211 }
1212 }
1213 }
1214 }
1215
1216 for (auto &ent_ptr : brother_ents_vec) {
1217 if (auto e = ent_ptr.lock()) {
1218 const EntityType type = e->getEntType();
1219 const int side = e->getSideNumberPtr()->side_number;
1220 const int brother_side = e->getSideNumberPtr()->brother_side_number;
1221 auto &dat = data.dataOnEntities[type][side];
1222 auto &dat_brother = data.dataOnEntities[type][brother_side];
1223 dat_brother.getBase() = dat.getBase();
1224 dat_brother.getSpace() = dat.getSpace();
1225 dat_brother.getOrder() = dat.getOrder();
1226 }
1227 }
1228 }
1230 };
1231
1232 for (auto &e : getDataFieldEnts()) {
1233 if (auto ent_data_ptr = e.lock()) {
1234 if (ent_data_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
1235 auto space = ent_data_ptr->getSpace();
1236 for (EntityType t = MBVERTEX; t != MBPOLYHEDRON; ++t) {
1237 for (auto &dat : (*dataOnElement[space]).dataOnEntities[t]) {
1238 for (auto &ptr : dat.getBBAlphaIndicesByOrderArray())
1239 ptr.reset();
1240 for (auto &ptr : dat.getBBNByOrderArray())
1241 ptr.reset();
1242 for (auto &ptr : dat.getBBDiffNByOrderArray())
1243 ptr.reset();
1244 }
1245 }
1246 }
1247 }
1248 }
1249
1250 auto check_space = [&](const auto space) {
1251 switch (space) {
1252 case H1:
1253 for (auto t = MBVERTEX; t <= ele_type; ++t) {
1254 if (dataOnElement[H1]->spacesOnEntities[t].test(H1))
1255 return true;
1256 }
1257 return false;
1258 case HCURL:
1259 for (auto t = MBEDGE; t <= ele_type; ++t) {
1260 if (dataOnElement[HCURL]->spacesOnEntities[t].test(HCURL))
1261 return true;
1262 }
1263 return false;
1264 case HDIV:
1265 for (auto t = MBTRI; t <= ele_type; ++t) {
1266 if (dataOnElement[HDIV]->spacesOnEntities[t].test(HDIV))
1267 return true;
1268 }
1269 return false;
1270 case L2:
1271 return dataOnElement[L2]->spacesOnEntities[ele_type].test(L2);
1272 break;
1273 default:
1274 THROW_MESSAGE("Not implemented");
1275 }
1276 };
1277
1278 std::set<string> fields_list;
1279 for (auto &e : getDataFieldEnts()) {
1280 if (auto ent_data_ptr = e.lock()) {
1281 if (ent_data_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
1282 auto field_name = ent_data_ptr->getName();
1283 if (fields_list.find(field_name) == fields_list.end()) {
1284 auto field_ptr = ent_data_ptr->getFieldRawPtr();
1285 auto space = ent_data_ptr->getSpace();
1286 CHKERR get_nodal_base_data(*dataOnElement[space], field_ptr);
1287 CHKERR get_entity_base_data(*dataOnElement[space], field_ptr, MBEDGE,
1288 ele_type);
1289 if (check_space(space)) {
1290 CHKERR getElementPolynomialBase()->getValue(
1291 gaussPts, boost::make_shared<EntPolynomialBaseCtx>(
1292 *dataOnElement[space], field_name,
1293 static_cast<FieldSpace>(space),
1295 fields_list.insert(field_name);
1296 }
1297 }
1298 }
1299 }
1300 }
1301
1303};
@ NOBASE
Definition: definitions.h:59
@ AINSWORTH_BERNSTEIN_BEZIER_BASE
Definition: definitions.h:64
FieldSpace
approximation spaces
Definition: definitions.h:82
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:346
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:416
#define CHKERR
Inline error check.
Definition: definitions.h:535
#define THROW_MESSAGE(msg)
Throw MoFEM exception.
Definition: definitions.h:561
BlockParamData * cache
static auto cmp_uid_lo(const boost::weak_ptr< FieldEntity > &a, const UId &b)
EntityHandle get_id_for_max_type()
static auto cmp_uid_hi(const UId &b, const boost::weak_ptr< FieldEntity > &a)
EntityHandle get_id_for_min_type()
constexpr double t
plate stiffness
Definition: plate.cpp:59
constexpr auto field_name
const FieldEntity_vector_view & getDataFieldEnts() const
auto getNumberOfNodes() const
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
auto & getElementPolynomialBase()
Get the Entity Polynomial Base object.
MatrixDouble gaussPts
Matrix of integration points.

◆ calHierarchicalBaseFunctionsOnElement() [1/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::calHierarchicalBaseFunctionsOnElement ( )
protected

Calculate base functions.

Returns
Error code

Use the some node base. Node base is usually used for construction other bases.

Definition at line 1086 of file ForcesAndSourcesCore.cpp.

1086 {
1088 /// Use the some node base. Node base is usually used for construction other
1089 /// bases.
1090 for (int space = HCURL; space != LASTSPACE; ++space) {
1091 dataOnElement[space]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
1092 dataOnElement[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
1093 dataOnElement[space]->dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(
1094 NOBASE) =
1095 dataOnElement[H1]->dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(
1096 NOBASE);
1097 }
1098 for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
1100 static_cast<FieldApproximationBase>(b));
1101 }
1103}
FieldApproximationBase
approximation base
Definition: definitions.h:58
@ LASTBASE
Definition: definitions.h:69
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:60
@ LASTSPACE
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:89
MoFEMErrorCode calHierarchicalBaseFunctionsOnElement()
Calculate base functions.

◆ calHierarchicalBaseFunctionsOnElement() [2/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::calHierarchicalBaseFunctionsOnElement ( const FieldApproximationBase  b)
protected

Calculate base functions.

Returns
Error code

Definition at line 1032 of file ForcesAndSourcesCore.cpp.

1033 {
1035 if (dataOnElement[H1]->bAse.test(b)) {
1036 switch (static_cast<FieldApproximationBase>(b)) {
1037 case NOBASE:
1038 break;
1040 // BERNSTEIN_BEZIER_BASE is not hierarchical base
1041 break;
1046 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1047 "Functions genrating approximation base not defined");
1048
1049 for (int space = H1; space != LASTSPACE; ++space) {
1050 if (dataOnElement[H1]->sPace.test(space) &&
1051 dataOnElement[H1]->bAse.test(b) &&
1052 dataOnElement[H1]->basesOnSpaces[space].test(b)) {
1053 CHKERR getElementPolynomialBase()->getValue(
1054 gaussPts,
1055 boost::make_shared<EntPolynomialBaseCtx>(
1056 *dataOnElement[space], static_cast<FieldSpace>(space),
1057 static_cast<FieldApproximationBase>(b), NOBASE));
1058 }
1059 }
1060 break;
1061 case USER_BASE:
1062 if (!getUserPolynomialBase())
1063 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1064 "Functions generating user approximation base not defined");
1065
1066 for (int space = H1; space != LASTSPACE; ++space)
1067 if (dataOnElement[H1]->sPace.test(space) &&
1068 dataOnElement[H1]->bAse.test(b) &&
1069 dataOnElement[H1]->basesOnSpaces[space].test(b)) {
1070 CHKERR getUserPolynomialBase()->getValue(
1071 gaussPts,
1072 boost::make_shared<EntPolynomialBaseCtx>(
1073 *dataOnElement[space], static_cast<FieldSpace>(space),
1074 static_cast<FieldApproximationBase>(b), NOBASE));
1075 }
1076 break;
1077 default:
1079 "Base <%s> not yet implemented",
1081 }
1082 }
1084}
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:62
@ USER_BASE
user implemented approximation base
Definition: definitions.h:68
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:66
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:31
static const char *const ApproximationBaseNames[]
Definition: definitions.h:72
virtual MPI_Comm & get_comm() const =0
auto & getUserPolynomialBase()
Get the User Polynomial Base object.

◆ createDataOnElement()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::createDataOnElement ( EntityType  type)
protected

Create a entity data on element object.

Returns
MoFEMErrorCode

Definition at line 1305 of file ForcesAndSourcesCore.cpp.

1305 {
1307
1308 // Data on elements for proper spaces
1309 for (int space = H1; space != LASTSPACE; ++space) {
1310 dataOnElement[space]->setElementType(type);
1311 derivedDataOnElement[space]->setElementType(type);
1312 }
1313
1315}

◆ getBitRefLevelOnData()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getBitRefLevelOnData ( )
protected

Get bit ref level in entities, and set it to data

Definition at line 595 of file ForcesAndSourcesCore.cpp.

595 {
597
598 for (auto &data : dataOnElement) {
599 if (data) {
600 for (auto &dat : data->dataOnEntities) {
601 for (auto &ent_dat : dat) {
602 ent_dat.getEntDataBitRefLevel().clear();
603 }
604 }
605 }
606 }
607
608 auto &field_ents = getDataFieldEnts();
609 for (auto it : field_ents) {
610 if (auto e = it.lock()) {
611 const FieldSpace space = e->getSpace();
612 if (space > NOFIELD) {
613 const EntityType type = e->getEntType();
614 const signed char side = e->getSideNumberPtr()->side_number;
615 if (side >= 0) {
616 if (auto &data = dataOnElement[space]) {
617 if (type == MBVERTEX) {
618 auto &dat = data->dataOnEntities[type][0];
619 dat.getEntDataBitRefLevel().resize(getNumberOfNodes(), false);
620 dat.getEntDataBitRefLevel()[side] = e->getBitRefLevel();
621 } else {
622 auto &dat = data->dataOnEntities[type][side];
623 dat.getEntDataBitRefLevel().resize(1, false);
624 dat.getEntDataBitRefLevel()[0] = e->getBitRefLevel();
625 }
626 }
627 }
628 } else {
629 if (auto &data = dataOnElement[NOFIELD]) {
630 auto &dat = data->dataOnEntities[MBENTITYSET][0];
631 dat.getEntDataBitRefLevel().resize(1, false);
632 dat.getEntDataBitRefLevel()[0] = e->getBitRefLevel();
633 }
634 }
635 }
636 }
637
639};

◆ getColNodesIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getColNodesIndices ( EntitiesFieldData data,
const int  bit_number 
) const
protected

get col node indices from FENumeredDofEntity_multiIndex

Definition at line 328 of file ForcesAndSourcesCore.cpp.

329 {
330
331 struct Extractor {
332 boost::weak_ptr<EntityCacheNumeredDofs>
333 operator()(boost::shared_ptr<FieldEntity> &e) {
334 return e->entityCacheColDofs;
335 }
336 };
337
338 return getNodesIndices(bit_number, getColFieldEnts(),
339 data.dataOnEntities[MBVERTEX][0].getIndices(),
340 data.dataOnEntities[MBVERTEX][0].getLocalIndices(),
341 Extractor());
342}
auto & getColFieldEnts() const
MoFEMErrorCode getNodesIndices(const int bit_number, FieldEntity_vector_view &ents_field, VectorInt &nodes_indices, VectorInt &local_nodes_indices, EXTRACTOR &&extractor) const
get node indices

◆ getDataOnElementBySpaceArray()

auto & MoFEM::ForcesAndSourcesCore::getDataOnElementBySpaceArray ( )
inline

Get data on entities and space.

Entities data are stored by space, by entity type, and entity side.

Returns
std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE>

Definition at line 156 of file ForcesAndSourcesCore.hpp.

156{ return dataOnElement; }

◆ getDerivedDataOnElementBySpaceArray()

auto & MoFEM::ForcesAndSourcesCore::getDerivedDataOnElementBySpaceArray ( )
inline

Get derived data on entities and space.

Entities data are stored by space, by entity type, and entity side. Derived data is used to store data on columns, so it shares infromatin about shape functions wih rows.

Returns
std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE>

Definition at line 167 of file ForcesAndSourcesCore.hpp.

167{ return derivedDataOnElement; }

◆ getElementPolynomialBase()

auto & MoFEM::ForcesAndSourcesCore::getElementPolynomialBase ( )
inline

Get the Entity Polynomial Base object.

Returns
boost::shared_ptr<BaseFunction>&&

Definition at line 90 of file ForcesAndSourcesCore.hpp.

boost::shared_ptr< BaseFunction > elementPolynomialBasePtr
Pointer to entity polynomial base.

◆ getEntData()

auto & MoFEM::ForcesAndSourcesCore::getEntData ( const FieldSpace  space,
const EntityType  type,
const int  side 
)
inline

Get the entity data.

Parameters
space
type
side
Returns
EntitiesFieldData::EntData&

Definition at line 144 of file ForcesAndSourcesCore.hpp.

145 {
146 return dataOnElement[space]->dataOnEntities[type][side];
147 }

◆ getEntityColIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntityColIndices ( EntitiesFieldData data,
const int  bit_number,
const EntityType  type_lo = MBVERTEX,
const EntityType  type_hi = MBPOLYHEDRON 
) const
protected

Definition at line 427 of file ForcesAndSourcesCore.cpp.

429 {
430
431 struct Extractor {
432 boost::weak_ptr<EntityCacheNumeredDofs>
433 operator()(boost::shared_ptr<FieldEntity> &e) {
434 return e->entityCacheColDofs;
435 }
436 };
437
438 return getEntityIndices(data, bit_number, getColFieldEnts(), type_lo, type_hi,
439 Extractor());
440}
MoFEMErrorCode getEntityIndices(EntitiesFieldData &data, const int bit_number, FieldEntity_vector_view &ents_field, const EntityType type_lo, const EntityType type_hi, EXTRACTOR &&extractor) const

◆ getEntityDataOrder() [1/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntityDataOrder ( const EntityType  type,
const FieldSpace  space,
boost::ptr_vector< EntitiesFieldData::EntData > &  data 
) const
protected

Get the entity data order.

Parameters
type
space
data
Returns
MoFEMErrorCode

Definition at line 139 of file ForcesAndSourcesCore.cpp.

141 {
143
144 auto set_order = [&]() {
146 auto &side_table = numeredEntFiniteElementPtr->getSideNumberTable();
147
148 for (unsigned int s = 0; s != data.size(); ++s)
149 data[s].getOrder() = 0;
150
151 const FieldEntity_vector_view &data_field_ent = getDataFieldEnts();
152
153 for (auto r = fieldsPtr->get<BitFieldId_space_mi_tag>().equal_range(space);
154 r.first != r.second; ++r.first) {
155
156 const auto field_bit_number = (*r.first)->getBitNumber();
157 const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
158 field_bit_number, get_id_for_min_type(type));
159 auto lo = std::lower_bound(data_field_ent.begin(), data_field_ent.end(),
160 lo_uid, cmp_uid_lo);
161 if (lo != data_field_ent.end()) {
162 const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
163 field_bit_number, get_id_for_max_type(type));
164 auto hi =
165 std::upper_bound(lo, data_field_ent.end(), hi_uid, cmp_uid_hi);
166 for (; lo != hi; ++lo) {
167
168 if (auto ptr = lo->lock()) {
169
170 auto &e = *ptr;
171 auto sit = side_table.find(e.getEnt());
172 if (sit != side_table.end()) {
173 auto &side = *sit;
174 if (const auto side_number = side->side_number;
175 side_number >= 0) {
176 ApproximationOrder ent_order = e.getMaxOrder();
177 auto &dat = data[side_number];
178 dat.getOrder() =
179 dat.getOrder() > ent_order ? dat.getOrder() : ent_order;
180 }
181 } else
182 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
183 "Entity on side of the element not found");
184 }
185 }
186 }
187 }
188
190 };
191
192 auto set_order_on_brother = [&]() {
194 auto &side_table =
195 numeredEntFiniteElementPtr->getSideNumberTable().get<0>();
196 auto sit = side_table.lower_bound(get_id_for_min_type(type));
197 if (sit != side_table.end()) {
198 auto hi_sit = side_table.upper_bound(get_id_for_max_type(type));
199 for (; sit != hi_sit; ++sit) {
200 const int brother_side_number = (*sit)->brother_side_number;
201 if (brother_side_number != -1) {
202 const int side_number = (*sit)->side_number;
203 data[brother_side_number].getOrder() = data[side_number].getOrder();
204 }
205 }
206 }
208 };
209
210 CHKERR set_order();
211 CHKERR set_order_on_brother();
212
214}
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:447
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:440
int ApproximationOrder
Approximation on the entity.
Definition: Types.hpp:26
std::vector< boost::weak_ptr< FieldEntity > > FieldEntity_vector_view
int r
Definition: sdf.py:8
const Field_multiIndex * fieldsPtr
raw pointer to fields container

◆ getEntityDataOrder() [2/2]

template<EntityType type>
MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntityDataOrder ( EntitiesFieldData data,
const FieldSpace  space 
) const
inlineprotected

Get the entity data order for given space.

Template Parameters
type
Parameters
data
space
Returns
MoFEMErrorCode

Definition at line 213 of file ForcesAndSourcesCore.hpp.

214 {
215 return getEntityDataOrder(type, space, data.dataOnEntities[type]);
216 }
MoFEMErrorCode getEntityDataOrder(const EntityType type, const FieldSpace space, boost::ptr_vector< EntitiesFieldData::EntData > &data) const
Get the entity data order.

◆ getEntityFieldData()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntityFieldData ( EntitiesFieldData data,
const int  bit_number,
const EntityType  type_lo = MBVERTEX,
const EntityType  type_hi = MBPOLYHEDRON 
) const
protected

Definition at line 770 of file ForcesAndSourcesCore.cpp.

772 {
774 for (EntityType t = type_lo; t != type_hi; ++t) {
775 for (auto &dat : data.dataOnEntities[t]) {
776 dat.getOrder() = 0;
777 dat.getBase() = NOBASE;
778 dat.getSpace() = NOSPACE;
779 dat.getFieldData().resize(0, false);
780 dat.getFieldDofs().resize(0, false);
781 dat.getFieldEntities().resize(0, false);
782 }
783 }
784
785 auto &field_ents = getDataFieldEnts();
786 const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
787 bit_number, get_id_for_min_type(type_lo));
788 auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
789 cmp_uid_lo);
790 if (lo != field_ents.end()) {
791 const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
792 bit_number, get_id_for_max_type(type_hi));
793 auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
794 if (lo != hi) {
795
796 std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
797
798 for (auto it = lo; it != hi; ++it)
799 if (auto e = it->lock()) {
800 auto side_ptr = e->getSideNumberPtr();
801 if (const auto side = side_ptr->side_number; side >= 0) {
802 const EntityType type = e->getEntType();
803 auto &dat = data.dataOnEntities[type][side];
804 auto &ent_field = dat.getFieldEntities();
805 auto &ent_field_dofs = dat.getFieldDofs();
806 auto &ent_field_data = dat.getFieldData();
807
808 const int brother_side = side_ptr->brother_side_number;
809 if (brother_side != -1)
810 brother_ents_vec.emplace_back(e);
811
812 dat.getBase() = e->getApproxBase();
813 dat.getSpace() = e->getSpace();
814 const int ent_order = e->getMaxOrder();
815 dat.getOrder() =
816 dat.getOrder() > ent_order ? dat.getOrder() : ent_order;
817
818 auto ent_data = e->getEntFieldData();
819 ent_field_data.resize(ent_data.size(), false);
820 noalias(ent_field_data) = ent_data;
821 ent_field_dofs.resize(ent_data.size(), false);
822 std::fill(ent_field_dofs.begin(), ent_field_dofs.end(), nullptr);
823 ent_field.resize(1, false);
824 ent_field[0] = e.get();
825 if (auto cache = e->entityCacheDataDofs.lock()) {
826 for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
827 ent_field_dofs[(*dit)->getEntDofIdx()] =
828 reinterpret_cast<FEDofEntity *>((*dit).get());
829 }
830 }
831 }
832 }
833
834 for (auto &it : brother_ents_vec) {
835 if (const auto e = it.lock()) {
836 const EntityType type = e->getEntType();
837 const int side = e->getSideNumberPtr()->side_number;
838 const int brother_side = e->getSideNumberPtr()->brother_side_number;
839 auto &dat = data.dataOnEntities[type][side];
840 auto &dat_brother = data.dataOnEntities[type][brother_side];
841 dat_brother.getBase() = dat.getBase();
842 dat_brother.getSpace() = dat.getSpace();
843 dat_brother.getOrder() = dat.getOrder();
844 dat_brother.getFieldData() = dat.getFieldData();
845 dat_brother.getFieldDofs() = dat.getFieldDofs();
846 dat_brother.getFieldEntities() = dat.getFieldEntities();
847 }
848 }
849 }
850 }
851
853}

◆ getEntityIndices()

template<typename EXTRACTOR >
MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntityIndices ( EntitiesFieldData data,
const int  bit_number,
FieldEntity_vector_view ents_field,
const EntityType  type_lo,
const EntityType  type_hi,
EXTRACTOR &&  extractor 
) const
protected

Definition at line 345 of file ForcesAndSourcesCore.cpp.

348 {
350
351 for (EntityType t = type_lo; t != type_hi; ++t) {
352 for (auto &dat : data.dataOnEntities[t]) {
353 dat.getIndices().resize(0, false);
354 dat.getLocalIndices().resize(0, false);
355 }
356 }
357
358 const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
359 bit_number, get_id_for_min_type(type_lo));
360 auto lo = std::lower_bound(ents_field.begin(), ents_field.end(), lo_uid,
361 cmp_uid_lo);
362 if (lo != ents_field.end()) {
363 const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
364 bit_number, get_id_for_max_type(type_hi));
365 auto hi = std::upper_bound(lo, ents_field.end(), hi_uid, cmp_uid_hi);
366 if (lo != hi) {
367
368 std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
369
370 for (auto it = lo; it != hi; ++it)
371 if (auto e = it->lock()) {
372
373 const EntityType type = e->getEntType();
374 auto side_ptr = e->getSideNumberPtr();
375 if (const auto side = side_ptr->side_number; side >= 0) {
376 const auto nb_dofs_on_ent = e->getNbDofsOnEnt();
377 const auto brother_side = side_ptr->brother_side_number;
378 auto &dat = data.dataOnEntities[type][side];
379 auto &ent_field_indices = dat.getIndices();
380 auto &ent_field_local_indices = dat.getLocalIndices();
381
382 ent_field_indices.resize(nb_dofs_on_ent, false);
383 ent_field_local_indices.resize(nb_dofs_on_ent, false);
384 std::fill(ent_field_indices.data().begin(),
385 ent_field_indices.data().end(), -1);
386 std::fill(ent_field_local_indices.data().begin(),
387 ent_field_local_indices.data().end(), -1);
388
389 if (auto cache = extractor(e).lock()) {
390 for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
391 const int idx = (*dit)->getEntDofIdx();
392 ent_field_indices[idx] = (*dit)->getPetscGlobalDofIdx();
393 ent_field_local_indices[idx] = (*dit)->getPetscLocalDofIdx();
394 }
395 }
396
397 if (brother_side != -1) {
398 auto &dat_brother = data.dataOnEntities[type][brother_side];
399 dat_brother.getIndices().resize(nb_dofs_on_ent, false);
400 dat_brother.getLocalIndices().resize(nb_dofs_on_ent, false);
401 noalias(dat_brother.getIndices()) = dat.getIndices();
402 noalias(dat_brother.getLocalIndices()) = dat.getLocalIndices();
403 }
404 }
405 }
406 }
407 }
408
410}

◆ getEntityRowIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntityRowIndices ( EntitiesFieldData data,
const int  bit_number,
const EntityType  type_lo = MBVERTEX,
const EntityType  type_hi = MBPOLYHEDRON 
) const
protected

Definition at line 412 of file ForcesAndSourcesCore.cpp.

414 {
415
416 struct Extractor {
417 boost::weak_ptr<EntityCacheNumeredDofs>
418 operator()(boost::shared_ptr<FieldEntity> &e) {
419 return e->entityCacheRowDofs;
420 }
421 };
422
423 return getEntityIndices(data, bit_number, getRowFieldEnts(), type_lo, type_hi,
424 Extractor());
425}
auto & getRowFieldEnts() const

◆ getEntitySense() [1/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntitySense ( const EntityType  type,
boost::ptr_vector< EntitiesFieldData::EntData > &  data 
) const
protected

get sense (orientation) of entity

Parameters
typetype of entity
dataentity data
Returns
error code

Definition at line 84 of file ForcesAndSourcesCore.cpp.

86 {
88
89 auto &side_table = numeredEntFiniteElementPtr->getSideNumberTable().get<0>();
90 auto sit = side_table.lower_bound(get_id_for_min_type(type));
91 if (sit != side_table.end()) {
92 auto hi_sit = side_table.upper_bound(get_id_for_max_type(type));
93 for (; sit != hi_sit; ++sit) {
94 if (const auto side_number = (*sit)->side_number; side_number >= 0) {
95 const int brother_side_number = (*sit)->brother_side_number;
96 const int sense = (*sit)->sense;
97
98 data[side_number].getSense() = sense;
99 if (brother_side_number != -1)
100 data[brother_side_number].getSense() = sense;
101 }
102 }
103 }
105}

◆ getEntitySense() [2/2]

template<EntityType type>
MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntitySense ( EntitiesFieldData data) const
inlineprotected

Get the entity sense (orientation)

Template Parameters
type
Parameters
data
Returns
MoFEMErrorCode

Definition at line 200 of file ForcesAndSourcesCore.hpp.

200 {
201 return getEntitySense(type, data.dataOnEntities[type]);
202 }
MoFEMErrorCode getEntitySense(const EntityType type, boost::ptr_vector< EntitiesFieldData::EntData > &data) const
get sense (orientation) of entity

◆ getFaceNodes()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getFaceNodes ( EntitiesFieldData data) const
protected

Get nodes on faces.

Definition at line 921 of file ForcesAndSourcesCore.cpp.

921 {
923 auto &face_nodes = data.facesNodes;
924 auto &face_nodes_order = data.facesNodesOrder;
925 auto &side_table = const_cast<SideNumber_multiIndex &>(
926 numeredEntFiniteElementPtr->getSideNumberTable());
927 const auto ent = numeredEntFiniteElementPtr->getEnt();
928 const auto type = numeredEntFiniteElementPtr->getEntType();
929 const auto nb_faces = CN::NumSubEntities(type, 2);
930 const EntityHandle *conn_ele;
931 int num_nodes_ele;
932 CHKERR mField.get_moab().get_connectivity(ent, conn_ele, num_nodes_ele, true);
933 auto side_ptr_it = side_table.get<1>().lower_bound(
934 boost::make_tuple(CN::TypeDimensionMap[2].first, 0));
935 auto hi_side_ptr_it = side_table.get<1>().upper_bound(
936 boost::make_tuple(CN::TypeDimensionMap[2].second, 100));
937
938 for (; side_ptr_it != hi_side_ptr_it; ++side_ptr_it) {
939 const auto side = (*side_ptr_it)->side_number;
940 const auto sense = (*side_ptr_it)->sense;
941 const auto offset = (*side_ptr_it)->offset;
942
943 EntityType face_type;
944 int nb_nodes_face;
945 auto face_indices =
946 CN::SubEntityVertexIndices(type, 2, side, face_type, nb_nodes_face);
947 face_nodes.resize(nb_faces, nb_nodes_face);
948 face_nodes_order.resize(nb_faces, nb_nodes_face);
949
950 if (sense == 1)
951 for (int n = 0; n != nb_nodes_face; ++n)
952 face_nodes_order(side, n) = (n + offset) % nb_nodes_face;
953 else
954 for (int n = 0; n != nb_nodes_face; ++n)
955 face_nodes_order(side, n) =
956 (nb_nodes_face - (n - offset) % nb_nodes_face) % nb_nodes_face;
957
958 for (int n = 0; n != nb_nodes_face; ++n)
959 face_nodes(side, n) = face_indices[face_nodes_order(side, n)];
960
961#ifndef NDEBUG
962 const auto face_ent = (*side_ptr_it)->ent;
963 auto check = [&]() {
965 const EntityHandle *conn_face;
966 // int nb_nodes_face;
967 CHKERR mField.get_moab().get_connectivity(face_ent, conn_face,
968 nb_nodes_face, true);
969 face_nodes.resize(nb_faces, nb_nodes_face);
970 for (int nn = 0; nn != nb_nodes_face; ++nn) {
971 if (face_nodes(side, nn) !=
972 std::distance(
973 conn_ele,
974 std::find(conn_ele, &conn_ele[num_nodes_ele], conn_face[nn]))) {
975 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
976 "Wrong face numeration");
977 }
978 }
980 };
981 CHKERR check();
982#endif
983 }
984
986}
FTensor::Index< 'n', SPACE_DIM > n
multi_index_container< boost::shared_ptr< SideNumber >, indexed_by< ordered_unique< member< SideNumber, EntityHandle, &SideNumber::ent > >, ordered_non_unique< composite_key< SideNumber, const_mem_fun< SideNumber, EntityType, &SideNumber::getEntType >, member< SideNumber, signed char, &SideNumber::side_number > > > > > SideNumber_multiIndex
SideNumber_multiIndex for SideNumber.
virtual moab::Interface & get_moab()=0

◆ getMaxColOrder()

int MoFEM::ForcesAndSourcesCore::getMaxColOrder ( ) const

Get max order of approximation for field in columns.

Definition at line 135 of file ForcesAndSourcesCore.cpp.

135 {
137}
static int getMaxOrder(const ENTMULTIINDEX &multi_index)

◆ getMaxDataOrder()

int MoFEM::ForcesAndSourcesCore::getMaxDataOrder ( ) const

Get max order of approximation for data fields.

getMeasure getMaxDataOrder () return maximal order on entities, for all data on the element. So for example if finite element is triangle, and triangle base function have order 4 and on edges base function have order 2, this function return 4.

If finite element has for example 2 or more approximated fields, for example Pressure (order 3) and displacement field (order 5), this function returns 5.

Definition at line 120 of file ForcesAndSourcesCore.cpp.

120 {
121 int max_order = 0;
122 for (auto e : getDataFieldEnts()) {
123 if (auto ptr = e.lock()) {
124 const int order = ptr->getMaxOrder();
125 max_order = (max_order < order) ? order : max_order;
126 }
127 }
128 return max_order;
129}
constexpr int order

◆ getMaxRowOrder()

int MoFEM::ForcesAndSourcesCore::getMaxRowOrder ( ) const

Get max order of approximation for field in rows.

Definition at line 131 of file ForcesAndSourcesCore.cpp.

131 {
133}

◆ getNodesFieldData()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNodesFieldData ( EntitiesFieldData data,
const int  bit_number 
) const
protected

Get data on nodes.

Parameters
dataData structure
field_nameField name
Returns
Error code

Definition at line 642 of file ForcesAndSourcesCore.cpp.

643 {
644
645 auto get_nodes_field_data = [&](VectorDouble &nodes_data,
646 VectorFieldEntities &field_entities,
647 VectorDofs &nodes_dofs, FieldSpace &space,
649 VectorInt &bb_node_order) {
651
652 nodes_data.resize(0, false);
653 nodes_dofs.resize(0, false);
654 field_entities.resize(0, false);
655
656 auto field_it =
657 fieldsPtr->get<BitFieldId_mi_tag>().find(BitFEId().set(bit_number - 1));
658 if (field_it != fieldsPtr->get<BitFieldId_mi_tag>().end()) {
659
660#ifndef NDEBUG
661 if ((*field_it)->getBitNumber() != bit_number)
662 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong bit number");
663#endif
664 const int nb_dofs_on_vert = (*field_it)->getNbOfCoeffs();
665 space = (*field_it)->getSpace();
666 base = (*field_it)->getApproxBase();
667
668 auto &field_ents = getDataFieldEnts();
669 const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
670 bit_number, get_id_for_min_type<MBVERTEX>());
671 auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
672 cmp_uid_lo);
673 if (lo != field_ents.end()) {
674 const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
675 bit_number, get_id_for_max_type<MBVERTEX>());
676 auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
677 if (lo != hi) {
678
679 int nb_dofs = 0;
680 for (auto it = lo; it != hi; ++it) {
681 if (auto e = it->lock()) {
682 if (auto cache = e->entityCacheDataDofs.lock()) {
683 if (cache->loHi[0] != cache->loHi[1]) {
684 nb_dofs += std::distance(cache->loHi[0], cache->loHi[1]);
685 break;
686 }
687 }
688 }
689 }
690
691 if (nb_dofs) {
692
693 const int num_nodes = getNumberOfNodes();
694 bb_node_order.resize(num_nodes, false);
695 bb_node_order.clear();
696 const int max_nb_dofs = nb_dofs_on_vert * num_nodes;
697 nodes_data.resize(max_nb_dofs, false);
698 nodes_dofs.resize(max_nb_dofs, false);
699 field_entities.resize(num_nodes, false);
700 std::fill(nodes_data.begin(), nodes_data.end(), 0);
701 std::fill(nodes_dofs.begin(), nodes_dofs.end(), nullptr);
702 std::fill(field_entities.begin(), field_entities.end(), nullptr);
703
704 std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
705
706 for (auto it = lo; it != hi; ++it) {
707 if (auto e = it->lock()) {
708 const auto &sn = e->getSideNumberPtr();
709 // Some field entities on skeleton can have negative side
710 // number
711 if (const auto side_number = sn->side_number;
712 side_number >= 0) {
713 const int brother_side_number = sn->brother_side_number;
714
715 field_entities[side_number] = e.get();
716 if (brother_side_number != -1) {
717 brother_ents_vec.emplace_back(e);
718 field_entities[side_number] = field_entities[side_number];
719 }
720
721 bb_node_order[side_number] = e->getMaxOrder();
722 int pos = side_number * nb_dofs_on_vert;
723 auto ent_filed_data_vec = e->getEntFieldData();
724 if (auto cache = e->entityCacheDataDofs.lock()) {
725 for (auto dit = cache->loHi[0]; dit != cache->loHi[1];
726 ++dit) {
727 const auto dof_idx = (*dit)->getEntDofIdx();
728 nodes_data[pos + dof_idx] = ent_filed_data_vec[dof_idx];
729 nodes_dofs[pos + dof_idx] =
730 reinterpret_cast<FEDofEntity *>((*dit).get());
731 }
732 }
733 }
734 }
735 }
736
737 for (auto &it : brother_ents_vec) {
738 if (const auto e = it.lock()) {
739 const auto &sn = e->getSideNumberPtr();
740 const int side_number = sn->side_number;
741 const int brother_side_number = sn->brother_side_number;
742 bb_node_order[brother_side_number] = bb_node_order[side_number];
743 int pos = side_number * nb_dofs_on_vert;
744 int brother_pos = brother_side_number * nb_dofs_on_vert;
745 for (int ii = 0; ii != nb_dofs_on_vert; ++ii) {
746 nodes_data[brother_pos] = nodes_data[pos];
747 nodes_dofs[brother_pos] = nodes_dofs[pos];
748 ++pos;
749 ++brother_pos;
750 }
751 }
752 }
753 }
754 }
755 }
756 }
757
759 };
760
761 return get_nodes_field_data(
762 data.dataOnEntities[MBVERTEX][0].getFieldData(),
763 data.dataOnEntities[MBVERTEX][0].getFieldEntities(),
764 data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
765 data.dataOnEntities[MBVERTEX][0].getSpace(),
766 data.dataOnEntities[MBVERTEX][0].getBase(),
767 data.dataOnEntities[MBVERTEX][0].getBBNodeOrder());
768}
std::bitset< BITFEID_SIZE > BitFEId
Finite element Id.
Definition: Types.hpp:43
UBlasVector< double > VectorDouble
Definition: Types.hpp:68
UBlasVector< int > VectorInt
Definition: Types.hpp:67
ublas::vector< FieldEntity *, FieldEntAllocator > VectorFieldEntities
Field_multiIndex::index< FieldName_mi_tag >::type::iterator field_it
ublas::vector< FEDofEntity *, DofsAllocator > VectorDofs

◆ getNodesIndices()

template<typename EXTRACTOR >
MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNodesIndices ( const int  bit_number,
FieldEntity_vector_view ents_field,
VectorInt nodes_indices,
VectorInt local_nodes_indices,
EXTRACTOR &&  extractor 
) const
protected

get node indices

Definition at line 219 of file ForcesAndSourcesCore.cpp.

222 {
224
225 auto field_it = fieldsPtr->get<BitFieldId_mi_tag>().find(
226 BitFieldId().set(bit_number - 1));
227 if (field_it != fieldsPtr->get<BitFieldId_mi_tag>().end()) {
228
229#ifndef NDEBUG
230 if ((*field_it)->getBitNumber() != bit_number)
231 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Wrong bit number");
232#endif
233 const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
234 bit_number, get_id_for_min_type<MBVERTEX>());
235 auto lo = std::lower_bound(ents_field.begin(), ents_field.end(), lo_uid,
236 cmp_uid_lo);
237 if (lo != ents_field.end()) {
238 const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
239 bit_number, get_id_for_max_type<MBVERTEX>());
240 auto hi = std::upper_bound(lo, ents_field.end(), hi_uid, cmp_uid_hi);
241
242 const int num_nodes = getNumberOfNodes();
243 const int nb_dofs_on_vert = (*field_it)->getNbOfCoeffs();
244 const int max_nb_dofs = nb_dofs_on_vert * num_nodes;
245
246 int nb_dofs = 0;
247 for (auto it = lo; it != hi; ++it) {
248 if (auto e = it->lock()) {
249 if (auto cache = extractor(e).lock()) {
250 if (cache->loHi[0] != cache->loHi[1]) {
251 nb_dofs += std::distance(cache->loHi[0], cache->loHi[1]);
252 break;
253 }
254 }
255 }
256 }
257
258 if (nb_dofs) {
259 nodes_indices.resize(max_nb_dofs, false);
260 local_nodes_indices.resize(max_nb_dofs, false);
261 } else {
262 nodes_indices.resize(0, false);
263 local_nodes_indices.resize(0, false);
264 }
265
266 if (nb_dofs != max_nb_dofs) {
267 std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
268 std::fill(local_nodes_indices.begin(), local_nodes_indices.end(), -1);
269 }
270
271 for (auto it = lo; it != hi; ++it) {
272 if (auto e = it->lock()) {
273 auto side_ptr = e->getSideNumberPtr();
274 if (const auto side_number = side_ptr->side_number;
275 side_number >= 0) {
276 const auto brother_side_number = side_ptr->brother_side_number;
277 if (auto cache = extractor(e).lock()) {
278 for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
279 auto &dof = **dit;
280 const int idx = dof.getPetscGlobalDofIdx();
281 const int local_idx = dof.getPetscLocalDofIdx();
282 const int pos =
283 side_number * nb_dofs_on_vert + dof.getDofCoeffIdx();
284 nodes_indices[pos] = idx;
285 local_nodes_indices[pos] = local_idx;
286 if (brother_side_number != -1) {
287 const int elem_idx = brother_side_number * nb_dofs_on_vert +
288 (*dit)->getDofCoeffIdx();
289 nodes_indices[elem_idx] = idx;
290 local_nodes_indices[elem_idx] = local_idx;
291 }
292 }
293 }
294 }
295 }
296 }
297 } else {
298 nodes_indices.resize(0, false);
299 local_nodes_indices.resize(0, false);
300 }
301
302 } else {
303 nodes_indices.resize(0, false);
304 local_nodes_indices.resize(0, false);
305 }
306
308}
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
Definition: Types.hpp:42

◆ getNoFieldColIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldColIndices ( EntitiesFieldData data,
const int  bit_number 
) const
protected

get col NoField indices

Definition at line 474 of file ForcesAndSourcesCore.cpp.

475 {
477#ifndef NDEBUG
478 if (data.dataOnEntities[MBENTITYSET].size() == 0) {
479 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
480 "data.dataOnEntities[MBENTITYSET] is empty");
481 }
482#endif
484 data.dataOnEntities[MBENTITYSET][0].getIndices());
486}
auto getColDofsPtr() const
MoFEMErrorCode getNoFieldIndices(const int bit_number, boost::shared_ptr< FENumeredDofEntity_multiIndex > dofs, VectorInt &nodes_indices) const
get NoField indices

◆ getNoFieldFieldData() [1/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldFieldData ( const int  bit_number,
VectorDouble ent_field_data,
VectorDofs ent_field_dofs,
VectorFieldEntities ent_field 
) const
protected

Get field data on nodes.

Definition at line 855 of file ForcesAndSourcesCore.cpp.

857 {
858
860
861 ent_field_data.resize(0, false);
862 ent_field_dofs.resize(0, false);
863 ent_field.resize(0, false);
864
865 auto &field_ents = getDataFieldEnts();
866 const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
867 bit_number, get_id_for_min_type<MBVERTEX>());
868 auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
869 cmp_uid_lo);
870 if (lo != field_ents.end()) {
871
872 ent_field.resize(field_ents.size(), false);
873 std::fill(ent_field.begin(), ent_field.end(), nullptr);
874
875 const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
876 bit_number, get_id_for_max_type<MBENTITYSET>());
877 auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
878 if (lo != hi) {
879
880 int side = 0;
881 for (auto it = lo; it != hi; ++it, ++side)
882 if (auto e = it->lock()) {
883
884 const auto size = e->getNbDofsOnEnt();
885 ent_field_data.resize(size, false);
886 ent_field_dofs.resize(size, false);
887 ent_field[side] = e.get();
888 noalias(ent_field_data) = e->getEntFieldData();
889
890 if (auto cache = e->entityCacheDataDofs.lock()) {
891 for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
892 ent_field_dofs[(*dit)->getEntDofIdx()] =
893 reinterpret_cast<FEDofEntity *>((*dit).get());
894 }
895 }
896 }
897 }
898 }
899
901}

◆ getNoFieldFieldData() [2/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldFieldData ( EntitiesFieldData data,
const int  bit_number 
) const
protected

Definition at line 904 of file ForcesAndSourcesCore.cpp.

905 {
907 if (data.dataOnEntities[MBENTITYSET].size() == 0)
908 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
909 "No space to insert data");
910
912 bit_number, data.dataOnEntities[MBENTITYSET][0].getFieldData(),
913 data.dataOnEntities[MBENTITYSET][0].getFieldDofs(),
914 data.dataOnEntities[MBENTITYSET][0].getFieldEntities());
916}
MoFEMErrorCode getNoFieldFieldData(const int bit_number, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs, VectorFieldEntities &ent_field) const
Get field data on nodes.

◆ getNoFieldIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldIndices ( const int  bit_number,
boost::shared_ptr< FENumeredDofEntity_multiIndex dofs,
VectorInt nodes_indices 
) const
protected

get NoField indices

Definition at line 442 of file ForcesAndSourcesCore.cpp.

444 {
446 auto dit = dofs->get<Unique_mi_tag>().lower_bound(
448 auto hi_dit = dofs->get<Unique_mi_tag>().upper_bound(
450 indices.resize(std::distance(dit, hi_dit));
451 for (; dit != hi_dit; dit++) {
452 int idx = (*dit)->getPetscGlobalDofIdx();
453 indices[(*dit)->getDofCoeffIdx()] = idx;
454 }
456}
static UId getHiBitNumberUId(const FieldBitNumber bit_number)
static UId getLoBitNumberUId(const FieldBitNumber bit_number)

◆ getNoFieldRowIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldRowIndices ( EntitiesFieldData data,
const int  bit_number 
) const
protected

get col NoField indices

Definition at line 459 of file ForcesAndSourcesCore.cpp.

460 {
462#ifndef NDEBUG
463 if (data.dataOnEntities[MBENTITYSET].size() == 0) {
464 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
465 "data.dataOnEntities[MBENTITYSET] is empty");
466 }
467#endif
469 data.dataOnEntities[MBENTITYSET][0].getIndices());
471}
auto getRowDofsPtr() const

◆ getOpPtrVector()

boost::ptr_deque< UserDataOperator > & MoFEM::ForcesAndSourcesCore::getOpPtrVector ( )
inline

Use to push back operator for row operator.

It can be used to calculate nodal forces or other quantities on the mesh.

Examples
MagneticElement.hpp, continuity_check_on_contact_prism_side_ele.cpp, continuity_check_on_skeleton_3d.cpp, continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, forces_and_sources_testing_flat_prism_element.cpp, and lorentz_force.cpp.

Definition at line 83 of file ForcesAndSourcesCore.hpp.

83{ return opPtrVector; }
boost::ptr_deque< UserDataOperator > opPtrVector
Vector of finite element users data operators.

◆ getProblemNodesColIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getProblemNodesColIndices ( const std::string &  field_name,
VectorInt nodes_indices 
) const
protected

Definition at line 579 of file ForcesAndSourcesCore.cpp.

580 {
582 nodes_indices);
583}
const Problem * problemPtr
raw pointer to problem
MoFEMErrorCode getProblemNodesIndices(const std::string &field_name, const NumeredDofEntity_multiIndex &dofs, VectorInt &nodes_indices) const
get indices of nodal indices which are declared for problem but not this particular element
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredColDofsPtr
store DOFs on columns for this problem

◆ getProblemNodesIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getProblemNodesIndices ( const std::string &  field_name,
const NumeredDofEntity_multiIndex dofs,
VectorInt nodes_indices 
) const
protected

get indices of nodal indices which are declared for problem but not this particular element

Definition at line 490 of file ForcesAndSourcesCore.cpp.

492 {
494
495 const Field *field_struture = mField.get_field_structure(field_name);
496 if (field_struture->getSpace() == H1) {
497
498 const int num_nodes = getNumberOfNodes();
499 nodes_indices.resize(field_struture->getNbOfCoeffs() * num_nodes, false);
500 std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
501
502 auto &side_table = const_cast<SideNumber_multiIndex &>(
503 numeredEntFiniteElementPtr->getSideNumberTable());
504 auto siit = side_table.get<1>().lower_bound(boost::make_tuple(MBVERTEX, 0));
505 auto hi_siit =
506 side_table.get<1>().lower_bound(boost::make_tuple(MBVERTEX, 10000));
507
508 int nn = 0;
509 for (; siit != hi_siit; siit++, nn++) {
510 if (siit->get()->side_number >= 0) {
511 auto bit_number = mField.get_field_bit_number(field_name);
512 const EntityHandle ent = siit->get()->ent;
513 auto dit = dofs.get<Unique_mi_tag>().lower_bound(
515 auto hi_dit = dofs.get<Unique_mi_tag>().upper_bound(
517 for (; dit != hi_dit; dit++) {
518 nodes_indices[siit->get()->side_number * (*dit)->getNbOfCoeffs() +
519 (*dit)->getDofCoeffIdx()] =
520 (*dit)->getPetscGlobalDofIdx();
521 }
522 }
523 }
524 } else {
525 nodes_indices.resize(0, false);
526 }
527
529}
virtual const Field * get_field_structure(const std::string &name, enum MoFEMTypes bh=MF_EXIST) const =0
get field structure
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
static UId getLoLocalEntityBitNumber(const char bit_number, const EntityHandle ent)
static UId getHiLocalEntityBitNumber(const char bit_number, const EntityHandle ent)

◆ getProblemNodesRowIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getProblemNodesRowIndices ( const std::string &  field_name,
VectorInt nodes_indices 
) const
protected

Definition at line 565 of file ForcesAndSourcesCore.cpp.

566 {
568 nodes_indices);
569}
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredRowDofsPtr
store DOFs on rows for this problem

◆ getProblemTypeColIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getProblemTypeColIndices ( const std::string &  field_name,
EntityType  type,
int  side_number,
VectorInt indices 
) const
protected

Definition at line 586 of file ForcesAndSourcesCore.cpp.

588 {
590 type, side_number, indices);
591}
MoFEMErrorCode getProblemTypeIndices(const std::string &field_name, const NumeredDofEntity_multiIndex &dofs, EntityType type, int side_number, VectorInt &indices) const
get indices by type (generic function) which are declared for problem but not this particular element

◆ getProblemTypeIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getProblemTypeIndices ( const std::string &  field_name,
const NumeredDofEntity_multiIndex dofs,
EntityType  type,
int  side_number,
VectorInt indices 
) const
protected

get indices by type (generic function) which are declared for problem but not this particular element

Definition at line 531 of file ForcesAndSourcesCore.cpp.

533 {
535
536 indices.resize(0);
537
538 auto &side_table = const_cast<SideNumber_multiIndex &>(
539 numeredEntFiniteElementPtr->getSideNumberTable());
540 auto siit =
541 side_table.get<1>().lower_bound(boost::make_tuple(type, side_number));
542 auto hi_siit =
543 side_table.get<1>().upper_bound(boost::make_tuple(type, side_number));
544
545 for (; siit != hi_siit; siit++) {
546 if (siit->get()->side_number >= 0) {
547
548 const EntityHandle ent = siit->get()->ent;
549 auto bit_number = mField.get_field_bit_number(field_name);
550 auto dit = dofs.get<Unique_mi_tag>().lower_bound(
552 auto hi_dit = dofs.get<Unique_mi_tag>().upper_bound(
554 indices.resize(std::distance(dit, hi_dit));
555 for (; dit != hi_dit; dit++) {
556
557 indices[(*dit)->getEntDofIdx()] = (*dit)->getPetscGlobalDofIdx();
558 }
559 }
560 }
561
563}

◆ getProblemTypeRowIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getProblemTypeRowIndices ( const std::string &  field_name,
EntityType  type,
int  side_number,
VectorInt indices 
) const
protected

Definition at line 572 of file ForcesAndSourcesCore.cpp.

574 {
576 type, side_number, indices);
577}

◆ getRowNodesIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getRowNodesIndices ( EntitiesFieldData data,
const int  bit_number 
) const
protected

get row node indices from FENumeredDofEntity_multiIndex

Definition at line 311 of file ForcesAndSourcesCore.cpp.

312 {
313
314 struct Extractor {
315 boost::weak_ptr<EntityCacheNumeredDofs>
316 operator()(boost::shared_ptr<FieldEntity> &e) {
317 return e->entityCacheRowDofs;
318 }
319 };
320
321 return getNodesIndices(bit_number, getRowFieldEnts(),
322 data.dataOnEntities[MBVERTEX][0].getIndices(),
323 data.dataOnEntities[MBVERTEX][0].getLocalIndices(),
324 Extractor());
325}

◆ getRule() [1/2]

int MoFEM::ForcesAndSourcesCore::getRule ( int  order)
protectedvirtual
Deprecated:
Use getRule(int row_order, int col_order, int data order)

Reimplemented in MoFEM::EdgeElementForcesAndSourcesCoreOnChildParent, MoFEM::FaceElementForcesAndSourcesCoreOnChildParent, MoFEM::FaceElementForcesAndSourcesCoreOnSide, MoFEM::VolumeElementForcesAndSourcesCoreOnContactPrismSide, MoFEM::VolumeElementForcesAndSourcesCoreOnSide, MoFEM::PostProcBrokenMeshInMoabBase< EdgeElementForcesAndSourcesCore >, MoFEM::PostProcBrokenMeshInMoabBase< FaceElementForcesAndSourcesCore >, MoFEM::PostProcBrokenMeshInMoabBase< VolumeElementForcesAndSourcesCore >, CohesiveElement::CohesiveInterfaceElement::MyPrism, AnalyticalDirichletBC::ApproxField::MyTriFE, BodyForceConstantField::MyVolumeFE, EdgeForce::MyFE, FluidPressure::MyTriangleFE, KelvinVoigtDamper::DamperFE, NonlinearElasticElement::MyVolumeFE, PostProcFaceOnRefinedMesh, PostProcEdgeOnRefinedMesh, NeumannForcesSurface::MyTriangleFE, NeumannForcesSurfaceComplexForLazy::MyTriangleSpatialFE, SurfaceSlidingConstrains::MyTriangleFE, EdgeSlidingConstrains::MyEdgeFE, ThermalElement::MyVolumeFE, ThermalElement::MyTriFE, ThermalStressElement::MyVolumeFE, BoneRemodeling::DensityMapFe, BoneRemodeling::Remodeling::Fe, CellEngineering::FaceElement, CellEngineering::DispMapFe, GelModule::Gel::GelFE, NitscheMethod::MyFace, NitscheMethod::MyVolumeFE, SmallStrainPlasticity::MyVolumeFE, MixTransport::MixTransportElement::MyVolumeFE, MixTransport::MixTransportElement::MyTriFE, MagneticElement::VolumeFE, and MagneticElement::TriFE.

Definition at line 1932 of file ForcesAndSourcesCore.cpp.

1932{ return 2 * order; }

◆ getRule() [2/2]

int MoFEM::ForcesAndSourcesCore::getRule ( int  order_row,
int  order_col,
int  order_data 
)
protectedvirtual

another variant of getRule

Parameters
order_roworder of base function on row
order_colorder of base function on columns
order_dataorder of base function approximating data
Returns
integration rule

This function is overloaded by the user. The integration rule is set such that specific operator implemented by the user is integrated accurately. For example if user implement bilinear operator

\[ b(u,v) = \int_\mathcal{T} \frac{\partial u_i}{\partial x_j}\frac{\partial v_i}{\partial x_j} \textrm{d}\mathcal{T} \]

then if \(u\) and \(v\) are polynomial of given order, then exact integral would be

int getRule(int order) { return
2*(order-1); };
virtual int getRule(int order_row, int order_col, int order_data)
another variant of getRule

The integration points and weights are set appropriately for given entity type and integration rule from quad.c

Method ForcesAndSourcesCore::getRule takes at argument takes maximal polynomial order set on the element on all fields defined on the element. If a user likes to have more control, another variant of this function can be called which distinguishing between field orders on rows, columns and data, the i.e. first argument of a bilinear form, the second argument of bilinear form and field coefficients on the element.

Note
If user set rule to -1 or any other negative integer, then method ForcesAndSourcesCore::setGaussPts is called. In that method user can implement own (specific) integration method.
Bug:
this function should be const

Reimplemented in TriFE, QuadFE, and EdgeFE.

Definition at line 1920 of file ForcesAndSourcesCore.cpp.

1921 {
1922 return getRuleHook ? getRuleHook(order_row, order_col, order_data)
1923 : getRule(order_data);
1924}

◆ getSpacesAndBaseOnEntities()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getSpacesAndBaseOnEntities ( EntitiesFieldData data) const
protected

Get field approximation space and base on entities.

Definition at line 990 of file ForcesAndSourcesCore.cpp.

991 {
993
994 if (nInTheLoop == 0) {
995 data.sPace.reset();
996 data.bAse.reset();
997 for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
998 data.spacesOnEntities[t].reset();
999 data.basesOnEntities[t].reset();
1000 }
1001 for (int s = 0; s != LASTSPACE; ++s) {
1002 data.basesOnSpaces[s].reset();
1003 }
1004 }
1005
1006 auto fe_type = numeredEntFiniteElementPtr->getEntType();
1007
1008 if (getDataFieldEntsPtr())
1009 for (auto e : getDataFieldEnts()) {
1010 if (auto ptr = e.lock()) {
1011 // get data from entity
1012 const EntityType type = ptr->getEntType();
1013 if (DefaultElementAdjacency::getDefTypeMap(fe_type, type)) {
1014 const FieldSpace space = ptr->getSpace();
1015 const FieldApproximationBase approx = ptr->getApproxBase();
1016 // set data
1017 data.sPace.set(space);
1018 data.bAse.set(approx);
1019 data.spacesOnEntities[type].set(space);
1020 data.basesOnEntities[type].set(approx);
1021 data.basesOnSpaces[space].set(approx);
1022 }
1023 }
1024 }
1025 else
1026 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1027 "data fields ents not allocated on element");
1028
1030}
int nInTheLoop
number currently of processed method
static bool getDefTypeMap(const EntityType fe_type, const EntityType ent_type)
boost::shared_ptr< FieldEntity_vector_view > & getDataFieldEntsPtr() const

◆ getUserPolynomialBase()

auto & MoFEM::ForcesAndSourcesCore::getUserPolynomialBase ( )
inline

Get the User Polynomial Base object.

Returns
boost::shared_ptr<BaseFunction>&

Definition at line 97 of file ForcesAndSourcesCore.hpp.

97{ return userPolynomialBasePtr; }
boost::shared_ptr< BaseFunction > userPolynomialBasePtr
Pointer to user polynomail base.

◆ loopOverOperators()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::loopOverOperators ( )
protected

Iterate user data operators.

Returns
MoFEMErrorCode

Definition at line 1347 of file ForcesAndSourcesCore.cpp.

1347 {
1349
1350 using UDO = UserDataOperator;
1351
1352 std::array<std::string, 2> field_name;
1353 std::array<const Field *, 3> field_struture;
1354 std::array<int, 2>
1355 field_id; // note the this is field bit number (nor field bit)
1356 std::array<FieldSpace, 2> space;
1357 std::array<FieldApproximationBase, 2> base;
1358
1359 constexpr std::array<UDO::OpType, 2> types = {UDO::OPROW, UDO::OPCOL};
1360 std::array<int, 2> last_eval_field_id = {0, 0};
1361
1362 std::array<boost::shared_ptr<EntitiesFieldData>, 2> op_data;
1363
1364 auto swap_bases = [&](auto &op) {
1366 for (size_t ss = 0; ss != 2; ++ss) {
1367 if (op.getOpType() & types[ss] || op.getOpType() & UDO::OPROWCOL) {
1368 switch (base[ss]) {
1370 CHKERR op_data[ss]->baseSwap(field_name[ss],
1372 default:
1373 break;
1374 }
1375 }
1376 }
1377
1379 };
1380
1381 const EntityType type = numeredEntFiniteElementPtr->getEntType();
1382
1383 // evaluate entity data only, no field specific data provided or known
1384 auto evaluate_op_space = [&](auto &op) {
1386
1387 // reseat all data which all field dependent
1388 dataOnElement[op.sPace]->resetFieldDependentData();
1389 std::fill(last_eval_field_id.begin(), last_eval_field_id.end(), 0);
1390
1391 switch (op.sPace) {
1392 case NOSPACE:
1393 try {
1394 CHKERR op.doWork(
1395 0, MBENTITYSET,
1396 dataOnElement[NOSPACE]->dataOnEntities[MBENTITYSET][0]);
1397 }
1398 CATCH_OP_ERRORS(op);
1399 break;
1400 case NOFIELD:
1401 case H1:
1402 case HCURL:
1403 case HDIV:
1404 case L2:
1405 try {
1406
1407 for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
1408 for (auto &e : dataOnElement[op.sPace]->dataOnEntities[t]) {
1409 e.getSpace() = op.sPace;
1410 }
1411 }
1412
1413 CHKERR op.opRhs(*dataOnElement[op.sPace], false);
1414 }
1415 CATCH_OP_ERRORS(op);
1416 break;
1417 default:
1418 SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1419 "Not implemented for this space", op.sPace);
1420 }
1421
1423 };
1424
1425 // set entity data
1426 auto set_op_entityties_data = [&](auto ss, auto &op) {
1428
1429#ifndef NDEBUG
1430 if ((op.getNumeredEntFiniteElementPtr()->getBitFieldIdData() &
1432 .none()) {
1433 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1434 "no data field < %s > on finite element < %s >",
1435 field_name[ss].c_str(), getFEName().c_str());
1436 }
1437#endif
1438
1439 op_data[ss] =
1440 !ss ? dataOnElement[space[ss]] : derivedDataOnElement[space[ss]];
1441
1442 for (auto &data : op_data[ss]->dataOnEntities[MBENTITYSET]) {
1443 CHKERR data.resetFieldDependentData();
1444 }
1445
1446 auto get_data_for_nodes = [&]() {
1448 if (!ss)
1449 CHKERR getRowNodesIndices(*op_data[ss], field_id[ss]);
1450 else
1451 CHKERR getColNodesIndices(*op_data[ss], field_id[ss]);
1452 CHKERR getNodesFieldData(*op_data[ss], field_id[ss]);
1454 };
1455
1456 // get data on entities but not nodes
1457 auto get_data_for_entities = [&]() {
1459 CHKERR getEntityFieldData(*op_data[ss], field_id[ss], MBEDGE);
1460 if (!ss)
1461 CHKERR getEntityRowIndices(*op_data[ss], field_id[ss], MBEDGE);
1462 else
1463 CHKERR getEntityColIndices(*op_data[ss], field_id[ss], MBEDGE);
1465 };
1466
1467 auto get_data_for_meshset = [&]() {
1469 if (!ss) {
1470 CHKERR getNoFieldRowIndices(*op_data[ss], field_id[ss]);
1471 } else {
1472 CHKERR getNoFieldColIndices(*op_data[ss], field_id[ss]);
1473 }
1474 CHKERR getNoFieldFieldData(*op_data[ss], field_id[ss]);
1476 };
1477
1478 switch (space[ss]) {
1479 case H1:
1480 CHKERR get_data_for_nodes();
1481 case HCURL:
1482 case HDIV:
1483 CHKERR get_data_for_entities();
1484 break;
1485 case L2:
1486 switch (type) {
1487 case MBVERTEX:
1488 CHKERR get_data_for_nodes();
1489 break;
1490 default:
1491 CHKERR get_data_for_entities();
1492 }
1493 break;
1494 case NOFIELD:
1495 // if (!getNinTheLoop()) {
1496 // NOFIELD data are the same for each element, can be
1497 // retrieved only once
1498 CHKERR get_data_for_meshset();
1499 // }
1500 break;
1501 default:
1502 SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1503 "not implemented for this space < %s >",
1504 FieldSpaceNames[space[ss]]);
1505 }
1507 };
1508
1509 // evalate operators with field data
1510 auto evaluate_op_for_fields = [&](auto &op) {
1512
1513 if (op.getOpType() & UDO::OPROW) {
1514 try {
1515 CHKERR op.opRhs(*op_data[0], false);
1516 }
1517 CATCH_OP_ERRORS(op);
1518 }
1519
1520 if (op.getOpType() & UDO::OPCOL) {
1521 try {
1522 CHKERR op.opRhs(*op_data[1], false);
1523 }
1524 CATCH_OP_ERRORS(op);
1525 }
1526
1527 if (op.getOpType() & UDO::OPROWCOL) {
1528 try {
1529 CHKERR op.opLhs(*op_data[0], *op_data[1]);
1530 }
1531 CATCH_OP_ERRORS(op);
1532 }
1534 };
1535
1536 // Collect bit ref level on all entities, fields and spaces
1538
1539 auto oit = opPtrVector.begin();
1540 auto hi_oit = opPtrVector.end();
1541
1542 // interate over all operators on element
1543 for (; oit != hi_oit; oit++) {
1544
1545 try {
1546
1547 CHKERR oit->setPtrFE(this);
1548
1549 if ((oit->opType & UDO::OPSPACE) == UDO::OPSPACE) {
1550
1551 // run operator for space but specific field
1552 CHKERR evaluate_op_space(*oit);
1553
1554 } else if (
1555
1556 (oit->opType & (UDO::OPROW | UDO::OPCOL | UDO::OPROWCOL)) ==
1557 oit->opType
1558
1559 ) {
1560
1561 field_name[0] = oit->rowFieldName;
1562 field_name[1] = oit->colFieldName;
1563
1564 for (size_t ss = 0; ss != 2; ss++) {
1565
1566#ifndef NDEBUG
1567 if (field_name[ss].empty()) {
1568 SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1569 "Not set Field name in operator %d (0-row, 1-column) in "
1570 "operator %s",
1571 ss,
1572 (boost::typeindex::type_id_runtime(*oit).pretty_name())
1573 .c_str());
1574 }
1575#endif
1576
1577 field_struture[ss] = mField.get_field_structure(field_name[ss]);
1578 field_id[ss] = field_struture[ss]->getBitNumber();
1579 space[ss] = field_struture[ss]->getSpace();
1580 base[ss] = field_struture[ss]->getApproxBase();
1581 }
1582
1583 // not that if field name do not change between operators, entity field
1584 // data are nor rebuild
1585 for (size_t ss = 0; ss != 2; ss++) {
1586
1587 if (oit->getOpType() & types[ss] ||
1588 oit->getOpType() & UDO::OPROWCOL) {
1589 if (last_eval_field_id[ss] != field_id[ss]) {
1590 CHKERR set_op_entityties_data(ss, *oit);
1591 last_eval_field_id[ss] = field_id[ss];
1592 }
1593 }
1594 }
1595
1596 CHKERR swap_bases(*oit);
1597
1598 // run operator for given field or row, column or both
1599 CHKERR evaluate_op_for_fields(*oit);
1600
1601 CHKERR swap_bases(*oit);
1602
1603 } else {
1604
1605 for (int i = 0; i != 5; ++i)
1606 if (oit->opType & (1 << i))
1607 MOFEM_LOG("SELF", Sev::error) << UDO::OpTypeNames[i];
1608 SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1609 "Impossible operator type");
1610 }
1611 }
1612 CATCH_OP_ERRORS(*oit);
1613 }
1615}
#define CATCH_OP_ERRORS(OP)
ForcesAndSourcesCore::UserDataOperator UserDataOperator
static const char *const FieldSpaceNames[]
Definition: definitions.h:92
@ MOFEM_NOT_IMPLEMENTED
Definition: definitions.h:32
#define MOFEM_LOG(channel, severity)
Log.
Definition: LogManager.hpp:308
FTensor::Index< 'i', SPACE_DIM > i
virtual BitFieldId get_field_id(const std::string &name) const =0
Get field Id.
auto getFEName() const
get finite element name
FieldBitNumber getBitNumber() const
Get number of set bit in Field ID. Each field has uid, get getBitNumber get number of bit set for giv...
MoFEMErrorCode getEntityRowIndices(EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
MoFEMErrorCode getRowNodesIndices(EntitiesFieldData &data, const int bit_number) const
get row node indices from FENumeredDofEntity_multiIndex
MoFEMErrorCode getNoFieldRowIndices(EntitiesFieldData &data, const int bit_number) const
get col NoField indices
MoFEMErrorCode getEntityFieldData(EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
MoFEMErrorCode getColNodesIndices(EntitiesFieldData &data, const int bit_number) const
get col node indices from FENumeredDofEntity_multiIndex
MoFEMErrorCode getNodesFieldData(EntitiesFieldData &data, const int bit_number) const
Get data on nodes.
MoFEMErrorCode getEntityColIndices(EntitiesFieldData &data, const int bit_number, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
MoFEMErrorCode getNoFieldColIndices(EntitiesFieldData &data, const int bit_number) const
get col NoField indices

◆ operator()()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::operator() ( )
virtual

function is run for every finite element

It is used to calculate element local matrices and assembly. It can be used for post-processing.

Reimplemented from MoFEM::BasicMethod.

Reimplemented in MoFEM::ContactPrismElementForcesAndSourcesCore, MoFEM::EdgeElementForcesAndSourcesCore, MoFEM::FaceElementForcesAndSourcesCore, MoFEM::FatPrismElementForcesAndSourcesCore, MoFEM::FlatPrismElementForcesAndSourcesCore, MoFEM::VertexElementForcesAndSourcesCore, MoFEM::VolumeElementForcesAndSourcesCore, and NeumannForcesSurfaceComplexForLazy::MyTriangleSpatialFE.

Definition at line 1966 of file ForcesAndSourcesCore.cpp.

1966 {
1968 if (operatorHook) {
1969 ierr = operatorHook();
1970 CHKERRG(ierr);
1971 } else {
1972#ifndef NDEBUG
1973 MOFEM_LOG("SELF", Sev::warning)
1974 << "No method operator() overloaded on element entity on finite "
1975 "element <"
1976 << boost::typeindex::type_id_runtime(*this).pretty_name() << ">";
1977#endif
1978 }
1980}
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:483
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:76
boost::function< MoFEMErrorCode()> operatorHook
Hook function for operator.

◆ postProcess()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::postProcess ( )
virtual

function is run at the end of loop

It is used to assembly matrices and vectors, calculating global variables, f.e. total internal energy, ect.

Iterating over dofs: Example1 iterating over dofs in row by name of the field for(IT_GET_FEROW_BY_NAME_DOFS_FOR_LOOP(this,"DISPLACEMENT",it)) { ... }

Reimplemented from MoFEM::BasicMethod.

Reimplemented in MoFEM::PostProcBrokenMeshInMoabBase< EdgeElementForcesAndSourcesCore >, MoFEM::PostProcBrokenMeshInMoabBase< FaceElementForcesAndSourcesCore >, MoFEM::PostProcBrokenMeshInMoabBase< VolumeElementForcesAndSourcesCore >, KelvinVoigtDamper::DamperFE, NonlinearElasticElement::MyVolumeFE, PostProcFatPrismOnRefinedMesh, PostProcFaceOnRefinedMesh, PostProcEdgeOnRefinedMesh, Smoother::MyVolumeFE, GelModule::Gel::GelFE, SolidShellModule::SolidShellPrismElement::SolidShell, SolidShellModule::SolidShellPrismElement::SolidShellError, and SolidShellModule::SolidShellPrismElement::PostProcFatPrismOnTriangleOnRefinedMesh.

Definition at line 1981 of file ForcesAndSourcesCore.cpp.

1981 {
1983 if (postProcessHook) {
1985 CHKERRG(ierr);
1986 }
1988}
boost::function< MoFEMErrorCode()> postProcessHook
Hook function for post-processing.

◆ preProcess()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::preProcess ( )
virtual

◆ setGaussPts() [1/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::setGaussPts ( int  order)
protectedvirtual

◆ setGaussPts() [2/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::setGaussPts ( int  order_row,
int  order_col,
int  order_data 
)
protectedvirtual

set user specific integration rule

This function allows for user defined integration rule. The key is to called matrix gaussPts, which is used by other MoFEM procedures. Matrix has number of rows equal to problem dimension plus one, where last index is used to store weight values. Number of columns is equal to number of integration points.

Note
This function is called if method ForcesAndSourcesCore::getRule is returning integer -1 or any other negative integer.

User sets

where

gaussPts.resize(dim+1,nb_gauss_pts);
const int dim

number rows represents local coordinates of integration points in reference element, where last index in row is for integration weight.

Reimplemented in TriFE, QuadFE, and EdgeFE.

Definition at line 1926 of file ForcesAndSourcesCore.cpp.

1927 {
1928 return setRuleHook ? setRuleHook(this, order_row, order_col, order_data)
1929 : setGaussPts(order_data);
1930}
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule

◆ setRefineFEPtr()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::setRefineFEPtr ( const ForcesAndSourcesCore refine_fe_ptr)
private

Set the pointer to face element refined.

Parameters
refine_fe_ptr
Returns
MoFEMErrorCode

Definition at line 1661 of file ForcesAndSourcesCore.cpp.

1662 {
1664 refinePtrFE = const_cast<ForcesAndSourcesCore *>(refine_fe_ptr);
1666}
ForcesAndSourcesCore(Interface &m_field)

◆ setSideFEPtr()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::setSideFEPtr ( const ForcesAndSourcesCore side_fe_ptr)
private

Set the pointer to face element on the side.

Note
Function is is used by face element, while it iterates over elements on the side
Parameters
side_fe_ptr
Returns
MoFEMErrorCode

Definition at line 1655 of file ForcesAndSourcesCore.cpp.

1655 {
1657 sidePtrFE = const_cast<ForcesAndSourcesCore *>(side_fe_ptr);
1659}

Friends And Related Function Documentation

◆ EdgeElementForcesAndSourcesCoreOnChildParent

Definition at line 537 of file ForcesAndSourcesCore.hpp.

◆ FaceElementForcesAndSourcesCoreOnChildParent

Definition at line 536 of file ForcesAndSourcesCore.hpp.

◆ FaceElementForcesAndSourcesCoreOnSide

Definition at line 535 of file ForcesAndSourcesCore.hpp.

◆ UserDataOperator

friend class UserDataOperator
friend

Definition at line 482 of file ForcesAndSourcesCore.hpp.

◆ VolumeElementForcesAndSourcesCoreOnContactPrismSide

Definition at line 538 of file ForcesAndSourcesCore.hpp.

◆ VolumeElementForcesAndSourcesCoreOnSide

Definition at line 534 of file ForcesAndSourcesCore.hpp.

Member Data Documentation

◆ coordsAtGaussPts

MatrixDouble MoFEM::ForcesAndSourcesCore::coordsAtGaussPts
protected

coordinated at gauss points

Definition at line 541 of file ForcesAndSourcesCore.hpp.

◆ dataH1

EntitiesFieldData& MoFEM::ForcesAndSourcesCore::dataH1
protected

Definition at line 471 of file ForcesAndSourcesCore.hpp.

◆ dataHcurl

EntitiesFieldData& MoFEM::ForcesAndSourcesCore::dataHcurl
protected

Definition at line 472 of file ForcesAndSourcesCore.hpp.

◆ dataHdiv

EntitiesFieldData& MoFEM::ForcesAndSourcesCore::dataHdiv
protected

Definition at line 473 of file ForcesAndSourcesCore.hpp.

◆ dataL2

EntitiesFieldData& MoFEM::ForcesAndSourcesCore::dataL2
protected

Definition at line 474 of file ForcesAndSourcesCore.hpp.

◆ dataNoField

EntitiesFieldData& MoFEM::ForcesAndSourcesCore::dataNoField
protected

Definition at line 470 of file ForcesAndSourcesCore.hpp.

◆ dataOnElement

const std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE> MoFEM::ForcesAndSourcesCore::dataOnElement
protected

Entity data on element entity rows fields.

Definition at line 461 of file ForcesAndSourcesCore.hpp.

◆ derivedDataOnElement

const std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE> MoFEM::ForcesAndSourcesCore::derivedDataOnElement
protected

Entity data on element entity columns fields.

Definition at line 468 of file ForcesAndSourcesCore.hpp.

◆ elementMeasure

double MoFEM::ForcesAndSourcesCore::elementMeasure
protected

Depending on dimension of elements, stores length, area or volume of element.

Definition at line 542 of file ForcesAndSourcesCore.hpp.

◆ elementPolynomialBasePtr

boost::shared_ptr<BaseFunction> MoFEM::ForcesAndSourcesCore::elementPolynomialBasePtr
private

Pointer to entity polynomial base.

Definition at line 496 of file ForcesAndSourcesCore.hpp.

◆ gaussPts

MatrixDouble MoFEM::ForcesAndSourcesCore::gaussPts

Matrix of integration points.

Columns is equal to number of integration points, numver of rows depends on dimension of finite element entity, for example for tetrahedron rows are x,y,z,weight. Last row is integration weight.

FIXME: that should be moved to private class data and acessed only by member function

Examples
hanging_node_approx.cpp, and prism_elements_from_surface.cpp.

Definition at line 109 of file ForcesAndSourcesCore.hpp.

◆ getRuleHook

RuleHookFun MoFEM::ForcesAndSourcesCore::getRuleHook

Hook to get rule.

Todo:
check preferred format how works with gcc and clang, see http://www.boost.org/doc/libs/1_64_0/doc/html/function/tutorial.html#idp247873024
Examples
field_evaluator.cpp, and lorentz_force.cpp.

Definition at line 42 of file ForcesAndSourcesCore.hpp.

◆ lastEvaluatedElementEntityType

EntityType MoFEM::ForcesAndSourcesCore::lastEvaluatedElementEntityType
protected

Last evaluated type of element entity.

Definition at line 489 of file ForcesAndSourcesCore.hpp.

◆ mField

Interface& MoFEM::ForcesAndSourcesCore::mField
Examples
prism_elements_from_surface.cpp.

Definition at line 24 of file ForcesAndSourcesCore.hpp.

◆ opPtrVector

boost::ptr_deque<UserDataOperator> MoFEM::ForcesAndSourcesCore::opPtrVector
protected

Vector of finite element users data operators.

Definition at line 480 of file ForcesAndSourcesCore.hpp.

◆ refinePtrFE

ForcesAndSourcesCore* MoFEM::ForcesAndSourcesCore::refinePtrFE
private

Element to integrate parent or child.

u

Definition at line 524 of file ForcesAndSourcesCore.hpp.

◆ setRuleHook

GaussHookFun MoFEM::ForcesAndSourcesCore::setRuleHook

Set function to calculate integration rule.

Definition at line 48 of file ForcesAndSourcesCore.hpp.

◆ sidePtrFE

ForcesAndSourcesCore* MoFEM::ForcesAndSourcesCore::sidePtrFE
private

Element to integrate on the sides.

Definition at line 507 of file ForcesAndSourcesCore.hpp.

◆ userPolynomialBasePtr

boost::shared_ptr<BaseFunction> MoFEM::ForcesAndSourcesCore::userPolynomialBasePtr
private

Pointer to user polynomail base.

Definition at line 501 of file ForcesAndSourcesCore.hpp.


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