v0.10.0
Classes | Public Types | Public Member Functions | Public Attributes | Protected Member Functions | List of all members
MoFEM::ForcesAndSourcesCore Struct Reference

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

#include <src/finite_elements/ForcesAndSourcesCore.hpp>

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

Classes

struct  UserDataOperator
 Data operator to do calculations at integration points.Is inherited and implemented by user to do calculations. It can be used in many different ways but typically is used to integrate matrices (f.e. stiffness matrix) and the right hand vector (f.e. force vector). More...
 

Public Types

typedef boost::function< int(int order_row, int order_col, int order_data)> RuleHookFun
 
- 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_vector< 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...
 
const DataForcesAndSourcesCore::EntDatagetEntData (const FieldSpace space, const EntityType type, const int side) const
 Get the entity data. More...
 
DataForcesAndSourcesCore::EntDatagetEntData (const FieldSpace space, const EntityType type, const int side)
 Get the entity data. More...
 
- Public Member Functions inherited from MoFEM::FEMethod
MoFEMErrorCode query_interface (const MOFEMuuid &uuid, UnknownInterface **iface) const
 
 FEMethod ()
 
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
 
MoFEMErrorCode getNumberOfNodes (int &num_nodes) const
 Get number of DOFs on element. More...
 
EntityHandle getFEEntityHandle () const
 
MoFEMErrorCode getNodeData (const std::string field_name, VectorDouble &data, const bool reset_dofs=true)
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_begin (const MULTIINDEX &index, const std::string &field_name, const EntityType type) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_end (const MULTIINDEX &index, const std::string &field_name, const EntityType type) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_begin (const MULTIINDEX &index, const std::string &field_name) const
 
template<class MULTIINDEX >
MULTIINDEX::iterator get_end (const MULTIINDEX &index, const std::string &field_name) const
 
- Public Member Functions inherited from MoFEM::BasicMethod
 BasicMethod ()
 
virtual ~BasicMethod ()=default
 
int getNinTheLoop () const
 get number of evaluated element in the loop More...
 
int getLoopSize () const
 get loop size 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...
 
- Public Member Functions inherited from MoFEM::KspMethod
 KspMethod ()
 
virtual ~KspMethod ()=default
 
MoFEMErrorCode copyKsp (const KspMethod &ksp)
 copy data form another method More...
 
- Public Member Functions inherited from MoFEM::PetscData
 PetscData ()
 
virtual ~PetscData ()=default
 
- Public Member Functions inherited from MoFEM::UnknownInterface
template<class IFACE >
MoFEMErrorCode registerInterface (const MOFEMuuid &uuid, bool error_if_registration_failed=true)
 Register interface. More...
 
template<class IFACE , bool VERIFY = false>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, IFACE *&iface) const
 Get interface by uuid and return reference to pointer of 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
 
virtual MoFEMErrorCode getLibVersion (Version &version) const
 Get library version. More...
 
virtual const MoFEMErrorCode getFileVersion (moab::Interface &moab, Version &version) const
 Get database major version. More...
 
virtual MoFEMErrorCode getInterfaceVersion (Version &version) const
 Get database major version. More...
 
template<>
MoFEMErrorCode getInterface (const MOFEMuuid &uuid, UnknownInterface *&iface) const
 
- Public Member Functions inherited from MoFEM::SnesMethod
 SnesMethod ()
 
virtual ~SnesMethod ()=default
 
MoFEMErrorCode copySnes (const SnesMethod &snes)
 Copy snes data. More...
 
DEPRECATED MoFEMErrorCode setSnesCtx (SNESContext ctx)
 
- Public Member Functions inherited from MoFEM::TSMethod
 TSMethod ()
 
virtual ~TSMethod ()=default
 
MoFEMErrorCode copyTs (const TSMethod &ts)
 Copy TS solver data. More...
 
DEPRECATED MoFEMErrorCode setTsCtx (TSContext ctx)
 

Public Attributes

InterfacemField
 
RuleHookFun getRuleHook
 Hook to get rule. More...
 
RuleHookFun 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
 
- Public Attributes inherited from MoFEM::BasicMethod
int nInTheLoop
 number currently of processed method More...
 
int loopSize
 local number oe methods to process 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
 
- 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 More...
 
PetscReal ts_a
 shift for U_tt (see PETSc Time Solver) More...
 
PetscReal ts_v
 shift for U_t shift for U_t More...
 
PetscReal ts_t
 time 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< DataForcesAndSourcesCore::EntData > &data) const
 get sense (orientation) of entity More...
 
MoFEMErrorCode getEntityDataOrder (const EntityType type, const FieldSpace space, boost::ptr_vector< DataForcesAndSourcesCore::EntData > &data) const
 Get the entity data order. More...
 
template<EntityType type>
MoFEMErrorCode getEntitySense (DataForcesAndSourcesCore &data) const
 Get the entity sense (orientation) More...
 
template<EntityType type>
MoFEMErrorCode getEntityDataOrder (DataForcesAndSourcesCore &data, const FieldSpace space) const
 Get the entity data order for given space. More...
 
MoFEMErrorCode getFaceTriNodes (DataForcesAndSourcesCore &data) const
 Get nodes on triangles. More...
 
