v0.11.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. 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_t (see PETSc Time Solver) More...
 
PetscReal ts_aa
 shift for U_tt shift for U_tt 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

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  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // NOSPACE,
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) {}
@ L2
field with C-1 continuity
Definition: definitions.h:180
@ NOFIELD
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
@ H1
continuous field
Definition: definitions.h:177
@ HCURL
field with continuous tangents
Definition: definitions.h:178
@ HDIV
field with continuous normal traction
Definition: definitions.h:179
DataForcesAndSourcesCore & dataNoField
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnElement
Entity data on element entity columns fields.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.
DataForcesAndSourcesCore & dataH1
DataForcesAndSourcesCore & dataHdiv
DataForcesAndSourcesCore & dataHcurl
RuleHookFun setRuleHook
Set function to calculate integration rule.
EntityType lastEvaluatedElementEntityType
Last evaluated type of element entity.
DataForcesAndSourcesCore & dataL2
RuleHookFun getRuleHook
Hook to get rule.

Member Function Documentation

◆ calBernsteinBezierBaseFunctionsOnElement()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::calBernsteinBezierBaseFunctionsOnElement ( )
protected

Calculate Bernstein-Bezier base.

Returns
MoFEMErrorCode

Definition at line 1087 of file ForcesAndSourcesCore.cpp.

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

◆ calHierarchicalBaseFunctionsOnElement() [1/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::calHierarchicalBaseFunctionsOnElement ( )
protected

Calculate base functions.

Returns
Error code

Definition at line 1067 of file ForcesAndSourcesCore.cpp.

1067  {
1069  /// Use the some node base. Node base is usually used for construction other
1070  /// bases.
1071  for (int space = HCURL; space != LASTSPACE; ++space) {
1072  dataOnElement[space]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
1073  dataOnElement[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
1074  dataOnElement[space]->dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(
1075  NOBASE) =
1076  dataOnElement[H1]->dataOnEntities[MBVERTEX][0].getDiffNSharedPtr(
1077  NOBASE);
1078  }
1079  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
1081  static_cast<FieldApproximationBase>(b));
1082  }
1084 }
FieldApproximationBase
approximation base
Definition: definitions.h:150
@ LASTBASE
Definition: definitions.h:161
@ AINSWORTH_LEGENDRE_BASE
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:152
@ LASTSPACE
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:181
MoFEMErrorCode calHierarchicalBaseFunctionsOnElement()
Calculate base functions.

◆ calHierarchicalBaseFunctionsOnElement() [2/2]

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

Calculate base functions.

Returns
Error code

Definition at line 1013 of file ForcesAndSourcesCore.cpp.

1014  {
1016  if (dataOnElement[H1]->bAse.test(b)) {
1017  switch (static_cast<FieldApproximationBase>(b)) {
1018  case NOBASE:
1019  break;
1021  // BERNSTEIN_BEZIER_BASE is not hierarchical base
1022  break;
1025  case DEMKOWICZ_JACOBI_BASE:
1026  if (!getElementPolynomialBase())
1027  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1028  "Functions genrating approximation base not defined");
1029 
1030  for (int space = H1; space != LASTSPACE; ++space) {
1031  if (dataOnElement[H1]->sPace.test(space) &&
1032  dataOnElement[H1]->bAse.test(b) &&
1033  dataOnElement[H1]->basesOnSpaces[space].test(b)) {
1034  CHKERR getElementPolynomialBase()->getValue(
1035  gaussPts,
1036  boost::make_shared<EntPolynomialBaseCtx>(
1037  *dataOnElement[space], static_cast<FieldSpace>(space),
1038  static_cast<FieldApproximationBase>(b), NOBASE));
1039  }
1040  }
1041  break;
1042  case USER_BASE:
1043  if (!getUserPolynomialBase())
1044  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1045  "Functions genrating user approximation base not defined");
1046 
1047  for (int space = H1; space != LASTSPACE; ++space)
1048  if (dataOnElement[H1]->sPace.test(space) &&
1049  dataOnElement[H1]->bAse.test(b) &&
1050  dataOnElement[H1]->basesOnSpaces[space].test(b)) {
1051  CHKERR getUserPolynomialBase()->getValue(
1052  gaussPts,
1053  boost::make_shared<EntPolynomialBaseCtx>(
1054  *dataOnElement[space], static_cast<FieldSpace>(space),
1055  static_cast<FieldApproximationBase>(b), NOBASE));
1056  }
1057  break;
1058  default:
1060  "Base <%s> not yet implemented",
1062  }
1063  }
1065 }
@ AINSWORTH_LOBATTO_BASE
Definition: definitions.h:154
@ USER_BASE
user implemented approximation base
Definition: definitions.h:160
@ DEMKOWICZ_JACOBI_BASE
Definition: definitions.h:158
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:123
static const char *const ApproximationBaseNames[]
Definition: definitions.h:164
virtual MPI_Comm & get_comm() const =0
auto & getUserPolynomialBase()
Get the User Polynomial Base object.

◆ createDataOnElement()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::createDataOnElement ( )
protected

Create a entity data on element object.

Returns
MoFEMErrorCode

Definition at line 1256 of file ForcesAndSourcesCore.cpp.

1256  {
1258 
1259  const EntityType type = numeredEntFiniteElementPtr->getEntType();
1262 
1263  // Data on elements for proper spaces
1264  for (int space = H1; space != LASTSPACE; ++space) {
1265  dataOnElement[space]->setElementType(type);
1266  derivedDataOnElement[space]->setElementType(type);
1267  }
1268 
1270 
1272 }
#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

◆ getColNodesIndices()

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

get col node indices from FENumeredDofEntity_multiIndex

Definition at line 334 of file ForcesAndSourcesCore.cpp.

335  {
336 
337  struct Extractor {
338  boost::weak_ptr<EntityCacheNumeredDofs>
339  operator()(boost::shared_ptr<FieldEntity> &e) {
340  return e->entityCacheColDofs;
341  }
342  };
343 
344  return getNodesIndices(field_name, getColFieldEnts(),
345  data.dataOnEntities[MBVERTEX][0].getIndices(),
346  data.dataOnEntities[MBVERTEX][0].getLocalIndices(),
347  Extractor());
348 }
auto & getColFieldEnts() const
virtual MoFEMErrorCode operator()()
function is run for every finite element
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

◆ getElementPolynomialBase()

auto& MoFEM::ForcesAndSourcesCore::getElementPolynomialBase ( )

Get the Entity Polynomial Base object.

Returns
boost::shared_ptr<BaseFunction>&&

Definition at line 432 of file ForcesAndSourcesCore.hpp.

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