MoFEMErrorCode getSpacesAndBaseOnEntities (DataForcesAndSourcesCore &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 ()
 Create a entity data on element object. More...
 
MoFEMErrorCode loopOverOperators ()
 Iterate user data operators. More...
 
Indices
template<typename EXTRACTOR >
MoFEMErrorCode getNodesIndices (const std::string &field_name, FieldEntity_vector_view &ents_field, VectorInt &nodes_indices, VectorInt &local_nodes_indices, EXTRACTOR &&extractor) const
 get node indices More...
 
MoFEMErrorCode getRowNodesIndices (DataForcesAndSourcesCore &data, const std::string &field_name) const
 get row node indices from FENumeredDofEntity_multiIndex More...
 
MoFEMErrorCode getColNodesIndices (DataForcesAndSourcesCore &data, const std::string &field_name) const
 get col node indices from FENumeredDofEntity_multiIndex More...
 
template<typename EXTRACTOR >
MoFEMErrorCode getEntityIndices (DataForcesAndSourcesCore &data, const std::string &field_name, FieldEntity_vector_view &ents_field, const EntityType type_lo, const EntityType type_hi, EXTRACTOR &&extractor) const
 
MoFEMErrorCode getEntityRowIndices (DataForcesAndSourcesCore &data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
 
MoFEMErrorCode getEntityColIndices (DataForcesAndSourcesCore &data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
 
MoFEMErrorCode getNoFieldIndices (const std::string &field_name, boost::shared_ptr< FENumeredDofEntity_multiIndex > dofs, VectorInt &nodes_indices) const
 get NoField indices More...
 
MoFEMErrorCode getNoFieldRowIndices (DataForcesAndSourcesCore &data, const std::string &field_name) const
 get col NoField indices More...
 
MoFEMErrorCode getNoFieldColIndices (DataForcesAndSourcesCore &data, const std::string &field_name) const
 get col NoField indices More...
 
Data
MoFEMErrorCode getNoFieldFieldData (const std::string field_name, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs, VectorFieldEntities &ent_field) const
 Get field data on nodes. More...
 
MoFEMErrorCode getNoFieldFieldData (DataForcesAndSourcesCore &data, const std::string field_name) const
 
MoFEMErrorCode getNodesFieldData (DataForcesAndSourcesCore &data, const std::string &field_name) const
 Get data on nodes. More...
 
MoFEMErrorCode getEntityFieldData (DataForcesAndSourcesCore &data, const std::string &field_name, 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
 
- Protected Member Functions inherited from MoFEM::UnknownInterface
boost::typeindex::type_index getClassIdx (const MOFEMuuid &uid) const
 Get type name for interface Id. More...
 
MOFEMuuid getUId (const boost::typeindex::type_index &class_idx) const
 Get interface Id for class name. More...
 

Deprecated (do not use)

const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACEdataOnElement
 Entity data on element entity rows fields. More...
 
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACEderivedDataOnElement
 Entity data on element entity columns fields. More...
 
DataForcesAndSourcesCoredataNoField
 
DataForcesAndSourcesCoredataH1
 
DataForcesAndSourcesCoredataHcurl
 
DataForcesAndSourcesCoredataHdiv
 
DataForcesAndSourcesCoredataL2
 
boost::ptr_vector< UserDataOperatoropPtrVector
 Vector of finite element users data operators. More...
 
EntityType lastEvaluatedElementEntityType
 Last evaluated type of element entity. More...
 
class UserDataOperator
 
class VolumeElementForcesAndSourcesCoreOnSideBase
 
class FaceElementForcesAndSourcesCoreOnSideBase
 
class VolumeElementForcesAndSourcesCoreOnContactPrismSideBase
 
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...
 
virtual int getRule (int order)
 
virtual MoFEMErrorCode setGaussPts (int order)
 
MoFEMErrorCode setSideFEPtr (const ForcesAndSourcesCore *side_fe_ptr)
 Set the pointer to face element on the side. More...
 

Additional Inherited Members

- 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 DataForcesAndSourcesCore

Examples
EshelbianPlasticity.cpp, hello_world.cpp, prism_elements_from_surface.cpp, prism_polynomial_approximation.cpp, and quad_polynomial_approximation.cpp.

Definition at line 34 of file ForcesAndSourcesCore.hpp.

Member Typedef Documentation

◆ RuleHookFun

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

Definition at line 41 of file ForcesAndSourcesCore.hpp.

Constructor & Destructor Documentation

◆ ForcesAndSourcesCore()

MoFEM::ForcesAndSourcesCore::ForcesAndSourcesCore ( Interface m_field)

Definition at line 54 of file ForcesAndSourcesCore.cpp.

55  :
56 
57  mField(m_field), getRuleHook(0), setRuleHook(0),
59 
60  nullptr,
61  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // NOFIELD
62  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // H1
63  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // HCURL
64  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // HDIV
65  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET) // L2
66 
67  },
69 
70  nullptr,
71  boost::make_shared<DerivedDataForcesAndSourcesCore>(
72  dataOnElement[NOFIELD]), // NOFIELD
73  boost::make_shared<DerivedDataForcesAndSourcesCore>(
74  dataOnElement[H1]), // H1
75  boost::make_shared<DerivedDataForcesAndSourcesCore>(
76  dataOnElement[HCURL]), // HCURL
77  boost::make_shared<DerivedDataForcesAndSourcesCore>(
78  dataOnElement[HDIV]), // HDIV
79  boost::make_shared<DerivedDataForcesAndSourcesCore>(
80  dataOnElement[L2]) // L2
81 
82  },
85  dataHdiv(*dataOnElement[HDIV].get()), dataL2(*dataOnElement[L2].get()),
86  lastEvaluatedElementEntityType(MBMAXTYPE), sidePtrFE(nullptr) {}
DataForcesAndSourcesCore & dataL2
field with continuous normal traction
Definition: definitions.h:179
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.
DataForcesAndSourcesCore & dataH1
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
RuleHookFun getRuleHook
Hook to get rule.
DataForcesAndSourcesCore & dataHcurl
EntityType lastEvaluatedElementEntityType
Last evaluated type of element entity.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnElement
Entity data on element entity columns fields.
field with continuous tangents
Definition: definitions.h:178
DataForcesAndSourcesCore & dataHdiv
continuous field
Definition: definitions.h:177
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
RuleHookFun setRuleHook
Set function to calculate integration rule.
field with C-1 continuity
Definition: definitions.h:180
DataForcesAndSourcesCore & dataNoField

Member Function Documentation

◆ calBernsteinBezierBaseFunctionsOnElement()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::calBernsteinBezierBaseFunctionsOnElement ( )
protected

Calculate Bernstein-Bezier base.

Returns
MoFEMErrorCode

Definition at line 1079 of file ForcesAndSourcesCore.cpp.

1079  {
1081 
1082  auto get_nodal_base_data = [&](DataForcesAndSourcesCore &data,
1083  auto field_ptr) {
1085  auto &space = data.dataOnEntities[MBVERTEX][0].getSpace();
1086  auto &base = data.dataOnEntities[MBVERTEX][0].getBase();
1087  auto &bb_node_order = data.dataOnEntities[MBVERTEX][0].getBBNodeOrder();
1088 
1089  auto &field_ents = getDataFieldEnts();
1090  auto bit_number = field_ptr->getBitNumber();
1091  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
1092  bit_number, get_id_for_min_type<MBVERTEX>());
1093  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
1094  cmp_uid_lo);
1095  if (lo != field_ents.end()) {
1096  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
1097  bit_number, get_id_for_max_type<MBVERTEX>());
1098  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
1099  if (lo != hi) {
1100 
1101  for (auto it = lo; it != hi; ++it)
1102  if (auto first_e = it->lock()) {
1103  space = first_e->getSpace();
1104  base = first_e->getApproxBase();
1105  int num_nodes;
1106  CHKERR getNumberOfNodes(num_nodes);
1107  bb_node_order.resize(num_nodes, false);
1108  bb_node_order.clear();
1109  const int nb_dof_idx = first_e->getNbOfCoeffs();
1110 
1111  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
1112 
1113  for (; it != hi; ++it) {
1114  if (auto e = it->lock()) {
1115  const auto &sn = e->getSideNumberPtr();
1116  const int side_number = sn->side_number;
1117  const int brother_side_number = sn->brother_side_number;
1118  if (brother_side_number != -1)
1119  brother_ents_vec.emplace_back(e);
1120  bb_node_order[side_number] = e->getMaxOrder();
1121  }
1122  }
1123 
1124  for (auto &it : brother_ents_vec) {
1125  if (const auto e = it.lock()) {
1126  const auto &sn = e->getSideNumberPtr();
1127  const int side_number = sn->side_number;
1128  const int brother_side_number = sn->brother_side_number;
1129  bb_node_order[brother_side_number] = bb_node_order[side_number];
1130  }
1131  }
1132 
1133  break;
1134  }
1135  }
1136  }
1137 
1139  };
1140 
1141  auto get_entity_base_data = [&](DataForcesAndSourcesCore &data,
1142  auto field_ptr, const EntityType type_lo,
1143  const EntityType type_hi) {
1145  for (EntityType t = type_lo; t != type_hi; ++t) {
1146  for (auto &dat : data.dataOnEntities[t]) {
1147  dat.getDataOrder() = 0;
1148  dat.getBase() = NOBASE;
1149  dat.getSpace() = NOSPACE;
1150  dat.getFieldData().resize(0, false);
1151  dat.getFieldDofs().resize(0, false);
1152  }
1153  }
1154 
1155  auto &field_ents = getDataFieldEnts();
1156  auto bit_number = field_ptr->getBitNumber();
1157  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
1158  bit_number, get_id_for_min_type(type_lo));
1159  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
1160  cmp_uid_lo);
1161  if (lo != field_ents.end()) {
1162  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
1163  bit_number, get_id_for_max_type(type_hi));
1164  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
1165 
1166  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
1167  for (; lo != hi; ++lo) {
1168  if (auto e = lo->lock()) {
1169  if (auto cache = e->entityCacheDataDofs.lock()) {
1170  if (cache->loHi[0] != cache->loHi[1]) {
1171  const EntityType type = e->getEntType();
1172  const int side = e->getSideNumberPtr()->side_number;
1173  auto &dat = data.dataOnEntities[type][side];
1174  const int brother_side =
1175  e->getSideNumberPtr()->brother_side_number;
1176  if (brother_side != -1)
1177  brother_ents_vec.emplace_back(e);
1178  dat.getBase() = e->getApproxBase();
1179  dat.getSpace() = e->getSpace();
1180  const auto ent_order = e->getMaxOrder();
1181  dat.getDataOrder() = dat.getDataOrder() > ent_order
1182  ? dat.getDataOrder()
1183  : ent_order;
1184  }
1185  }
1186  }
1187  }
1188 
1189  for (auto &ent_ptr : brother_ents_vec) {
1190  if (auto e = ent_ptr.lock()) {
1191  const EntityType type = e->getEntType();
1192  const int side = e->getSideNumberPtr()->side_number;
1193  const int brother_side = e->getSideNumberPtr()->brother_side_number;
1194  auto &dat = data.dataOnEntities[type][side];
1195  auto &dat_brother = data.dataOnEntities[type][brother_side];
1196  dat_brother.getBase() = dat.getBase();
1197  dat_brother.getSpace() = dat.getSpace();
1198  dat_brother.getDataOrder() = dat.getDataOrder();
1199  }
1200  }
1201  }
1203  };
1204 
1205  for (auto &e : getDataFieldEnts()) {
1206  if (auto ent_data_ptr = e.lock()) {
1207  if (ent_data_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
1208  auto space = ent_data_ptr->getSpace();
1209  for (EntityType t = MBVERTEX; t != MBPOLYHEDRON; ++t) {
1210  for (auto &dat : (*dataOnElement[space]).dataOnEntities[t]) {
1211  for (auto &ptr : dat.getBBAlphaIndicesByOrderArray())
1212  ptr.reset();
1213  for (auto &ptr : dat.getBBNByOrderArray())
1214  ptr.reset();
1215  for (auto &ptr : dat.getBBDiffNByOrderArray())
1216  ptr.reset();
1217  }
1218  }
1219  }
1220  }
1221  }
1222 
1223  std::set<string> fields_list;
1224  for (auto &e : getDataFieldEnts()) {
1225  if (auto ent_data_ptr = e.lock()) {
1226  if (ent_data_ptr->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
1227  auto field_name = ent_data_ptr->getName();
1228  if (fields_list.find(field_name) == fields_list.end()) {
1229  auto field_ptr = ent_data_ptr->getFieldRawPtr();
1230  auto space = ent_data_ptr->getSpace();
1231  CHKERR get_nodal_base_data(*dataOnElement[space], field_ptr);
1232  CHKERR get_entity_base_data(*dataOnElement[space], field_ptr, MBEDGE,
1233  MBPOLYHEDRON);
1234  CHKERR getElementPolynomialBase()->getValue(
1235  gaussPts, boost::make_shared<EntPolynomialBaseCtx>(
1236  *dataOnElement[space], field_name,
1237  static_cast<FieldSpace>(space),
1239  fields_list.insert(field_name);
1240  }
1241  }
1242  }
1243  }
1244 
1246 };
EntityHandle get_id_for_max_type()
EntityHandle get_id_for_min_type()
MatrixDouble gaussPts
Matrix of integration points.
BlockParamData * cache
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
const FieldEntity_vector_view & getDataFieldEnts() const
static auto cmp_uid_lo(const boost::weak_ptr< FieldEntity > &a, const UId &b)
static auto cmp_uid_hi(const UId &b, const boost::weak_ptr< FieldEntity > &a)
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
MoFEMErrorCode getNumberOfNodes(int &num_nodes) const
Get number of DOFs on element.
auto & getElementPolynomialBase()
Get the Entity Polynomial Base object.

◆ calHierarchicalBaseFunctionsOnElement() [1/2]

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

Calculate base functions.

Returns
Error code

Definition at line 1009 of file ForcesAndSourcesCore.cpp.

1010  {
1012  if (dataOnElement[H1]->bAse.test(b)) {
1013  switch (static_cast<FieldApproximationBase>(b)) {
1014  case NOBASE:
1015  break;
1017  // BERNSTEIN_BEZIER_BASE is not hierarchical base
1018  break;
1021  case DEMKOWICZ_JACOBI_BASE:
1022  if (!getElementPolynomialBase())
1023  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1024  "Functions genrating approximation base not defined");
1025 
1026  for (int space = H1; space != LASTSPACE; ++space) {
1027  if (dataOnElement[H1]->sPace.test(space) &&
1028  dataOnElement[H1]->bAse.test(b) &&
1029  dataOnElement[H1]->basesOnSpaces[space].test(b)) {
1030  CHKERR getElementPolynomialBase()->getValue(
1031  gaussPts,
1032  boost::make_shared<EntPolynomialBaseCtx>(
1033  *dataOnElement[space], static_cast<FieldSpace>(space),
1034  static_cast<FieldApproximationBase>(b), NOBASE));
1035  }
1036  }
1037  break;
1038  case USER_BASE:
1039  if (!getUserPolynomialBase())
1040  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1041  "Functions genrating user approximation base not defined");
1042 
1043  for (int space = H1; space != LASTSPACE; ++space)
1044  if (dataOnElement[H1]->sPace.test(space) &&
1045  dataOnElement[H1]->bAse.test(b) &&
1046  dataOnElement[H1]->basesOnSpaces[space].test(b)) {
1047  CHKERR getUserPolynomialBase()->getValue(
1048  gaussPts,
1049  boost::make_shared<EntPolynomialBaseCtx>(
1050  *dataOnElement[space], static_cast<FieldSpace>(space),
1051  static_cast<FieldApproximationBase>(b), NOBASE));
1052  }
1053  break;
1054  default:
1056  "Base <%s> not yet implemented",
1057  ApproximationBaseNames[static_cast<FieldApproximationBase>(b)]);
1058  }
1059  }
1061 }
MatrixDouble gaussPts
Matrix of integration points.
user implemented approximation base
Definition: definitions.h:160
auto & getUserPolynomialBase()
Get the User Polynomial Base object.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:181
static const char *const ApproximationBaseNames[]
Definition: definitions.h:164
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:152
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
virtual MPI_Comm & get_comm() const =0
continuous field
Definition: definitions.h:177
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
auto & getElementPolynomialBase()
Get the Entity Polynomial Base object.