◆ getEntData() [1/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 501 of file ForcesAndSourcesCore.hpp.

501  {
502  return dataOnElement[space]->dataOnEntities[type][side];
503  }

◆ getEntData() [2/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 486 of file ForcesAndSourcesCore.hpp.

488  {
489  return dataOnElement[space]->dataOnEntities[type][side];
490  }

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

437  {
438 
439  struct Extractor {
440  boost::weak_ptr<EntityCacheNumeredDofs>
441  operator()(boost::shared_ptr<FieldEntity> &e) {
442  return e->entityCacheColDofs;
443  }
444  };
445 
446  return getEntityIndices(data, field_name, getColFieldEnts(), type_lo, type_hi,
447  Extractor());
448 }
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

◆ 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 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:509
int ApproximationOrder
Approximation on the entity.
Definition: Types.hpp:37
std::vector< boost::weak_ptr< FieldEntity > > FieldEntity_vector_view
const double r
rate factor
const Field_multiIndex * fieldsPtr
raw pointer to fields container

◆ 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 549 of file ForcesAndSourcesCore.hpp.

550  {
551  return getEntityDataOrder(type, space, data.dataOnEntities[type]);
552  }
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 722 of file ForcesAndSourcesCore.cpp.

724  {
726  for (EntityType t = type_lo; t != type_hi; ++t) {
727  for (auto &dat : data.dataOnEntities[t]) {
728  dat.getDataOrder() = 0;
729  dat.getBase() = NOBASE;
730  dat.getSpace() = NOSPACE;
731  dat.getFieldData().resize(0, false);
732  dat.getFieldDofs().resize(0, false);
733  dat.getFieldEntities().resize(0, false);
734  }
735  }
736 
737  auto &field_ents = getDataFieldEnts();
738  auto bit_number = mField.get_field_bit_number(field_name);
739  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
740  bit_number, get_id_for_min_type(type_lo));
741  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
742  cmp_uid_lo);
743  if (lo != field_ents.end()) {
744  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
745  bit_number, get_id_for_max_type(type_hi));
746  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
747  if (lo != hi) {
748 
749  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
750 
751  for (auto it = lo; it != hi; ++it)
752  if (auto e = it->lock()) {
753 
754  const EntityType type = e->getEntType();
755  auto side_ptr = e->getSideNumberPtr();
756  const int side = side_ptr->side_number;
757  if (side >= 0) {
758 
759  auto &dat = data.dataOnEntities[type][side];
760  auto &ent_field = dat.getFieldEntities();
761  auto &ent_field_dofs = dat.getFieldDofs();
762  auto &ent_field_data = dat.getFieldData();
763 
764  const int brother_side = side_ptr->brother_side_number;
765  if (brother_side != -1)
766  brother_ents_vec.emplace_back(e);
767 
768  dat.getBase() = e->getApproxBase();
769  dat.getSpace() = e->getSpace();
770  const int ent_order = e->getMaxOrder();
771  dat.getDataOrder() =
772  dat.getDataOrder() > ent_order ? dat.getDataOrder() : ent_order;
773 
774  auto ent_data = e->getEntFieldData();
775  ent_field_data.resize(ent_data.size(), false);
776  noalias(ent_field_data) = ent_data;
777  ent_field_dofs.resize(ent_data.size(), false);
778  std::fill(ent_field_dofs.begin(), ent_field_dofs.end(), nullptr);
779  ent_field.resize(1, false);
780  ent_field[0] = e.get();
781  if (auto cache = e->entityCacheDataDofs.lock()) {
782  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
783  ent_field_dofs[(*dit)->getEntDofIdx()] =
784  reinterpret_cast<FEDofEntity *>((*dit).get());
785  }
786  }
787  }
788  }
789 
790  for (auto &it : brother_ents_vec) {
791  if (const auto e = it.lock()) {
792  const EntityType type = e->getEntType();
793  const int side = e->getSideNumberPtr()->side_number;
794  const int brother_side = e->getSideNumberPtr()->brother_side_number;
795  auto &dat = data.dataOnEntities[type][side];
796  auto &dat_brother = data.dataOnEntities[type][brother_side];
797  dat_brother.getBase() = dat.getBase();
798  dat_brother.getSpace() = dat.getSpace();
799  dat_brother.getDataOrder() = dat.getDataOrder();
800  dat_brother.getFieldData() = dat.getFieldData();
801  dat_brother.getFieldDofs() = dat.getFieldDofs();
802  dat_brother.getFieldEntities() = dat.getFieldEntities();
803  }
804  }
805  }
806  }
807 
809 }
virtual FieldBitNumber get_field_bit_number(const std::string name) const =0
get field bit number

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

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

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

422  {
423 
424  struct Extractor {
425  boost::weak_ptr<EntityCacheNumeredDofs>
426  operator()(boost::shared_ptr<FieldEntity> &e) {
427  return e->entityCacheRowDofs;
428  }
429  };
430 
431  return getEntityIndices(data, field_name, getRowFieldEnts(), type_lo, type_hi,
432  Extractor());
433 }
auto & getRowFieldEnts() const

◆ 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 }

◆ 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 536 of file ForcesAndSourcesCore.hpp.

536  {
537  return getEntitySense(type, data.dataOnEntities[type]);
538  }
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 878 of file ForcesAndSourcesCore.cpp.

878  {
880  // PetscAttachDebugger();
881  data.facesNodes.resize(4, 3, false);
882  auto &side_table = const_cast<SideNumber_multiIndex &>(
883  numeredEntFiniteElementPtr->getSideNumberTable());
884  auto siit = side_table.get<1>().lower_bound(boost::make_tuple(MBTRI, 0));
885  auto hi_siit = side_table.get<1>().upper_bound(boost::make_tuple(MBTRI, 4));
886  if (std::distance(siit, hi_siit) != 4) {
887  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
888  "Should be 4 triangles on tet, side_table not initialized");
889  }
890  const int canonical_face_sense_p1[4][3] = {
891  {0, 1, 3},
892  {1, 2, 3},
893  {0, 3, 2} /**/,
894  {0, 2, 1} /**/}; // second index is offset (positive sense)
895  const int canonical_face_sense_m1[4][3] = {
896  {0, 3, 1},
897  {1, 3, 2},
898  {0, 2, 3},
899  {0, 1, 2}}; // second index is offset (negative sense
900  for (; siit != hi_siit; siit++) {
901  const boost::shared_ptr<SideNumber> side = *siit;
902  int face_conn[3] = {-1, -1, -1};
903  if (side->offset == 0) {
904  face_conn[0] = side->sense == 1
905  ? canonical_face_sense_p1[(int)side->side_number][0]
906  : canonical_face_sense_m1[(int)side->side_number][0];
907  face_conn[1] = side->sense == 1
908  ? canonical_face_sense_p1[(int)side->side_number][1]
909  : canonical_face_sense_m1[(int)side->side_number][1];
910  face_conn[2] = side->sense == 1
911  ? canonical_face_sense_p1[(int)side->side_number][2]
912  : canonical_face_sense_m1[(int)side->side_number][2];
913  }
914  if (side->offset == 1) {
915  face_conn[0] =
916  side->sense == 1
917  ? canonical_face_sense_p1[(int)side->side_number][1]
918  : canonical_face_sense_m1[(int)side->side_number][2] /**/;
919  face_conn[1] = side->sense == 1
920  ? canonical_face_sense_p1[(int)side->side_number][2]
921  : canonical_face_sense_m1[(int)side->side_number][0];
922  face_conn[2] = side->sense == 1
923  ? canonical_face_sense_p1[(int)side->side_number][0]
924  : canonical_face_sense_m1[(int)side->side_number][1];
925  }
926  if (side->offset == 2) {
927  face_conn[0] =
928  side->sense == 1
929  ? canonical_face_sense_p1[(int)side->side_number][2]
930  : canonical_face_sense_m1[(int)side->side_number][1] /**/;
931  face_conn[1] = side->sense == 1
932  ? canonical_face_sense_p1[(int)side->side_number][0]
933  : canonical_face_sense_m1[(int)side->side_number][2];
934  face_conn[2] = side->sense == 1
935  ? canonical_face_sense_p1[(int)side->side_number][1]
936  : canonical_face_sense_m1[(int)side->side_number][0];
937  }
938  for (int nn = 0; nn < 3; nn++)
939  data.facesNodes(side->side_number, nn) = face_conn[nn];
940  {
941  const EntityHandle *conn_tet;
942  int num_nodes_tet;
943  EntityHandle ent = numeredEntFiniteElementPtr->getEnt();
944  CHKERR mField.get_moab().get_connectivity(ent, conn_tet, num_nodes_tet,
945  true);
946  if (num_nodes_tet != 4)
947  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
948  "data inconsistency");
949  int num_nodes_face;
950  const EntityHandle *conn_face;
951  CHKERR mField.get_moab().get_connectivity(side->ent, conn_face,
952  num_nodes_face, true);
953  if (num_nodes_face != 3)
954  SETERRQ(PETSC_COMM_SELF, 1, "data inconsistency");
955  if (conn_face[0] != conn_tet[data.facesNodes(side->side_number, 0)])
956  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
957  "data inconsistency");
958  if (conn_face[1] != conn_tet[data.facesNodes(side->side_number, 1)])
959  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
960  "data inconsistency");
961  if (conn_face[2] != conn_tet[data.facesNodes(side->side_number, 2)])
962  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
963  "data inconsistency");
964  }
965  }
967 }
multi_index_container< boost::shared_ptr< SideNumber >, indexed_by< ordered_unique< member< SideNumber, EntityHandle, &SideNumber::ent > >, ordered_non_unique< composite_key< SideNumber, const_mem_fun< SideNumber, EntityType, &SideNumber::getEntType >, member< SideNumber, char, &SideNumber::side_number > > > > > SideNumber_multiIndex
SideNumber_multiIndex for SideNumber.
virtual moab::Interface & get_moab()=0

◆ getMaxColOrder()

int MoFEM::ForcesAndSourcesCore::getMaxColOrder ( ) const

Get max order of approximation for field in columns.

Definition at line 142 of file ForcesAndSourcesCore.cpp.

142  {
143  return getMaxOrder(getColFieldEnts());
144 }
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 }

◆ 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 }

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