◆ calHierarchicalBaseFunctionsOnElement() [2/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::calHierarchicalBaseFunctionsOnElement ( )
protected

Calculate base functions.

Returns
Error code

Definition at line 1063 of file ForcesAndSourcesCore.cpp.

1063  {
1065  /// Use the some node base. Node base is usually used for construction other
1066  /// bases.
1067  for (int space = HCURL; space != LASTSPACE; ++space) {
1068  dataOnElement[space]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
1069  dataOnElement[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
1070  }
1071  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
1073  static_cast<FieldApproximationBase>(b));
1074  }
1076 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:181
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:152
field with continuous tangents
Definition: definitions.h:178
#define CHKERR
Inline error check.
Definition: definitions.h:604
MoFEMErrorCode calHierarchicalBaseFunctionsOnElement()
Calculate base functions.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
continuous field
Definition: definitions.h:177
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.

◆ createDataOnElement()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::createDataOnElement ( )
protected

Create a entity data on element object.

Returns
MoFEMErrorCode

Definition at line 1248 of file ForcesAndSourcesCore.cpp.

1248  {
1250 
1251  const EntityType type = numeredEntFiniteElementPtr->getEntType();
1254 
1255  // Data on elements for proper spaces
1256  for (int space = H1; space != LASTSPACE; ++space) {
1257  dataOnElement[space]->setElementType(type);
1258  derivedDataOnElement[space]->setElementType(type);
1259  }
1260 
1262 
1264 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:181
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
EntityType lastEvaluatedElementEntityType
Last evaluated type of element entity.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnElement
Entity data on element entity columns fields.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
continuous field
Definition: definitions.h:177
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.

◆ getColNodesIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getColNodesIndices ( DataForcesAndSourcesCore data,
const std::string &  field_name 
) const
protected

get col node indices from FENumeredDofEntity_multiIndex

Definition at line 332 of file ForcesAndSourcesCore.cpp.

333  {
334 
335  struct Extractor {
336  boost::weak_ptr<EntityCacheNumeredDofs>
337  operator()(boost::shared_ptr<FieldEntity> &e) {
338  return e->entityCacheColDofs;
339  }
340  };
341 
342  return getNodesIndices(field_name, getColFieldEnts(),
343  data.dataOnEntities[MBVERTEX][0].getIndices(),
344  data.dataOnEntities[MBVERTEX][0].getLocalIndices(),
345  Extractor());
346 }
MoFEMErrorCode getNodesIndices(const std::string &field_name, FieldEntity_vector_view &ents_field, VectorInt &nodes_indices, VectorInt &local_nodes_indices, EXTRACTOR &&extractor) const
get node indices
auto & getColFieldEnts() const
virtual MoFEMErrorCode operator()()
function is run for every finite element

◆ getElementPolynomialBase()

auto& MoFEM::ForcesAndSourcesCore::getElementPolynomialBase ( )

Get the Entity Polynomial Base object.

Returns
boost::shared_ptr<BaseFunction>&&

Definition at line 430 of file ForcesAndSourcesCore.hpp.

430 { return elementPolynomialBasePtr; }
boost::shared_ptr< BaseFunction > elementPolynomialBasePtr
Pointer to entity polynomial base.

◆ getEntData() [1/2]

const DataForcesAndSourcesCore::EntData& MoFEM::ForcesAndSourcesCore::getEntData ( const FieldSpace  space,
const EntityType  type,
const int  side 
) const

Get the entity data.

Parameters
space
type
side
Returns
const DataForcesAndSourcesCore::EntData&

Definition at line 484 of file ForcesAndSourcesCore.hpp.

486  {
487  return dataOnElement[space]->dataOnEntities[type][side];
488  }
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.

◆ getEntData() [2/2]

DataForcesAndSourcesCore::EntData& MoFEM::ForcesAndSourcesCore::getEntData ( const FieldSpace  space,
const EntityType  type,
const int  side 
)

Get the entity data.

Parameters
space
type
side
Returns
DataForcesAndSourcesCore::EntData&

Definition at line 499 of file ForcesAndSourcesCore.hpp.

499  {
500  return dataOnElement[space]->dataOnEntities[type][side];
501  }
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.

◆ getEntityColIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntityColIndices ( DataForcesAndSourcesCore data,
const std::string &  field_name,
const EntityType  type_lo = MBVERTEX,
const EntityType  type_hi = MBPOLYHEDRON 
) const
protected

Definition at line 433 of file ForcesAndSourcesCore.cpp.

435  {
436 
437  struct Extractor {
438  boost::weak_ptr<EntityCacheNumeredDofs>
439  operator()(boost::shared_ptr<FieldEntity> &e) {
440  return e->entityCacheColDofs;
441  }
442  };
443 
444  return getEntityIndices(data, field_name, getColFieldEnts(), type_lo, type_hi,
445  Extractor());
446 }
MoFEMErrorCode getEntityIndices(DataForcesAndSourcesCore &data, const std::string &field_name, FieldEntity_vector_view &ents_field, const EntityType type_lo, const EntityType type_hi, EXTRACTOR &&extractor) const
auto & getColFieldEnts() const
virtual MoFEMErrorCode operator()()
function is run for every finite element

◆ getEntityDataOrder() [1/2]

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

Get the entity data order.

Parameters
type
space
data
Returns
MoFEMErrorCode

Definition at line 146 of file ForcesAndSourcesCore.cpp.

148  {
150 
151  auto set_order = [&]() {
153  auto &side_table = numeredEntFiniteElementPtr->getSideNumberTable();
154 
155  for (unsigned int s = 0; s != data.size(); ++s)
156  data[s].getDataOrder() = 0;
157 
158  const FieldEntity_vector_view &data_field_ent = getDataFieldEnts();
159 
160  for (auto r = fieldsPtr->get<BitFieldId_space_mi_tag>().equal_range(space);
161  r.first != r.second; ++r.first) {
162 
163  const auto field_bit_number = (*r.first)->getBitNumber();
164  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
165  field_bit_number, get_id_for_min_type(type));
166  auto lo = std::lower_bound(data_field_ent.begin(), data_field_ent.end(),
167  lo_uid, cmp_uid_lo);
168  if (lo != data_field_ent.end()) {
169  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
170  field_bit_number, get_id_for_max_type(type));
171  auto hi =
172  std::upper_bound(lo, data_field_ent.end(), hi_uid, cmp_uid_hi);
173  for (; lo != hi; ++lo) {
174 
175  if (auto ptr = lo->lock()) {
176 
177  auto &e = *ptr;
178  auto sit = side_table.find(e.getEnt());
179  if (sit != side_table.end()) {
180  auto &side = *sit;
181  const int side_number = side->side_number;
182  if (side_number >= 0) {
183  ApproximationOrder ent_order = e.getMaxOrder();
184  auto &dat = data[side_number];
185  dat.getDataOrder() = dat.getDataOrder() > ent_order
186  ? dat.getDataOrder()
187  : ent_order;
188  }
189  } else
190  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
191  "Entity on side of the element not found");
192  }
193  }
194  }
195  }
196 
198  };
199 
200  auto set_order_on_brother = [&]() {
202  auto &side_table =
203  numeredEntFiniteElementPtr->getSideNumberTable().get<0>();
204  auto sit = side_table.lower_bound(get_id_for_min_type(type));
205  if (sit != side_table.end()) {
206  auto hi_sit = side_table.upper_bound(get_id_for_max_type(type));
207  for (; sit != hi_sit; ++sit) {
208  const int brother_side_number = (*sit)->brother_side_number;
209  if (brother_side_number != -1) {
210  const int side_number = (*sit)->side_number;
211  data[brother_side_number].getDataOrder() =
212  data[side_number].getDataOrder();
213  }
214  }
215  }
217  };
218 
219  CHKERR set_order();
220  CHKERR set_order_on_brother();
221 
223 }
EntityHandle get_id_for_max_type()
EntityHandle get_id_for_min_type()
const Field_multiIndex * fieldsPtr
raw pointer to fields container
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
int ApproximationOrder
Approximation on the entity.
Definition: Types.hpp:37
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
std::vector< boost::weak_ptr< FieldEntity > > FieldEntity_vector_view
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
const FieldEntity_vector_view & getDataFieldEnts() const
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
static auto cmp_uid_lo(const boost::weak_ptr< FieldEntity > &a, const UId &b)
static auto cmp_uid_hi(const UId &b, const boost::weak_ptr< FieldEntity > &a)
#define CHKERR
Inline error check.
Definition: definitions.h:604
const double r
rate factor
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ getEntityDataOrder() [2/2]

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

Get the entity data order for given space.

Template Parameters
type
Parameters
data
space
Returns
MoFEMErrorCode

Definition at line 547 of file ForcesAndSourcesCore.hpp.

548  {
549  return getEntityDataOrder(type, space, data.dataOnEntities[type]);
550  }
MoFEMErrorCode getEntityDataOrder(const EntityType type, const FieldSpace space, boost::ptr_vector< DataForcesAndSourcesCore::EntData > &data) const
Get the entity data order.

◆ getEntityFieldData()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntityFieldData ( DataForcesAndSourcesCore data,
const std::string &  field_name,
const EntityType  type_lo = MBVERTEX,
const EntityType  type_hi = MBPOLYHEDRON 
) const
protected

Definition at line 716 of file ForcesAndSourcesCore.cpp.