603  {
604 
605  auto get_nodes_field_data = [&](VectorDouble &nodes_data,
606  VectorFieldEntities &field_entities,
607  VectorDofs &nodes_dofs, FieldSpace &space,
609  VectorInt &bb_node_order) {
611 
612  nodes_data.resize(0, false);
613  nodes_dofs.resize(0, false);
614  field_entities.resize(0, false);
615 
616  auto field_it = fieldsPtr->get<FieldName_mi_tag>().find(field_name);
617  if (field_it != fieldsPtr->get<FieldName_mi_tag>().end()) {
618 
619  auto bit_number = (*field_it)->getBitNumber();
620  const int nb_dofs_on_vert = (*field_it)->getNbOfCoeffs();
621  space = (*field_it)->getSpace();
622  base = (*field_it)->getApproxBase();
623 
624  auto &field_ents = getDataFieldEnts();
625  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
626  bit_number, get_id_for_min_type<MBVERTEX>());
627  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
628  cmp_uid_lo);
629  if (lo != field_ents.end()) {
630  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
631  bit_number, get_id_for_max_type<MBVERTEX>());
632  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
633  if (lo != hi) {
634 
635  int nb_dofs = 0;
636  for (auto it = lo; it != hi; ++it) {
637  if (auto e = it->lock()) {
638  nb_dofs += e->getNbDofsOnEnt();
639  }
640  }
641 
642  if (nb_dofs) {
643 
644  int num_nodes;
645  CHKERR getNumberOfNodes(num_nodes);
646  bb_node_order.resize(num_nodes, false);
647  bb_node_order.clear();
648  const int max_nb_dofs = nb_dofs_on_vert * num_nodes;
649  nodes_data.resize(max_nb_dofs, false);
650  nodes_dofs.resize(max_nb_dofs, false);
651  field_entities.resize(num_nodes, false);
652  std::fill(nodes_data.begin(), nodes_data.end(), 0);
653  std::fill(nodes_dofs.begin(), nodes_dofs.end(), nullptr);
654  std::fill(field_entities.begin(), field_entities.end(), nullptr);
655 
656  std::vector<boost::weak_ptr<FieldEntity>> brother_ents_vec;
657 
658  for (auto it = lo; it != hi; ++it) {
659  if (auto e = it->lock()) {
660  const auto &sn = e->getSideNumberPtr();
661  const int side_number = sn->side_number;
662  // Some field entities on skeleton can have negative side
663  // numbeer
664  if (side_number >= 0) {
665  const int brother_side_number = sn->brother_side_number;
666 
667  field_entities[side_number] = e.get();
668  if (brother_side_number != -1) {
669  brother_ents_vec.emplace_back(e);
670  field_entities[side_number] = field_entities[side_number];
671  }
672 
673  bb_node_order[side_number] = e->getMaxOrder();
674  int pos = side_number * nb_dofs_on_vert;
675  auto ent_filed_data_vec = e->getEntFieldData();
676  if (auto cache = e->entityCacheDataDofs.lock()) {
677  for (auto dit = cache->loHi[0]; dit != cache->loHi[1];
678  ++dit) {
679  const auto dof_idx = (*dit)->getEntDofIdx();
680  nodes_data[pos + dof_idx] = ent_filed_data_vec[dof_idx];
681  nodes_dofs[pos + dof_idx] =
682  reinterpret_cast<FEDofEntity *>((*dit).get());
683  }
684  }
685  }
686  }
687  }
688 
689  for (auto &it : brother_ents_vec) {
690  if (const auto e = it.lock()) {
691  const auto &sn = e->getSideNumberPtr();
692  const int side_number = sn->side_number;
693  const int brother_side_number = sn->brother_side_number;
694  bb_node_order[brother_side_number] = bb_node_order[side_number];
695  int pos = side_number * nb_dofs_on_vert;
696  int brother_pos = brother_side_number * nb_dofs_on_vert;
697  for (int ii = 0; ii != nb_dofs_on_vert; ++ii) {
698  nodes_data[brother_pos] = nodes_data[pos];
699  nodes_dofs[brother_pos] = nodes_dofs[pos];
700  ++pos;
701  ++brother_pos;
702  }
703  }
704  }
705  }
706  }
707  }
708  }
709 
711  };
712 
713  return get_nodes_field_data(
714  data.dataOnEntities[MBVERTEX][0].getFieldData(),
715  data.dataOnEntities[MBVERTEX][0].getFieldEntities(),
716  data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
717  data.dataOnEntities[MBVERTEX][0].getSpace(),
718  data.dataOnEntities[MBVERTEX][0].getBase(),
719  data.dataOnEntities[MBVERTEX][0].getBBNodeOrder());
720 }
ublas::vector< int, IntAllocator > VectorInt
Definition: Types.hpp:73
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74
ublas::vector< FieldEntity *, FieldEntAllocator > VectorFieldEntities
ublas::vector< FEDofEntity *, DofsAllocator > VectorDofs

◆ 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  if (side_number >= 0) {
282  const auto brother_side_number = side_ptr->brother_side_number;
283  if (auto cache = extractor(e).lock()) {
284  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
285  auto &dof = **dit;
286  const int idx = dof.getPetscGlobalDofIdx();
287  const int local_idx = dof.getPetscLocalDofIdx();
288  const int pos =
289  side_number * nb_dofs_on_vert + dof.getDofCoeffIdx();
290  nodes_indices[pos] = idx;
291  local_nodes_indices[pos] = local_idx;
292  if (brother_side_number != -1) {
293  const int elem_idx = brother_side_number * nb_dofs_on_vert +
294  (*dit)->getDofCoeffIdx();
295  nodes_indices[elem_idx] = idx;
296  local_nodes_indices[elem_idx] = local_idx;
297  }
298  }
299  }
300  }
301  }
302  }
303  } else {
304  nodes_indices.resize(0, false);
305  local_nodes_indices.resize(0, false);
306  }
307 
308  } else {
309  nodes_indices.resize(0, false);
310  local_nodes_indices.resize(0, false);
311  }
312 
314 }

◆ getNoFieldColIndices()

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

get col NoField indices

Definition at line 479 of file ForcesAndSourcesCore.cpp.

480  {
482  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
483  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
484  }
485  CHKERR getNoFieldIndices(field_name, getColDofsPtr(),
486  data.dataOnEntities[MBENTITYSET][0].getIndices());
488 }
auto getColDofsPtr() const
MoFEMErrorCode getNoFieldIndices(const std::string &field_name, boost::shared_ptr< FENumeredDofEntity_multiIndex > dofs, VectorInt &nodes_indices) const
get NoField indices

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

813  {
814 
816 
817  ent_field_data.resize(0, false);
818  ent_field_dofs.resize(0, false);
819  ent_field.resize(0, false);
820 
821  auto &field_ents = getDataFieldEnts();
822  auto bit_number = mField.get_field_bit_number(field_name);
823  const auto lo_uid = FieldEntity::getLocalUniqueIdCalculate(
824  bit_number, get_id_for_min_type<MBVERTEX>());
825  auto lo = std::lower_bound(field_ents.begin(), field_ents.end(), lo_uid,
826  cmp_uid_lo);
827  if (lo != field_ents.end()) {
828 
829  ent_field.resize(field_ents.size(), false);
830  std::fill(ent_field.begin(), ent_field.end(), nullptr);
831 
832  const auto hi_uid = FieldEntity::getLocalUniqueIdCalculate(
833  bit_number, get_id_for_max_type<MBVERTEX>());
834  auto hi = std::upper_bound(lo, field_ents.end(), hi_uid, cmp_uid_hi);
835  if (lo != hi) {
836 
837  int side = 0;
838  for (auto it = lo; it != hi; ++it, ++side)
839  if (auto e = it->lock()) {
840 
841  const auto size = e->getNbDofsOnEnt();
842  ent_field_data.resize(size, false);
843  ent_field_dofs.resize(size, false);
844  ent_field[side] = e.get();
845  noalias(ent_field_data) = e->getEntFieldData();
846 
847  if (auto cache = e->entityCacheDataDofs.lock()) {
848  for (auto dit = cache->loHi[0]; dit != cache->loHi[1]; ++dit) {
849  ent_field_dofs[(*dit)->getEntDofIdx()] =
850  reinterpret_cast<FEDofEntity *>((*dit).get());
851  }
852  }
853  }
854  }
855  }
856 
858 }

◆ getNoFieldFieldData() [2/2]

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

Definition at line 861 of file ForcesAndSourcesCore.cpp.

862  {
864  if (data.dataOnEntities[MBENTITYSET].size() == 0)
865  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
866  "No space to insert data");
867 
869  field_name, data.dataOnEntities[MBENTITYSET][0].getFieldData(),
870  data.dataOnEntities[MBENTITYSET][0].getFieldDofs(),
871  data.dataOnEntities[MBENTITYSET][0].getFieldEntities());
873 }
MoFEMErrorCode getNoFieldFieldData(const std::string field_name, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs, VectorFieldEntities &ent_field) const
Get field data on nodes.

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