718  {
720  for (EntityType t = type_lo; t != type_hi; ++t) {
721  for (auto &dat : data.dataOnEntities[t]) {
722  dat.getDataOrder() = 0;
723  dat.getBase() = NOBASE;
724  dat.getSpace() = NOSPACE;
725  dat.getFieldData().resize(0, false);
726  dat.getFieldDofs().resize(0, false);
727  dat.getFieldEntities().resize(0, false);
728  }
729  }
730 
731  auto &field_ents = getDataFieldEnts();
732  auto bit_number = mField.get_field_bit_number(field_name);
733  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
734  bit_number, get_id_for_min_type(type_lo));
735  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
736  cmp_uid_lo);
737  if (lo != field_ents.end()) {
738  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
739  bit_number, get_id_for_max_type(type_hi));
740  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
741  if (lo != hi) {
742 
743  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
744 
745  for (auto it = lo; it != hi; ++it)
746  if (auto e = it->lock()) {
747 
748  const EntityType type = e->getEntType();
749  auto side_ptr = e->getSideNumberPtr();
750  const int side = side_ptr->side_number;
751  if (side >= 0) {
752 
753  auto &dat = data.dataOnEntities[type][side];
754  auto &ent_field = dat.getFieldEntities();
755  auto &ent_field_dofs = dat.getFieldDofs();
756  auto &ent_field_data = dat.getFieldData();
757 
758  const int brother_side = side_ptr->brother_side_number;
759  if (brother_side != -1)
760  brother_ents_vec.emplace_back(e);
761 
762  dat.getBase() = e->getApproxBase();
763  dat.getSpace() = e->getSpace();
764  const int ent_order = e->getMaxOrder();
765  dat.getDataOrder() =
766  dat.getDataOrder() > ent_order ? dat.getDataOrder() : ent_order;
767 
768  auto ent_data = e->getEntFieldData();
769  ent_field_data.resize(ent_data.size(), false);
770  noalias(ent_field_data) = ent_data;
771  ent_field_dofs.resize(ent_data.size(), false);
772  std::fill(ent_field_dofs.begin(), ent_field_dofs.end(), nullptr);
773  ent_field.resize(1, false);
774  ent_field[0] = e.get();
775  if (auto cache = e->entityCacheDataDofs.lock()) {
776  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
777  ent_field_dofs[(*dit)->getEntDofIdx()] =
778  reinterpret_cast<FEDofEntity *>((*dit).get());
779  }
780  }
781  }
782  }
783 
784  for (auto &it : brother_ents_vec) {
785  if (const auto e = it.lock()) {
786  const EntityType type = e->getEntType();
787  const int side = e->getSideNumberPtr()->side_number;
788  const int brother_side = e->getSideNumberPtr()->brother_side_number;
789  auto &dat = data.dataOnEntities[type][side];
790  auto &dat_brother = data.dataOnEntities[type][brother_side];
791  dat_brother.getBase() = dat.getBase();
792  dat_brother.getSpace() = dat.getSpace();
793  dat_brother.getDataOrder() = dat.getDataOrder();
794  dat_brother.getFieldData() = dat.getFieldData();
795  dat_brother.getFieldDofs() = dat.getFieldDofs();
796  dat_brother.getFieldEntities() = dat.getFieldEntities();
797  }
798  }
799  }
800  }
801 
803 }
EntityHandle get_id_for_max_type()
EntityHandle get_id_for_min_type()
BlockParamData * cache
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
const FieldEntity_vector_view & getDataFieldEnts() const
static auto cmp_uid_lo(const boost::weak_ptr< FieldEntity > &a, const UId &b)
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
static auto cmp_uid_hi(const UId &b, const boost::weak_ptr< FieldEntity > &a)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ getEntityIndices()

template<typename EXTRACTOR >
MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntityIndices ( DataForcesAndSourcesCore data,
const std::string &  field_name,
FieldEntity_vector_view ents_field,
const EntityType  type_lo,
const EntityType  type_hi,
EXTRACTOR &&  extractor 
) const
protected

Definition at line 349 of file ForcesAndSourcesCore.cpp.

352  {
354 
355  for (EntityType t = type_lo; t != type_hi; ++t) {
356  for (auto &dat : data.dataOnEntities[t]) {
357  dat.getIndices().resize(0, false);
358  dat.getLocalIndices().resize(0, false);
359  }
360  }
361 
362  auto bit_number = mField.get_field_bit_number(field_name);
363  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
364  bit_number, get_id_for_min_type(type_lo));
365  auto lo = std::lower_bound(ents_field.begin(), ents_field.end(), lo_uid,
366  cmp_uid_lo);
367  if (lo != ents_field.end()) {
368  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
369  bit_number, get_id_for_max_type(type_hi));
370  auto hi = std::upper_bound(lo, ents_field.end(), hi_uid, cmp_uid_hi);
371  if (lo != hi) {
372 
373  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
374 
375  for (auto it = lo; it != hi; ++it)
376  if (auto e = it->lock()) {
377 
378  const EntityType type = e->getEntType();
379  auto side_ptr = e->getSideNumberPtr();
380  const int side = side_ptr->side_number;
381  if (side >= 0) {
382  const int nb_dofs_on_ent = e->getNbDofsOnEnt();
383  const int brother_side = side_ptr->brother_side_number;
384  auto &dat = data.dataOnEntities[type][side];
385  auto &ent_field_indices = dat.getIndices();
386  auto &ent_field_local_indices = dat.getLocalIndices();
387 
388  ent_field_indices.resize(nb_dofs_on_ent, false);
389  ent_field_local_indices.resize(nb_dofs_on_ent, false);
390  std::fill(ent_field_indices.data().begin(),
391  ent_field_indices.data().end(), -1);
392  std::fill(ent_field_local_indices.data().begin(),
393  ent_field_local_indices.data().end(), -1);
394 
395  if (auto cache = extractor(e).lock()) {
396  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
397  const int idx = (*dit)->getEntDofIdx();
398  ent_field_indices[idx] = (*dit)->getPetscGlobalDofIdx();
399  ent_field_local_indices[idx] = (*dit)->getPetscLocalDofIdx();
400  }
401  }
402 
403  if (brother_side != -1) {
404  auto &dat_brother = data.dataOnEntities[type][brother_side];
405  dat_brother.getIndices().resize(nb_dofs_on_ent, false);
406  dat_brother.getLocalIndices().resize(nb_dofs_on_ent, false);
407  noalias(dat_brother.getIndices()) = dat.getIndices();
408  noalias(dat_brother.getLocalIndices()) = dat.getLocalIndices();
409  }
410  }
411  }
412  }
413  }
414 
416 }
EntityHandle get_id_for_max_type()
EntityHandle get_id_for_min_type()
BlockParamData * cache
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
static auto cmp_uid_lo(const boost::weak_ptr< FieldEntity > &a, const UId &b)
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
static auto cmp_uid_hi(const UId &b, const boost::weak_ptr< FieldEntity > &a)
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ getEntityRowIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntityRowIndices ( DataForcesAndSourcesCore data,
const std::string &  field_name,
const EntityType  type_lo = MBVERTEX,
const EntityType  type_hi = MBPOLYHEDRON 
) const
protected

Definition at line 418 of file ForcesAndSourcesCore.cpp.

420  {
421 
422  struct Extractor {
423  boost::weak_ptr<EntityCacheNumeredDofs>
424  operator()(boost::shared_ptr<FieldEntity> &e) {
425  return e->entityCacheRowDofs;
426  }
427  };
428 
429  return getEntityIndices(data, field_name, getRowFieldEnts(), type_lo, type_hi,
430  Extractor());
431 }
auto & getRowFieldEnts() const
MoFEMErrorCode getEntityIndices(DataForcesAndSourcesCore &data, const std::string &field_name, FieldEntity_vector_view &ents_field, const EntityType type_lo, const EntityType type_hi, EXTRACTOR &&extractor) const
virtual MoFEMErrorCode operator()()
function is run for every finite element

◆ getEntitySense() [1/2]

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

get sense (orientation) of entity

Parameters
typetype of entity
dataentity data
Returns
error code

Definition at line 90 of file ForcesAndSourcesCore.cpp.

92  {
94 
95  auto &side_table = numeredEntFiniteElementPtr->getSideNumberTable().get<0>();
96  auto sit = side_table.lower_bound(get_id_for_min_type(type));
97  if (sit != side_table.end()) {
98  auto hi_sit = side_table.upper_bound(get_id_for_max_type(type));
99  for (; sit != hi_sit; ++sit) {
100  const int side_number = (*sit)->side_number;
101  if (side_number >= 0) {
102  const int brother_side_number = (*sit)->brother_side_number;
103  const int sense = (*sit)->sense;
104 
105  data[side_number].getSense() = sense;
106  if (brother_side_number != -1)
107  data[brother_side_number].getSense() = sense;
108  }
109  }
110  }
112 }
EntityHandle get_id_for_max_type()
EntityHandle get_id_for_min_type()
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ getEntitySense() [2/2]

template<EntityType type>
MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getEntitySense ( DataForcesAndSourcesCore data) const
protected

Get the entity sense (orientation)

Template Parameters
type
Parameters
data
Returns
MoFEMErrorCode

Definition at line 534 of file ForcesAndSourcesCore.hpp.

534  {
535  return getEntitySense(type, data.dataOnEntities[type]);
536  }
MoFEMErrorCode getEntitySense(const EntityType type, boost::ptr_vector< DataForcesAndSourcesCore::EntData > &data) const
get sense (orientation) of entity

◆ getFaceTriNodes()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getFaceTriNodes ( DataForcesAndSourcesCore data) const
protected

Get nodes on triangles.

Definition at line 873 of file ForcesAndSourcesCore.cpp.

873  {
875  // PetscAttachDebugger();
876  data.facesNodes.resize(4, 3, false);
877  auto &side_table = const_cast<SideNumber_multiIndex &>(
878  numeredEntFiniteElementPtr->getSideNumberTable());
879  auto siit = side_table.get<1>().lower_bound(boost::make_tuple(MBTRI, 0));
880  auto hi_siit = side_table.get<1>().upper_bound(boost::make_tuple(MBTRI, 4));
881  if (std::distance(siit, hi_siit) != 4) {
882  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
883  "Should be 4 triangles on tet, side_table not initialized");
884  }
885  const int canonical_face_sense_p1[4][3] = {
886  {0, 1, 3},
887  {1, 2, 3},
888  {0, 3, 2} /**/,
889  {0, 2, 1} /**/}; // second index is offset (positive sense)
890  const int canonical_face_sense_m1[4][3] = {
891  {0, 3, 1},
892  {1, 3, 2},
893  {0, 2, 3},
894  {0, 1, 2}}; // second index is offset (negative sense
895  for (; siit != hi_siit; siit++) {
896  const boost::shared_ptr<SideNumber> side = *siit;
897  int face_conn[3] = {-1, -1, -1};
898  if (side->offset == 0) {
899  face_conn[0] = side->sense == 1
900  ? canonical_face_sense_p1[(int)side->side_number][0]
901  : canonical_face_sense_m1[(int)side->side_number][0];
902  face_conn[1] = side->sense == 1
903  ? canonical_face_sense_p1[(int)side->side_number][1]
904  : canonical_face_sense_m1[(int)side->side_number][1];
905  face_conn[2] = side->sense == 1
906  ? canonical_face_sense_p1[(int)side->side_number][2]
907  : canonical_face_sense_m1[(int)side->side_number][2];
908  }
909  if (side->offset == 1) {
910  face_conn[0] =
911  side->sense == 1
912  ? canonical_face_sense_p1[(int)side->side_number][1]
913  : canonical_face_sense_m1[(int)side->side_number][2] /**/;
914  face_conn[1] = side->sense == 1
915  ? canonical_face_sense_p1[(int)side->side_number][2]
916  : canonical_face_sense_m1[(int)side->side_number][0];
917  face_conn[2] = side->sense == 1
918  ? canonical_face_sense_p1[(int)side->side_number][0]
919  : canonical_face_sense_m1[(int)side->side_number][1];
920  }
921  if (side->offset == 2) {
922  face_conn[0] =
923  side->sense == 1
924  ? canonical_face_sense_p1[(int)side->side_number][2]
925  : canonical_face_sense_m1[(int)side->side_number][1] /**/;
926  face_conn[1] = side->sense == 1
927  ? canonical_face_sense_p1[(int)side->side_number][0]
928  : canonical_face_sense_m1[(int)side->side_number][2];
929  face_conn[2] = side->sense == 1
930  ? canonical_face_sense_p1[(int)side->side_number][1]
931  : canonical_face_sense_m1[(int)side->side_number][0];
932  }
933  for (int nn = 0; nn < 3; nn++)
934  data.facesNodes(side->side_number, nn) = face_conn[nn];
935  {
936  const EntityHandle *conn_tet;
937  int num_nodes_tet;
939  CHKERR mField.get_moab().get_connectivity(ent, conn_tet, num_nodes_tet,
940  true);
941  if (num_nodes_tet != 4)
942  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
943  "data inconsistency");
944  int num_nodes_face;
945  const EntityHandle *conn_face;
946  CHKERR mField.get_moab().get_connectivity(side->ent, conn_face,
947  num_nodes_face, true);
948  if (num_nodes_face != 3)
949  SETERRQ(PETSC_COMM_SELF, 1, "data inconsistency");
950  if (conn_face[0] != conn_tet[data.facesNodes(side->side_number, 0)])
951  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
952  "data inconsistency");
953  if (conn_face[1] != conn_tet[data.facesNodes(side->side_number, 1)])
954  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
955  "data inconsistency");
956  if (conn_face[2] != conn_tet[data.facesNodes(side->side_number, 2)])
957  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
958  "data inconsistency");
959  }
960  }
962 }
virtual moab::Interface & get_moab()=0
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ getMaxColOrder()

int MoFEM::ForcesAndSourcesCore::getMaxColOrder ( ) const

Get max order of approximation for field in columns.

Definition at line 142 of file ForcesAndSourcesCore.cpp.

142  {
143  return getMaxOrder(getColFieldEnts());
144 }
auto & getColFieldEnts() const
static int getMaxOrder(const ENTMULTIINDEX &multi_index)

◆ getMaxDataOrder()

int MoFEM::ForcesAndSourcesCore::getMaxDataOrder ( ) const

Get max order of approximation for data fields.

Method 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 127 of file ForcesAndSourcesCore.cpp.

127  {
128  int max_order = 0;
129  for (auto e : getDataFieldEnts()) {
130  if (auto ptr = e.lock()) {
131  const int order = ptr->getMaxOrder();
132  max_order = (max_order < order) ? order : max_order;
133  }
134  }
135  return max_order;
136 }
const FieldEntity_vector_view & getDataFieldEnts() const
constexpr int order

◆ getMaxRowOrder()

int MoFEM::ForcesAndSourcesCore::getMaxRowOrder ( ) const

Get max order of approximation for field in rows.

Definition at line 138 of file ForcesAndSourcesCore.cpp.

138  {
139  return getMaxOrder(getRowFieldEnts());
140 }
auto & getRowFieldEnts() const
static int getMaxOrder(const ENTMULTIINDEX &multi_index)

◆ getNodesFieldData()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNodesFieldData ( DataForcesAndSourcesCore data,
const std::string &  field_name 
) const
protected

Get data on nodes.

Parameters
dataData structure
field_nameField name
Returns
Error code

Definition at line 600 of file ForcesAndSourcesCore.cpp.

601  {
602 
603  auto get_nodes_field_data = [&](VectorDouble &nodes_data,
604  VectorFieldEntities &field_entities,
605  VectorDofs &nodes_dofs, FieldSpace &space,
607  VectorInt &bb_node_order) {
609 
610  nodes_data.resize(0, false);
611  nodes_dofs.resize(0, false);
612  field_entities.resize(0, false);
613 
614  auto field_it = fieldsPtr->get<FieldName_mi_tag>().find(field_name);
615  if (field_it != fieldsPtr->get<FieldName_mi_tag>().end()) {
616 
617  auto bit_number = (*field_it)->getBitNumber();
618  const int nb_dofs_on_vert = (*field_it)->getNbOfCoeffs();
619  space = (*field_it)->getSpace();
620  base = (*field_it)->getApproxBase();
621 
622  auto &field_ents = getDataFieldEnts();
623  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
624  bit_number, get_id_for_min_type<MBVERTEX>());
625  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
626  cmp_uid_lo);
627  if (lo != field_ents.end()) {
628  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
629  bit_number, get_id_for_max_type<MBVERTEX>());
630  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
631  if (lo != hi) {
632 
633  int nb_dofs = 0;
634  for (auto it = lo; it != hi; ++it) {
635  if (auto e = it->lock()) {
636  nb_dofs += e->getNbDofsOnEnt();
637  }
638  }
639 
640  if (nb_dofs) {
641 
642  int num_nodes;
643  CHKERR getNumberOfNodes(num_nodes);
644  bb_node_order.resize(num_nodes, false);
645  bb_node_order.clear();
646  const int max_nb_dofs = nb_dofs_on_vert * num_nodes;
647  nodes_data.resize(max_nb_dofs, false);
648  nodes_dofs.resize(max_nb_dofs, false);
649  field_entities.resize(num_nodes, false);
650  std::fill(nodes_data.begin(), nodes_data.end(), 0);
651  std::fill(nodes_dofs.begin(), nodes_dofs.end(), nullptr);
652  std::fill(field_entities.begin(), field_entities.end(), nullptr);
653 
654  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
655 
656  for (auto it = lo; it != hi; ++it) {
657  if (auto e = it->lock()) {
658  const auto &sn = e->getSideNumberPtr();
659  const int side_number = sn->side_number;
660  const int brother_side_number = sn->brother_side_number;
661 
662  field_entities[side_number] = e.get();
663  if (brother_side_number != -1) {
664  brother_ents_vec.emplace_back(e);
665  field_entities[side_number] = field_entities[side_number];
666  }
667 
668  bb_node_order[side_number] = e->getMaxOrder();
669  int pos = side_number * nb_dofs_on_vert;
670  auto ent_filed_data_vec = e->getEntFieldData();
671  if (auto cache = e->entityCacheDataDofs.lock()) {
672  for (auto dit = cache->loHi[0]; dit != cache->loHi[1];
673  ++dit) {
674  const auto dof_idx = (*dit)->getEntDofIdx();
675  nodes_data[pos + dof_idx] = ent_filed_data_vec[dof_idx];
676  nodes_dofs[pos + dof_idx] =
677  reinterpret_cast<FEDofEntity *>((*dit).get());
678  }
679  }
680  }
681  }
682 
683  for (auto &it : brother_ents_vec) {
684  if (const auto e = it.lock()) {
685  const auto &sn = e->getSideNumberPtr();
686  const int side_number = sn->side_number;
687  const int brother_side_number = sn->brother_side_number;
688  bb_node_order[brother_side_number] = bb_node_order[side_number];
689  int pos = side_number * nb_dofs_on_vert;
690  int brother_pos = brother_side_number * nb_dofs_on_vert;
691  for (int ii = 0; ii != nb_dofs_on_vert; ++ii) {
692  nodes_data[brother_pos] = nodes_data[pos];
693  nodes_dofs[brother_pos] = nodes_dofs[pos];
694  ++pos;
695  ++brother_pos;
696  }
697  }
698  }
699  }
700  }
701  }
702  }
703 
705  };
706 
707  return get_nodes_field_data(
708  data.dataOnEntities[MBVERTEX][0].getFieldData(),
709  data.dataOnEntities[MBVERTEX][0].getFieldEntities(),
710  data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
711  data.dataOnEntities[MBVERTEX][0].getSpace(),
712  data.dataOnEntities[MBVERTEX][0].getBase(),
713  data.dataOnEntities[MBVERTEX][0].getBBNodeOrder());
714 }
const Field_multiIndex * fieldsPtr
raw pointer to fields container
BlockParamData * cache
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
const FieldEntity_vector_view & getDataFieldEnts() const
ublas::vector< FieldEntity *, FieldEntAllocator > VectorFieldEntities
static auto cmp_uid_lo(const boost::weak_ptr< FieldEntity > &a, const UId &b)
FieldApproximationBase
approximation base
Definition: definitions.h:150
ublas::vector< FEDofEntity *, DofsAllocator > VectorDofs
static auto cmp_uid_hi(const UId &b, const boost::weak_ptr< FieldEntity > &a)
FieldSpace
approximation spaces
Definition: definitions.h:174
#define CHKERR
Inline error check.
Definition: definitions.h:604
ublas::vector< int, IntAllocator > VectorInt
Definition: Types.hpp:73
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74
MoFEMErrorCode getNumberOfNodes(int &num_nodes) const
Get number of DOFs on element.

◆ getNodesIndices()

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

get node indices

Definition at line 228 of file ForcesAndSourcesCore.cpp.

231  {
233 
234  auto field_it = fieldsPtr->get<FieldName_mi_tag>().find(field_name);
235  if (field_it != fieldsPtr->get<FieldName_mi_tag>().end()) {
236 
237  auto bit_number = (*field_it)->getBitNumber();
238  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
239  bit_number, get_id_for_min_type<MBVERTEX>());
240  auto lo = std::lower_bound(ents_field.begin(), ents_field.end(), lo_uid,
241  cmp_uid_lo);
242  if (lo != ents_field.end()) {
243  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
244  bit_number, get_id_for_max_type<MBVERTEX>());
245  auto hi = std::upper_bound(lo, ents_field.end(), hi_uid, cmp_uid_hi);
246 
247  int num_nodes;
248  CHKERR getNumberOfNodes(num_nodes);
249  const int nb_dofs_on_vert = (*field_it)->getNbOfCoeffs();
250  const int max_nb_dofs = nb_dofs_on_vert * num_nodes;
251 
252  int nb_dofs = 0;
253  for (auto it = lo; it != hi; ++it) {
254  if (auto e = it->lock()) {
255  if (auto cache = extractor(e).lock()) {
256  if (cache->loHi[0] != cache->loHi[1]) {
257  nb_dofs += std::distance(cache->loHi[0], cache->loHi[1]);
258  break;
259  }
260  }
261  }
262  }
263 
264  if (nb_dofs) {
265  nodes_indices.resize(max_nb_dofs, false);
266  local_nodes_indices.resize(max_nb_dofs, false);
267  } else {
268  nodes_indices.resize(0, false);
269  local_nodes_indices.resize(0, false);
270  }
271 
272  if (nb_dofs != max_nb_dofs) {
273  std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
274  std::fill(local_nodes_indices.begin(), local_nodes_indices.end(), -1);
275  }
276 
277  for (auto it = lo; it != hi; ++it) {
278  if (auto e = it->lock()) {
279  auto side_ptr = e->getSideNumberPtr();
280  const auto side_number = side_ptr->side_number;
281  const auto brother_side_number = side_ptr->brother_side_number;
282  if (auto cache = extractor(e).lock()) {
283  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
284  auto &dof = **dit;
285  const int idx = dof.getPetscGlobalDofIdx();
286  const int local_idx = dof.getPetscLocalDofIdx();
287  const int pos =
288  side_number * nb_dofs_on_vert + dof.getDofCoeffIdx();
289  nodes_indices[pos] = idx;
290  local_nodes_indices[pos] = local_idx;
291  if (brother_side_number != -1) {
292  const int elem_idx = brother_side_number * nb_dofs_on_vert +
293  (*dit)->getDofCoeffIdx();
294  nodes_indices[elem_idx] = idx;
295  local_nodes_indices[elem_idx] = local_idx;
296  }
297  }
298  }
299  }
300  }
301  } else {
302  nodes_indices.resize(0, false);
303  local_nodes_indices.resize(0, false);
304  }
305 
306  } else {
307  nodes_indices.resize(0, false);
308  local_nodes_indices.resize(0, false);
309  }
310 
312 }
const Field_multiIndex * fieldsPtr
raw pointer to fields container
BlockParamData * cache
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
static auto cmp_uid_lo(const boost::weak_ptr< FieldEntity > &a, const UId &b)
static auto cmp_uid_hi(const UId &b, const boost::weak_ptr< FieldEntity > &a)
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
MoFEMErrorCode getNumberOfNodes(int &num_nodes) const
Get number of DOFs on element.