453  {
455  auto field_it = fieldsPtr->get<FieldName_mi_tag>().find(field_name);
456  auto dit = dofs->get<Unique_mi_tag>().lower_bound(
457  FieldEntity::getLoBitNumberUId((*field_it)->getBitNumber()));
458  auto hi_dit = dofs->get<Unique_mi_tag>().upper_bound(
459  FieldEntity::getHiBitNumberUId((*field_it)->getBitNumber()));
460  indices.resize(std::distance(dit, hi_dit));
461  for (; dit != hi_dit; dit++) {
462  int idx = (*dit)->getPetscGlobalDofIdx();
463  indices[(*dit)->getDofCoeffIdx()] = idx;
464  }
466 }
static UId getHiBitNumberUId(const FieldBitNumber bit_number)
static UId getLoBitNumberUId(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 468 of file ForcesAndSourcesCore.cpp.

469  {
471  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
472  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
473  }
474  CHKERR getNoFieldIndices(field_name, getRowDofsPtr(),
475  data.dataOnEntities[MBENTITYSET][0].getIndices());
477 }
auto getRowDofsPtr() const

◆ 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
EshelbianPlasticity.cpp, MagneticElement.hpp, continuity_check_on_skeleton_with_simple_2d_for_h1.cpp, continuity_check_on_skeleton_with_simple_2d_for_hcurl.cpp, continuity_check_on_skeleton_with_simple_2d_for_hdiv.cpp, elasticity_mixed_formulation.cpp, lorentz_force.cpp, nonlinear_dynamics.cpp, and simple_contact.cpp.

Definition at line 425 of file ForcesAndSourcesCore.hpp.

425 { 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 585 of file ForcesAndSourcesCore.cpp.

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

◆ getProblemNodesIndices()

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

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

Definition at line 492 of file ForcesAndSourcesCore.cpp.

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

◆ getProblemNodesRowIndices()

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

Definition at line 571 of file ForcesAndSourcesCore.cpp.

572  {
573  return getProblemNodesIndices(field_name, *(problemPtr->numeredRowDofsPtr),
574  nodes_indices);
575 }
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredRowDofsPtr
store DOFs on rows for this problem

◆ getProblemTypeColIndices()

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

Definition at line 592 of file ForcesAndSourcesCore.cpp.

594  {
595  return getProblemTypeIndices(field_name, *(problemPtr->numeredColDofsPtr),
596  type, side_number, indices);
597 }
MoFEMErrorCode getProblemTypeIndices(const std::string &field_name, const NumeredDofEntity_multiIndex &dofs, EntityType type, int side_number, VectorInt &indices) const
get indices by type (generic function) which are declared for problem but not this particular element

◆ getProblemTypeIndices()

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

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

Definition at line 536 of file ForcesAndSourcesCore.cpp.

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

◆ getProblemTypeRowIndices()

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

Definition at line 578 of file ForcesAndSourcesCore.cpp.

580  {
581  return getProblemTypeIndices(field_name, *(problemPtr->numeredRowDofsPtr),
582  type, side_number, indices);
583 }

◆ getRowNodesIndices()

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

get row node indices from FENumeredDofEntity_multiIndex

Definition at line 317 of file ForcesAndSourcesCore.cpp.

318  {
319 
320  struct Extractor {
321  boost::weak_ptr<EntityCacheNumeredDofs>
322  operator()(boost::shared_ptr<FieldEntity> &e) {
323  return e->entityCacheRowDofs;
324  }
325  };
326 
327  return getNodesIndices(field_name, getRowFieldEnts(),
328  data.dataOnEntities[MBVERTEX][0].getIndices(),
329  data.dataOnEntities[MBVERTEX][0].getLocalIndices(),
330  Extractor());
331 }

◆ getRule() [1/2]

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

◆ getRule() [2/2]

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

another variant of getRule

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

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

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

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

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

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

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

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

Reimplemented in EdgeFE.

Definition at line 1591 of file ForcesAndSourcesCore.cpp.

1592  {
1593  return getRuleHook ? getRuleHook(order_row, order_col, order_data)
1594  : getRule(order_data);
1595 }

◆ getSpacesAndBaseOnEntities()

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

Get field approximation space and base on entities.

Definition at line 971 of file ForcesAndSourcesCore.cpp.

972  {
974 
975  if (nInTheLoop == 0) {
976  data.sPace.reset();
977  data.bAse.reset();
978  for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
979  data.spacesOnEntities[t].reset();
980  data.basesOnEntities[t].reset();
981  }
982  for (int s = 0; s != LASTSPACE; ++s) {
983  data.basesOnSpaces[s].reset();
984  }
985  }
986 
987  auto fe_type = numeredEntFiniteElementPtr->getEntType();
988 
989  if (getDataFieldEntsPtr())
990  for (auto e : getDataFieldEnts()) {
991  if (auto ptr = e.lock()) {
992  // get data from entity
993  const EntityType type = ptr->getEntType();
995  const FieldSpace space = ptr->getSpace();
996  const FieldApproximationBase approx = ptr->getApproxBase();
997  // set data
998  data.sPace.set(space);
999  data.bAse.set(approx);
1000  data.spacesOnEntities[type].set(space);
1001  data.basesOnEntities[type].set(approx);
1002  data.basesOnSpaces[space].set(approx);
1003  }
1004  }
1005  }
1006  else
1007  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1008  "data fields ents not allocated on element");
1009 
1011 }
int nInTheLoop
number currently of processed method
static bool getDefTypeMap(const EntityType fe_type, const EntityType ent_type)
boost::shared_ptr< FieldEntity_vector_view > & getDataFieldEntsPtr() const

◆ getUserPolynomialBase()

auto& MoFEM::ForcesAndSourcesCore::getUserPolynomialBase ( )

Get the User Polynomial Base object.

Returns
boost::shared_ptr<BaseFunction>&
Examples
EshelbianPlasticity.cpp.

Definition at line 439 of file ForcesAndSourcesCore.hpp.

439 { 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 1304 of file ForcesAndSourcesCore.cpp.

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

◆ operator()()

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

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

Definition at line 1645 of file ForcesAndSourcesCore.cpp.

1645  {
1647  if (postProcessHook) {
1648  ierr = postProcessHook();
1649  CHKERRG(ierr);
1650  }
1652 }
boost::function< MoFEMErrorCode()> postProcessHook
Hook function for post-processing.

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

Definition at line 1629 of file ForcesAndSourcesCore.cpp.

1629  {
1631  if (preProcessHook) {
1632  ierr = preProcessHook();
1633  CHKERRG(ierr);
1634  }
1636 }
boost::function< MoFEMErrorCode()> preProcessHook
Hook function for pre-processing.

◆ setGaussPts() [1/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 NitscheMethod::MyVolumeFE, PostProcEdgeOnRefinedMesh, PostProcTemplateVolumeOnRefinedMesh< MoFEM::VolumeElementForcesAndSourcesCore >, MoFEM::VolumeElementForcesAndSourcesCoreOnSideBase, MoFEM::VolumeElementForcesAndSourcesCoreOnContactPrismSideBase, and MoFEM::FaceElementForcesAndSourcesCoreOnSideBase.

Definition at line 1607 of file ForcesAndSourcesCore.cpp.

1607  {
1609  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Sorry, not implemented");
1611 }

◆ setGaussPts() [2/2]

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

set user specific integration rule

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

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

User sets

ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76

where

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

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

Reimplemented in EdgeFE.

Definition at line 1597 of file ForcesAndSourcesCore.cpp.

1598  {
1599  return setRuleHook ? setRuleHook(order_row, order_col, order_data)
1600  : setGaussPts(order_data);
1601 }
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule

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

1542  {
1544  sidePtrFE = const_cast<ForcesAndSourcesCore *>(side_fe_ptr);
1546 }
ForcesAndSourcesCore(Interface &m_field)

Friends And Related Function Documentation

◆ FaceElementForcesAndSourcesCoreOnSideBase

Definition at line 857 of file ForcesAndSourcesCore.hpp.

◆ UserDataOperator

friend class UserDataOperator
friend

◆ VolumeElementForcesAndSourcesCoreOnContactPrismSideBase

Definition at line 858 of file ForcesAndSourcesCore.hpp.

◆ VolumeElementForcesAndSourcesCoreOnSideBase

Definition at line 856 of file ForcesAndSourcesCore.hpp.

Member Data Documentation

◆ dataH1

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataH1
protected

Definition at line 807 of file ForcesAndSourcesCore.hpp.

◆ dataHcurl

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataHcurl
protected

Definition at line 808 of file ForcesAndSourcesCore.hpp.

◆ dataHdiv

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataHdiv
protected

Definition at line 809 of file ForcesAndSourcesCore.hpp.

◆ dataL2

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataL2
protected

Definition at line 810 of file ForcesAndSourcesCore.hpp.

◆ dataNoField

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataNoField
protected

Definition at line 806 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 797 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 804 of file ForcesAndSourcesCore.hpp.

◆ elementPolynomialBasePtr

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

Pointer to entity polynomial base.

Definition at line 832 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 451 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 825 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 816 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 843 of file ForcesAndSourcesCore.hpp.

◆ userPolynomialBasePtr

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

Pointer to user polynomail base.

Definition at line 837 of file ForcesAndSourcesCore.hpp.


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