◆ getNoFieldColIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldColIndices ( DataForcesAndSourcesCore data,
const std::string &  field_name 
) const
protected

get col NoField indices

Definition at line 477 of file ForcesAndSourcesCore.cpp.

478  {
480  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
481  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
482  }
483  CHKERR getNoFieldIndices(field_name, getColDofsPtr(),
484  data.dataOnEntities[MBENTITYSET][0].getIndices());
486 }
MoFEMErrorCode getNoFieldIndices(const std::string &field_name, boost::shared_ptr< FENumeredDofEntity_multiIndex > dofs, VectorInt &nodes_indices) const
get NoField indices
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define CHKERR
Inline error check.
Definition: definitions.h:604
auto getColDofsPtr() const
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ getNoFieldFieldData() [1/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldFieldData ( const std::string  field_name,
VectorDouble ent_field_data,
VectorDofs ent_field_dofs,
VectorFieldEntities ent_field 
) const
protected

Get field data on nodes.

Definition at line 805 of file ForcesAndSourcesCore.cpp.

807  {
808 
810 
811  ent_field_data.resize(0, false);
812  ent_field_dofs.resize(0, false);
813  ent_field.resize(0, false);
814 
815  auto &field_ents = getDataFieldEnts();
816  auto bit_number = mField.get_field_bit_number(field_name);
817  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
818  bit_number, get_id_for_min_type<MBVERTEX>());
819  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
820  cmp_uid_lo);
821  if (lo != field_ents.end()) {
822 
823  ent_field.resize(field_ents.size(), false);
824  std::fill(ent_field.begin(), ent_field.end(), nullptr);
825 
826  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
827  bit_number, get_id_for_max_type<MBVERTEX>());
828  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
829  if (lo != hi) {
830 
831  int side = 0;
832  for (auto it = lo; it != hi; ++it, ++side)
833  if (auto e = it->lock()) {
834 
835  const auto size = e->getNbDofsOnEnt();
836  ent_field_data.resize(size, false);
837  ent_field_dofs.resize(size, false);
838  ent_field[side] = e.get();
839  noalias(ent_field_data) = e->getEntFieldData();
840 
841 
842  if (auto cache = e->entityCacheDataDofs.lock()) {
843  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
844  ent_field_dofs[(*dit)->getEntDofIdx()] =
845  reinterpret_cast<FEDofEntity *>((*dit).get());
846  }
847  }
848  }
849  }
850  }
851 
853 }
BlockParamData * cache
UId getLocalUniqueIdCalculate()
Get the Local Unique Id Calculate object.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
const FieldEntity_vector_view & getDataFieldEnts() const
static auto cmp_uid_lo(const boost::weak_ptr< FieldEntity > &a, const UId &b)
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
static auto cmp_uid_hi(const UId &b, const boost::weak_ptr< FieldEntity > &a)

◆ getNoFieldFieldData() [2/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldFieldData ( DataForcesAndSourcesCore data,
const std::string  field_name 
) const
protected

Definition at line 856 of file ForcesAndSourcesCore.cpp.

857  {
859  if (data.dataOnEntities[MBENTITYSET].size() == 0)
860  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
861  "No space to insert data");
862 
864  field_name, data.dataOnEntities[MBENTITYSET][0].getFieldData(),
865  data.dataOnEntities[MBENTITYSET][0].getFieldDofs(),
866  data.dataOnEntities[MBENTITYSET][0].getFieldEntities());
868 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
MoFEMErrorCode getNoFieldFieldData(const std::string field_name, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs, VectorFieldEntities &ent_field) const
Get field data on nodes.
#define CHKERR
Inline error check.
Definition: definitions.h:604
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ getNoFieldIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldIndices ( const std::string &  field_name,
boost::shared_ptr< FENumeredDofEntity_multiIndex dofs,
VectorInt nodes_indices 
) const
protected

get NoField indices

Definition at line 448 of file ForcesAndSourcesCore.cpp.

451  {
453  auto field_it = fieldsPtr->get<FieldName_mi_tag>().find(field_name);
454  auto dit = dofs->get<Unique_mi_tag>().lower_bound(
455  FieldEntity::getLoBitNumberUId((*field_it)->getBitNumber()));
456  auto hi_dit = dofs->get<Unique_mi_tag>().upper_bound(
457  FieldEntity::getHiBitNumberUId((*field_it)->getBitNumber()));
458  indices.resize(std::distance(dit, hi_dit));
459  for (; dit != hi_dit; dit++) {
460  int idx = (*dit)->getPetscGlobalDofIdx();
461  indices[(*dit)->getDofCoeffIdx()] = idx;
462  }
464 }
const Field_multiIndex * fieldsPtr
raw pointer to fields container
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
static UId getLoBitNumberUId(const FieldBitNumber bit_number)
static UId getHiBitNumberUId(const FieldBitNumber bit_number)

◆ getNoFieldRowIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldRowIndices ( DataForcesAndSourcesCore data,
const std::string &  field_name 
) const
protected

get col NoField indices

Definition at line 466 of file ForcesAndSourcesCore.cpp.

467  {
469  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
470  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
471  }
472  CHKERR getNoFieldIndices(field_name, getRowDofsPtr(),
473  data.dataOnEntities[MBENTITYSET][0].getIndices());
475 }
MoFEMErrorCode getNoFieldIndices(const std::string &field_name, boost::shared_ptr< FENumeredDofEntity_multiIndex > dofs, VectorInt &nodes_indices) const
get NoField indices
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
#define CHKERR
Inline error check.
Definition: definitions.h:604
auto getRowDofsPtr() const
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415

◆ getOpPtrVector()

boost::ptr_vector<UserDataOperator>& MoFEM::ForcesAndSourcesCore::getOpPtrVector ( )

Use to push back operator for row operator.

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

Examples
continuity_check_on_contact_prism_side_ele.cpp, continuity_check_on_skeleton_3d.cpp, continuity_check_on_skeleton_with_simple_2d.cpp, forces_and_sources_testing_edge_element.cpp, forces_and_sources_testing_flat_prism_element.cpp, and lorentz_force.cpp.

Definition at line 423 of file ForcesAndSourcesCore.hpp.

423 { return opPtrVector; }
boost::ptr_vector< 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 583 of file ForcesAndSourcesCore.cpp.

584  {
585  return getProblemNodesIndices(field_name, *(problemPtr->numeredColDofsPtr),
586  nodes_indices);
587 }
const Problem * problemPtr
raw pointer to problem
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredColDofsPtr
store DOFs on columns for this 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

◆ 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  int num_nodes;
499  CHKERR getNumberOfNodes(num_nodes);
500  nodes_indices.resize(field_struture->getNbOfCoeffs() * num_nodes, false);
501  std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
502 
503  auto &side_table = const_cast<SideNumber_multiIndex &>(
504  numeredEntFiniteElementPtr->getSideNumberTable());
505  auto siit = side_table.get<1>().lower_bound(boost::make_tuple(MBVERTEX, 0));
506  auto hi_siit =
507  side_table.get<1>().lower_bound(boost::make_tuple(MBVERTEX, 10000));
508 
509  int nn = 0;
510  for (; siit != hi_siit; siit++, nn++) {
511 
512  if (siit->get()->side_number == -1)
513  continue;
514 
515  auto bit_number = mField.get_field_bit_number(field_name);
516  const EntityHandle ent = siit->get()->ent;
517  auto dit = dofs.get<Unique_mi_tag>().lower_bound(
518  FieldEntity::getLoLocalEntityBitNumber(bit_number, ent));
519  auto hi_dit = dofs.get<Unique_mi_tag>().upper_bound(
520  FieldEntity::getHiLocalEntityBitNumber(bit_number, ent));
521  for (; dit != hi_dit; dit++) {
522  nodes_indices[siit->get()->side_number * (*dit)->getNbOfCoeffs() +
523  (*dit)->getDofCoeffIdx()] =
524  (*dit)->getPetscGlobalDofIdx();
525  }
526  }
527  } else {
528  nodes_indices.resize(0, false);
529  }
530 
532 }
static UId getHiLocalEntityBitNumber(const char bit_number, const EntityHandle ent)
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
static UId getLoLocalEntityBitNumber(const char bit_number, const EntityHandle ent)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number
#define CHKERR
Inline error check.
Definition: definitions.h:604
continuous field
Definition: definitions.h:177
virtual const Field * get_field_structure(const std::string &name)=0
get field structure
MoFEMErrorCode getNumberOfNodes(int &num_nodes) const
Get number of DOFs on element.

◆ getProblemNodesRowIndices()

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

Definition at line 569 of file ForcesAndSourcesCore.cpp.

570  {
571  return getProblemNodesIndices(field_name, *(problemPtr->numeredRowDofsPtr),
572  nodes_indices);
573 }
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredRowDofsPtr
store DOFs on rows for this problem
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

◆ getProblemTypeColIndices()

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

Definition at line 590 of file ForcesAndSourcesCore.cpp.

592  {
593  return getProblemTypeIndices(field_name, *(problemPtr->numeredColDofsPtr), type,
594  side_number, indices);
595 }
const Problem * problemPtr
raw pointer to problem
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
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredColDofsPtr
store DOFs on columns for this problem

◆ 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 534 of file ForcesAndSourcesCore.cpp.

536  {
538 
539  indices.resize(0);
540 
541  auto &side_table = const_cast<SideNumber_multiIndex &>(
542  numeredEntFiniteElementPtr->getSideNumberTable());
543  auto siit =
544  side_table.get<1>().lower_bound(boost::make_tuple(type, side_number));
545  auto hi_siit =
546  side_table.get<1>().upper_bound(boost::make_tuple(type, side_number));
547 
548  for (; siit != hi_siit; siit++) {
549 
550  if (siit->get()->side_number == -1)
551  continue;
552 
553  const EntityHandle ent = siit->get()->ent;
554  auto bit_number = mField.get_field_bit_number(field_name);
555  auto dit = dofs.get<Unique_mi_tag>().lower_bound(
556  FieldEntity::getLoLocalEntityBitNumber(bit_number, ent));
557  auto hi_dit = dofs.get<Unique_mi_tag>().upper_bound(
558  FieldEntity::getHiLocalEntityBitNumber(bit_number, ent));
559  indices.resize(std::distance(dit, hi_dit));
560  for (; dit != hi_dit; dit++) {
561 
562  indices[(*dit)->getEntDofIdx()] = (*dit)->getPetscGlobalDofIdx();
563  }
564  }
565 
567 }
static UId getHiLocalEntityBitNumber(const char bit_number, const EntityHandle ent)
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
static UId getLoLocalEntityBitNumber(const char bit_number, const EntityHandle ent)
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number

◆ getProblemTypeRowIndices()

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

Definition at line 576 of file ForcesAndSourcesCore.cpp.

578  {
579  return getProblemTypeIndices(field_name, *(problemPtr->numeredRowDofsPtr), type,
580  side_number, indices);
581 }
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredRowDofsPtr
store DOFs on rows for this problem
const Problem * problemPtr
raw pointer to problem
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

◆ getRowNodesIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getRowNodesIndices ( DataForcesAndSourcesCore data,
const std::string &  field_name 
) const
protected

get row node indices from FENumeredDofEntity_multiIndex

Definition at line 315 of file ForcesAndSourcesCore.cpp.

316  {
317 
318  struct Extractor {
319  boost::weak_ptr<EntityCacheNumeredDofs>
320  operator()(boost::shared_ptr<FieldEntity> &e) {
321  return e->entityCacheRowDofs;
322  }
323  };
324 
325  return getNodesIndices(field_name, getRowFieldEnts(),
326  data.dataOnEntities[MBVERTEX][0].getIndices(),
327  data.dataOnEntities[MBVERTEX][0].getLocalIndices(),
328  Extractor());
329 }
auto & getRowFieldEnts() const
MoFEMErrorCode getNodesIndices(const std::string &field_name, FieldEntity_vector_view &ents_field, VectorInt &nodes_indices, VectorInt &local_nodes_indices, EXTRACTOR &&extractor) const
get node indices
virtual MoFEMErrorCode operator()()
function is run for every finite element

◆ getRule() [1/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); };

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 EdgeFE.

Definition at line 1580 of file ForcesAndSourcesCore.cpp.

1581  {
1582  return getRuleHook ? getRuleHook(order_row, order_col, order_data)
1583  : getRule(order_data);
1584 }
virtual int getRule(int order_row, int order_col, int order_data)
another variant of getRule
RuleHookFun getRuleHook
Hook to get rule.

◆ getRule() [2/2]

int MoFEM::ForcesAndSourcesCore::getRule ( int  order)
protectedvirtual

◆ getSpacesAndBaseOnEntities()

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

Get field approximation space and base on entities.

Definition at line 967 of file ForcesAndSourcesCore.cpp.

968  {
970 
971  if (nInTheLoop == 0) {
972  data.sPace.reset();
973  data.bAse.reset();
974  for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
975  data.spacesOnEntities[t].reset();
976  data.basesOnEntities[t].reset();
977  }
978  for (int s = 0; s != LASTSPACE; ++s) {
979  data.basesOnSpaces[s].reset();
980  }
981  }
982 
983  auto fe_type = numeredEntFiniteElementPtr->getEntType();
984 
985  if (getDataFieldEntsPtr())
986  for (auto e : getDataFieldEnts()) {
987  if (auto ptr = e.lock()) {
988  // get data from entity
989  const EntityType type = ptr->getEntType();
991  const FieldSpace space = ptr->getSpace();
992  const FieldApproximationBase approx = ptr->getApproxBase();
993  // set data
994  data.sPace.set(space);
995  data.bAse.set(approx);
996  data.spacesOnEntities[type].set(space);
997  data.basesOnEntities[type].set(approx);
998  data.basesOnSpaces[space].set(approx);
999  }
1000  }
1001  }
1002  else
1003  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1004  "data fields ents not allocated on element");
1005 
1007 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:181
const FieldEntity_vector_view & getDataFieldEnts() const
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
FieldApproximationBase
approximation base
Definition: definitions.h:150
FieldSpace
approximation spaces
Definition: definitions.h:174
int nInTheLoop
number currently of processed method
boost::shared_ptr< FieldEntity_vector_view > & getDataFieldEntsPtr() const
static bool getDefTypeMap(const EntityType fe_type, const EntityType ent_type)

◆ getUserPolynomialBase()

auto& MoFEM::ForcesAndSourcesCore::getUserPolynomialBase ( )

Get the User Polynomial Base object.

Returns
boost::shared_ptr<BaseFunction>&

Definition at line 437 of file ForcesAndSourcesCore.hpp.

437 { 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 1296 of file ForcesAndSourcesCore.cpp.

1296  {
1298 
1299  const EntityType type = numeredEntFiniteElementPtr->getEntType();
1302  std::vector<std::string> last_eval_field_name(2);
1303 
1304  boost::ptr_vector<UserDataOperator>::iterator oit, hi_oit;
1305  oit = opPtrVector.begin();
1306  hi_oit = opPtrVector.end();
1307 
1308  for (; oit != hi_oit; oit++) {
1309 
1310  try {
1311 
1312  CHKERR oit->setPtrFE(this);
1313 
1314  if (oit->opType == UserDataOperator::OPLAST) {
1315 
1316  // Set field
1317  switch (oit->sPace) {
1318  case NOSPACE:
1319  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Unknown space");
1320  case NOFIELD:
1321  case H1:
1322  case HCURL:
1323  case HDIV:
1324  case L2:
1325  break;
1326  default:
1327  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1328  "Not implemented for this space", oit->sPace);
1329  }
1330 
1331  // Reseat all data which all field dependent
1332  dataOnElement[oit->sPace]->resetFieldDependentData();
1333  last_eval_field_name[0] = "";
1334 
1335  // Run operator
1336  try {
1337  CHKERR oit->opRhs(*dataOnElement[oit->sPace], false);
1338  }
1339  CATCH_OP_ERRORS(*oit);
1340 
1341  } else {
1342 
1343  boost::shared_ptr<DataForcesAndSourcesCore> op_data[2];
1344  std::array<bool, 2> base_swap;
1345  std::array<std::pair<std::string, FieldApproximationBase>, 2>
1346  base_swap_data;
1347  auto swap_bases = [&]() {
1349  for (size_t ss = 0; ss != 2; ++ss)
1350  if (base_swap[ss])
1351  CHKERR op_data[ss]->baseSwap(base_swap_data[ss].first,
1352  base_swap_data[ss].second);
1354  };
1355 
1356  for (size_t ss = 0; ss != 2; ss++) {
1357 
1358  const std::string field_name =
1359  !ss ? oit->rowFieldName : oit->colFieldName;
1360  if (field_name.empty()) {
1361  SETERRQ2(
1362  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1363  "No field name in operator %d (0-row, 1-column) in operator %s",
1364  ss,
1365  (boost::typeindex::type_id_runtime(*oit).pretty_name())
1366  .c_str());
1367  }
1368  const Field *field_struture = mField.get_field_structure(field_name);
1369  const BitFieldId data_id = field_struture->getId();
1370  const FieldSpace space = field_struture->getSpace();
1371  const FieldApproximationBase base = field_struture->getApproxBase();
1372  op_data[ss] =
1373  !ss ? dataOnElement[space] : derivedDataOnElement[space];
1374 
1375  switch (base) {
1377  base_swap_data[ss] = std::pair<std::string, FieldApproximationBase>(
1378  field_name, AINSWORTH_BERNSTEIN_BEZIER_BASE);
1379  base_swap[ss] = true;
1380  break;
1381  default:
1382  base_swap[ss] = false;
1383  };
1384 
1385  if ((oit->getNumeredEntFiniteElementPtr()->getBitFieldIdData() &
1386  data_id)
1387  .none()) {
1388  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1389  "no data field < %s > on finite element < %s >",
1390  field_name.c_str(), feName.c_str());
1391  }
1392 
1393  if (oit->getOpType() & types[ss] ||
1394  oit->getOpType() & UserDataOperator::OPROWCOL) {
1395 
1396  switch (space) {
1397  case NOSPACE:
1398  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1399  "unknown space");
1400  break;
1401  case NOFIELD:
1402  case H1:
1403  case HCURL:
1404  case HDIV:
1405  case L2:
1406  break;
1407  default:
1408  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1409  "Not implemented for this space", space);
1410  }
1411 
1412  if (last_eval_field_name[ss] != field_name) {
1413 
1414  CHKERR getEntityFieldData(*op_data[ss], field_name, MBEDGE);
1415  if (!ss)
1416  CHKERR getEntityRowIndices(*op_data[ss], field_name, MBEDGE);
1417  else
1418  CHKERR getEntityColIndices(*op_data[ss], field_name, MBEDGE);
1419 
1420  switch (space) {
1421  case NOSPACE:
1422  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1423  "unknown space");
1424  break;
1425  case H1:
1426  if (!ss)
1427  CHKERR getRowNodesIndices(*op_data[ss], field_name);
1428  else
1429  CHKERR getColNodesIndices(*op_data[ss], field_name);
1430  CHKERR getNodesFieldData(*op_data[ss], field_name);
1431  break;
1432  case HCURL:
1433  case HDIV:
1434  break;
1435  case L2:
1436  switch (type) {
1437  case MBVERTEX:
1438  CHKERR getNodesFieldData(*op_data[ss], field_name);
1439  break;
1440  default:
1441  break;
1442  }
1443  break;
1444  case NOFIELD:
1445  if (!getNinTheLoop()) {
1446  // NOFIELD data are the same for each element, can be
1447  // retrieved only once
1448  if (!ss) {
1449  CHKERR getNoFieldRowIndices(*op_data[ss], field_name);
1450  } else {
1451  CHKERR getNoFieldColIndices(*op_data[ss], field_name);
1452  }
1453  CHKERR getNoFieldFieldData(*op_data[ss], field_name);
1454  }
1455  break;
1456  default:
1457  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1458  "Not implemented for this space", space);
1459  }
1460  last_eval_field_name[ss] = field_name;
1461  }
1462  }
1463  }
1464 
1465  CHKERR swap_bases();
1466 
1467  if (oit->getOpType() & UserDataOperator::OPROW) {
1468  try {
1469  CHKERR oit->opRhs(*op_data[0], false);
1470  }
1471  CATCH_OP_ERRORS(*oit);
1472  }
1473 
1474  if (oit->getOpType() & UserDataOperator::OPCOL) {
1475  try {
1476  CHKERR oit->opRhs(*op_data[1], false);
1477  }
1478  CATCH_OP_ERRORS(*oit);
1479  }
1480 
1481  if (oit->getOpType() & UserDataOperator::OPROWCOL) {
1482  try {
1483  CHKERR oit->opLhs(*op_data[0], *op_data[1]);
1484  }
1485  CATCH_OP_ERRORS(*oit);
1486  }
1487 
1488  CHKERR swap_bases();
1489  }
1490  }
1491  CATCH_OP_ERRORS(*oit);
1492  }
1494 }
boost::ptr_vector< UserDataOperator > opPtrVector
Vector of finite element users data operators.
field with continuous normal traction
Definition: definitions.h:179
MoFEMErrorCode getNoFieldColIndices(DataForcesAndSourcesCore &data, const std::string &field_name) const
get col NoField indices
MoFEMErrorCode getColNodesIndices(DataForcesAndSourcesCore &data, const std::string &field_name) const
get col node indices from FENumeredDofEntity_multiIndex
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:485
std::string feName
Name of finite element.
int getNinTheLoop() const
get number of evaluated element in the loop
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
MoFEMErrorCode getNoFieldRowIndices(DataForcesAndSourcesCore &data, const std::string &field_name) const
get col NoField indices
OpType
Controls loop over entities on element.
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
MoFEMErrorCode getNodesFieldData(DataForcesAndSourcesCore &data, const std::string &field_name) const
Get data on nodes.
#define CATCH_OP_ERRORS(OP)
MoFEMErrorCode getNoFieldFieldData(const std::string field_name, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs, VectorFieldEntities &ent_field) const
Get field data on nodes.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnElement
Entity data on element entity columns fields.
FieldApproximationBase
approximation base
Definition: definitions.h:150
field with continuous tangents
Definition: definitions.h:178
FieldSpace
approximation spaces
Definition: definitions.h:174
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
Definition: Types.hpp:53
#define CHKERR
Inline error check.
Definition: definitions.h:604
MoFEMErrorCode getEntityFieldData(DataForcesAndSourcesCore &data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
MoFEMErrorCode getEntityColIndices(DataForcesAndSourcesCore &data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:415
continuous field
Definition: definitions.h:177
virtual const Field * get_field_structure(const std::string &name)=0
get field structure
MoFEMErrorCode getRowNodesIndices(DataForcesAndSourcesCore &data, const std::string &field_name) const
get row node indices from FENumeredDofEntity_multiIndex
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
MoFEMErrorCode getEntityRowIndices(DataForcesAndSourcesCore &data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
field with C-1 continuity
Definition: definitions.h:180

◆ 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::VolumeElementForcesAndSourcesCoreSwitch< SWITCH >, MoFEM::FaceElementForcesAndSourcesCoreSwitch< SWITCH >, MoFEM::ContactPrismElementForcesAndSourcesCore, MoFEM::VolumeElementForcesAndSourcesCoreOnSideSwitch< SWITCH >, MoFEM::FatPrismElementForcesAndSourcesCore, MoFEM::EdgeElementForcesAndSourcesCoreSwitch< SWITCH >, MoFEM::FaceElementForcesAndSourcesCoreOnSideSwitch< SWITCH >, MoFEM::FaceElementForcesAndSourcesCoreOnSideSwitch< FaceElementForcesAndSourcesCore::NO_HO_GEOMETRY >, MoFEM::FaceElementForcesAndSourcesCoreOnSideSwitch< FaceElementForcesAndSourcesCore::NO_CONTRAVARIANT_TRANSFORM_HDIV|FaceElementForcesAndSourcesCore::NO_COVARIANT_TRANSFORM_HCURL >, MoFEM::FlatPrismElementForcesAndSourcesCore, MoFEM::VolumeElementForcesAndSourcesCoreOnContactPrismSideSwitch< SWITCH >, and MoFEM::VertexElementForcesAndSourcesCore.

Definition at line 1626 of file ForcesAndSourcesCore.cpp.

1626  {
1628  if (operatorHook) {
1629  ierr = operatorHook();
1630  CHKERRG(ierr);
1631  }
1633 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:552
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
boost::function< MoFEMErrorCode()> operatorHook
Hook function for operator.
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87

◆ 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 PostProcEdgeOnRefinedMesh, PostProcFatPrismOnRefinedMesh, SolidShellModule::SolidShellPrismElement::PostProcFatPrismOnTriangleOnRefinedMesh, PostProcTemplateVolumeOnRefinedMesh< MoFEM::VolumeElementForcesAndSourcesCore >, GelModule::Gel::GelFE, KelvinVoigtDamper::DamperFE, SolidShellModule::SolidShellPrismElement::SolidShellError, SolidShellModule::SolidShellPrismElement::SolidShell, Smoother::MyVolumeFE, and NonlinearElasticElement::MyVolumeFE.

Definition at line 1634 of file ForcesAndSourcesCore.cpp.

1634  {
1636  if (postProcessHook) {
1637  ierr = postProcessHook();
1638  CHKERRG(ierr);
1639  }
1641 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:552
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
boost::function< MoFEMErrorCode()> postProcessHook
Hook function for post-processing.
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87

◆ preProcess()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::preProcess ( )
virtual

function is run at the beginning of loop

It is used to zeroing matrices and vectors, calculation of shape functions on reference element, preprocessing boundary conditions, etc.

Reimplemented from MoFEM::BasicMethod.

Reimplemented in EdgeSlidingConstrains::MyEdgeFE, PostProcEdgeOnRefinedMesh, PostProcFatPrismOnRefinedMesh, PostProcTemplateVolumeOnRefinedMesh< MoFEM::VolumeElementForcesAndSourcesCore >, GelModule::Gel::GelFE, SurfaceSlidingConstrains::MyTriangleFE, KelvinVoigtDamper::DamperFE, SolidShellModule::SolidShellPrismElement::SolidShellError, SolidShellModule::SolidShellPrismElement::SolidShell, Smoother::MyVolumeFE, and NonlinearElasticElement::MyVolumeFE.

Definition at line 1618 of file ForcesAndSourcesCore.cpp.

1618  {
1620  if (preProcessHook) {
1621  ierr = preProcessHook();
1622  CHKERRG(ierr);
1623  }
1625 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:552
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
boost::function< MoFEMErrorCode()> preProcessHook
Hook function for pre-processing.

◆ setGaussPts() [1/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);

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

Reimplemented in EdgeFE.

Definition at line 1586 of file ForcesAndSourcesCore.cpp.

1587  {
1588  return setRuleHook ? setRuleHook(order_row, order_col, order_data)
1589  : setGaussPts(order_data);
1590 }
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule
RuleHookFun setRuleHook
Set function to calculate integration rule.

◆ setGaussPts() [2/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::setGaussPts ( int  order)
protectedvirtual
Deprecated:
setGaussPts(int row_order, int col_order, int data order);
Deprecated:
setGaussPts(int row_order, int col_order, int data order);

Reimplemented in PostProcEdgeOnRefinedMesh, PostProcTemplateVolumeOnRefinedMesh< MoFEM::VolumeElementForcesAndSourcesCore >, NitscheMethod::MyVolumeFE, MoFEM::VolumeElementForcesAndSourcesCoreOnSideBase, MoFEM::FaceElementForcesAndSourcesCoreOnSideBase, and MoFEM::VolumeElementForcesAndSourcesCoreOnContactPrismSideBase.

Definition at line 1596 of file ForcesAndSourcesCore.cpp.

1596  {
1598  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Sorry, not implemented");
1600 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516

◆ 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 1531 of file ForcesAndSourcesCore.cpp.

1531  {
1533  sidePtrFE = const_cast<ForcesAndSourcesCore *>(side_fe_ptr);
1535 }
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:516

Friends And Related Function Documentation

◆ FaceElementForcesAndSourcesCoreOnSideBase

Definition at line 855 of file ForcesAndSourcesCore.hpp.

◆ UserDataOperator

friend class UserDataOperator
friend

◆ VolumeElementForcesAndSourcesCoreOnContactPrismSideBase

Definition at line 856 of file ForcesAndSourcesCore.hpp.

◆ VolumeElementForcesAndSourcesCoreOnSideBase

Definition at line 854 of file ForcesAndSourcesCore.hpp.

Member Data Documentation

◆ dataH1

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataH1
protected

Definition at line 805 of file ForcesAndSourcesCore.hpp.

◆ dataHcurl

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataHcurl
protected

Definition at line 806 of file ForcesAndSourcesCore.hpp.

◆ dataHdiv

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataHdiv
protected

Definition at line 807 of file ForcesAndSourcesCore.hpp.

◆ dataL2

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataL2
protected

Definition at line 808 of file ForcesAndSourcesCore.hpp.

◆ dataNoField

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataNoField
protected

Definition at line 804 of file ForcesAndSourcesCore.hpp.

◆ dataOnElement

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

Entity data on element entity rows fields.

@/}

Definition at line 795 of file ForcesAndSourcesCore.hpp.

◆ derivedDataOnElement

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

Entity data on element entity columns fields.

Definition at line 802 of file ForcesAndSourcesCore.hpp.

◆ elementPolynomialBasePtr

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

Pointer to entity polynomial base.

Definition at line 830 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
prism_elements_from_surface.cpp.

Definition at line 449 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 50 of file ForcesAndSourcesCore.hpp.

◆ lastEvaluatedElementEntityType

EntityType MoFEM::ForcesAndSourcesCore::lastEvaluatedElementEntityType
protected

Last evaluated type of element entity.

Definition at line 823 of file ForcesAndSourcesCore.hpp.

◆ mField

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

Definition at line 36 of file ForcesAndSourcesCore.hpp.

◆ opPtrVector

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

Vector of finite element users data operators.

Definition at line 814 of file ForcesAndSourcesCore.hpp.

◆ setRuleHook

RuleHookFun MoFEM::ForcesAndSourcesCore::setRuleHook

Set function to calculate integration rule.

Definition at line 56 of file ForcesAndSourcesCore.hpp.

◆ sidePtrFE

ForcesAndSourcesCore* MoFEM::ForcesAndSourcesCore::sidePtrFE
private

Element to integrate on the sides.

Definition at line 841 of file ForcesAndSourcesCore.hpp.

◆ userPolynomialBasePtr

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

Pointer to user polynomail base.

Definition at line 835 of file ForcesAndSourcesCore.hpp.


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