v0.9.1
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 ()
 
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...
 
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
 
boost::shared_ptr< const FENumeredDofEntity_multiIndexrowPtr
 Pointer to finite element rows dofs view. More...
 
boost::shared_ptr< const FENumeredDofEntity_multiIndexcolPtr
 Pointer to finite element columns dofs view. More...
 
boost::shared_ptr< const FEDofEntity_multiIndexdataPtr
 Pointer to finite element data dofs. More...
 
boost::shared_ptr< const FieldEntity_vector_viewrowFieldEntsPtr
 Pointer to finite element field entities row view. More...
 
boost::shared_ptr< const FieldEntity_vector_viewcolFieldEntsPtr
 Pointer to finite element field entities column view. More...
 
boost::shared_ptr< const FieldEntity_multiIndex_spaceType_viewdataFieldEntsPtr
 Pointer to finite element field entities data view. More...
 
- 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
MoFEMErrorCode getNodesIndices (const boost::string_ref field_name, FENumeredDofEntity_multiIndex &dofs, VectorInt &nodes_indices, VectorInt &local_nodes_indices) 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...
 
MoFEMErrorCode getEntityIndices (DataForcesAndSourcesCore &data, const std::string &field_name, FENumeredDofEntity_multiIndex &dofs, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) 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, 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 boost::string_ref field_name, FEDofEntity_multiIndex &dofs, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs) const
 Get field data on nodes. More...
 
MoFEMErrorCode getNoFieldFieldData (DataForcesAndSourcesCore &data, const boost::string_ref 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
 
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 33 of file ForcesAndSourcesCore.cpp.

34  :
35 
36  mField(m_field), getRuleHook(0), setRuleHook(0),
38 
39  nullptr,
40  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // NOFIELD
41  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // H1
42  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // HCURL
43  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // HDIV
44  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET) // L2
45 
46  },
48 
49  nullptr,
50  boost::make_shared<DerivedDataForcesAndSourcesCore>(
51  dataOnElement[NOFIELD]), // NOFIELD
52  boost::make_shared<DerivedDataForcesAndSourcesCore>(
53  dataOnElement[H1]), // H1
54  boost::make_shared<DerivedDataForcesAndSourcesCore>(
55  dataOnElement[HCURL]), // HCURL
56  boost::make_shared<DerivedDataForcesAndSourcesCore>(
57  dataOnElement[HDIV]), // HDIV
58  boost::make_shared<DerivedDataForcesAndSourcesCore>(
59  dataOnElement[L2]) // L2
60 
61  },
64  dataHdiv(*dataOnElement[HDIV].get()), dataL2(*dataOnElement[L2].get()),
65  lastEvaluatedElementEntityType(MBMAXTYPE), sidePtrFE(nullptr) {}
DataForcesAndSourcesCore & dataL2
field with continuous normal traction
Definition: definitions.h:178
ForcesAndSourcesCore * sidePtrFE
Element to integrate on the sides.
DataForcesAndSourcesCore & dataH1
scalar or vector of scalars describe (no true field)
Definition: definitions.h:175
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:177
DataForcesAndSourcesCore & dataHdiv
continuous field
Definition: definitions.h:176
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:179
DataForcesAndSourcesCore & dataNoField

Member Function Documentation

◆ calBernsteinBezierBaseFunctionsOnElement()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::calBernsteinBezierBaseFunctionsOnElement ( )
protected

Calculate Bernstein-Bezier base.

Returns
MoFEMErrorCode

Definition at line 861 of file ForcesAndSourcesCore.cpp.

861  {
863 
864  auto get_nodal_base_data = [&](DataForcesAndSourcesCore &data,
865  const std::string &field_name) {
867  auto &space = data.dataOnEntities[MBVERTEX][0].getSpace();
868  auto &base = data.dataOnEntities[MBVERTEX][0].getBase();
869  auto &bb_node_order = data.dataOnEntities[MBVERTEX][0].getBBNodeOrder();
870 
871  auto &dofs_by_name_and_type =
872  dataPtr->get<Composite_Name_And_Type_mi_tag>();
873  auto tuple = boost::make_tuple(field_name, MBVERTEX);
874  auto dit = dofs_by_name_and_type.lower_bound(tuple);
875  if (dit == dofs_by_name_and_type.end())
876  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
877  "No nodal dofs on element");
878  auto hi_dit =
879  dataPtr->get<Composite_Name_And_Type_mi_tag>().upper_bound(tuple);
880 
881  if (dit != hi_dit) {
882  auto &first_dof = **dit;
883  space = first_dof.getSpace();
884  base = first_dof.getApproxBase();
885  int num_nodes;
886  CHKERR getNumberOfNodes(num_nodes);
887  bb_node_order.resize(num_nodes, false);
888  bb_node_order.clear();
889  const int nb_dof_idx = first_dof.getNbOfCoeffs();
890 
891  std::vector<boost::weak_ptr<FEDofEntity>> brother_dofs_vec;
892  for (; dit != hi_dit;) {
893  const auto &dof_ptr = *dit;
894  const auto &dof = *dof_ptr;
895  const auto &sn = *dof.sideNumberPtr;
896  const int side_number = sn.side_number;
897  const int brother_side_number = sn.brother_side_number;
898  if (brother_side_number != -1)
899  brother_dofs_vec.emplace_back(dof_ptr);
900  bb_node_order[side_number] = dof.getMaxOrder();
901  for (int ii = 0; ii != nb_dof_idx; ++ii)
902  ++dit;
903  }
904 
905  for (auto &dof_ptr : brother_dofs_vec) {
906  if (const auto d = dof_ptr.lock()) {
907  const auto &sn = d->sideNumberPtr;
908  const int side_number = sn->side_number;
909  const int brother_side_number = sn->brother_side_number;
910  bb_node_order[brother_side_number] = bb_node_order[side_number];
911  }
912  }
913  }
915  };
916 
917  auto get_entity_base_data = [&](DataForcesAndSourcesCore &data,
918  const std::string &field_name,
919  const EntityType type_lo,
920  const EntityType type_hi) {
922  for (EntityType t = type_lo; t != type_hi; ++t) {
923  for (auto &dat : data.dataOnEntities[t]) {
924  dat.getDataOrder() = 0;
925  dat.getBase() = NOBASE;
926  dat.getSpace() = NOSPACE;
927  dat.getFieldData().resize(0, false);
928  dat.getFieldDofs().resize(0, false);
929  }
930  }
931 
932  auto &dofs = const_cast<FEDofEntity_multiIndex &>(
933  numeredEntFiniteElementPtr->getDataDofs());
934  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
935  auto dit = dofs_by_type.lower_bound(boost::make_tuple(field_name, type_lo));
936  if (dit == dofs_by_type.end())
938  auto hi_dit =
939  dofs_by_type.lower_bound(boost::make_tuple(field_name, type_hi));
940  std::vector<boost::weak_ptr<FEDofEntity>> brother_dofs_vec;
941  for (; dit != hi_dit;) {
942 
943  auto &dof = **dit;
944  const int nb_dofs_on_ent = dof.getNbDofsOnEnt();
945  if (nb_dofs_on_ent) {
946  const EntityType type = dof.getEntType();
947  const int side = dof.sideNumberPtr->side_number;
948  auto &dat = data.dataOnEntities[type][side];
949 
950  const int brother_side = dof.sideNumberPtr->brother_side_number;
951  if (brother_side != -1)
952  brother_dofs_vec.emplace_back(*dit);
953 
954  dat.getBase() = dof.getApproxBase();
955  dat.getSpace() = dof.getSpace();
956  const int ent_order = dof.getMaxOrder();
957  dat.getDataOrder() =
958  dat.getDataOrder() > ent_order ? dat.getDataOrder() : ent_order;
959  for (int ii = 0; ii != nb_dofs_on_ent; ++ii) {
960  ++dit;
961  }
962  }
963  }
964 
965  for (auto &dof_ptr : brother_dofs_vec) {
966  if (auto d = dof_ptr.lock()) {
967  const EntityType type = d->getEntType();
968  const int side = d->sideNumberPtr->side_number;
969  const int brother_side = d->sideNumberPtr->brother_side_number;
970  auto &dat = data.dataOnEntities[type][side];
971  auto &dat_brother = data.dataOnEntities[type][brother_side];
972  dat_brother.getBase() = dat.getBase();
973  dat_brother.getSpace() = dat.getSpace();
974  dat_brother.getDataOrder() = dat.getDataOrder();
975  }
976  }
978  };
979 
980  auto &ents_data = *dataFieldEntsPtr;
981  for (auto &e : ents_data) {
982  if (e->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
983  auto space = e->getSpace();
984  for (EntityType t = MBVERTEX; t != MBPOLYHEDRON; ++t) {
985  for (auto &dat : (*dataOnElement[space]).dataOnEntities[t]) {
986  for (auto &ptr : dat.getBBAlphaIndicesByOrderArray())
987  ptr.reset();
988  for (auto &ptr : dat.getBBNByOrderArray())
989  ptr.reset();
990  for (auto &ptr : dat.getBBDiffNByOrderArray())
991  ptr.reset();
992  }
993  }
994  }
995  }
996 
997  for (auto &e : ents_data) {
998  if (e->getApproxBase() == AINSWORTH_BERNSTEIN_BEZIER_BASE) {
999  auto field_name = e->getName();
1000  auto space = e->getSpace();
1001  CHKERR get_nodal_base_data(*dataOnElement[space], field_name);
1002  CHKERR get_entity_base_data(*dataOnElement[space], field_name, MBEDGE,
1003  MBPOLYHEDRON);
1004  CHKERR getElementPolynomialBase()->getValue(
1005  gaussPts,
1006  boost::make_shared<EntPolynomialBaseCtx>(
1007  *dataOnElement[space], field_name, static_cast<FieldSpace>(space),
1009  }
1010  }
1011 
1013 };
MatrixDouble gaussPts
Matrix of integration points.
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
boost::shared_ptr< const FEDofEntity_multiIndex > dataPtr
Pointer to finite element data dofs.
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
#define CHKERR
Inline error check.
Definition: definitions.h:601
boost::shared_ptr< const FieldEntity_multiIndex_spaceType_view > dataFieldEntsPtr
Pointer to finite element field entities data view.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
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 791 of file ForcesAndSourcesCore.cpp.

792  {
794  if (dataOnElement[H1]->bAse.test(b)) {
795  switch (static_cast<FieldApproximationBase>(b)) {
796  case NOBASE:
797  break;
799  // BERNSTEIN_BEZIER_BASE is not hierarchical base
800  break;
805  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
806  "Functions genrating approximation base not defined");
807 
808  for (int space = H1; space != LASTSPACE; ++space) {
809  if (dataOnElement[H1]->sPace.test(space) &&
810  dataOnElement[H1]->bAse.test(b) &&
811  dataOnElement[H1]->basesOnSpaces[space].test(b)) {
812  CHKERR getElementPolynomialBase()->getValue(
813  gaussPts,
814  boost::make_shared<EntPolynomialBaseCtx>(
815  *dataOnElement[space], static_cast<FieldSpace>(space),
816  static_cast<FieldApproximationBase>(b), NOBASE));
817  }
818  }
819  break;
820  case USER_BASE:
821  if (!getUserPolynomialBase())
822  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
823  "Functions genrating user approximation base not defined");
824 
825  for (int space = H1; space != LASTSPACE; ++space)
826  if (dataOnElement[H1]->sPace.test(space) &&
827  dataOnElement[H1]->bAse.test(b) &&
828  dataOnElement[H1]->basesOnSpaces[space].test(b)) {
829  CHKERR getUserPolynomialBase()->getValue(
830  gaussPts,
831  boost::make_shared<EntPolynomialBaseCtx>(
832  *dataOnElement[space], static_cast<FieldSpace>(space),
833  static_cast<FieldApproximationBase>(b), NOBASE));
834  }
835  break;
836  default:
838  "Base <%s> not yet implemented",
839  ApproximationBaseNames[static_cast<FieldApproximationBase>(b)]);
840  }
841  }
843 }
MatrixDouble gaussPts
Matrix of integration points.
user implemented approximation base
Definition: definitions.h:159
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:482
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:180
static const char *const ApproximationBaseNames[]
Definition: definitions.h:163
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:151
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
virtual MPI_Comm & get_comm() const =0
continuous field
Definition: definitions.h:176
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 845 of file ForcesAndSourcesCore.cpp.

845  {
847  /// Use the some node base. Node base is usually used for construction other
848  /// bases.
849  for (int space = HCURL; space != LASTSPACE; ++space) {
850  dataOnElement[space]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
851  dataOnElement[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
852  }
853  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
855  static_cast<FieldApproximationBase>(b));
856  }
858 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:180
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:151
field with continuous tangents
Definition: definitions.h:177
#define CHKERR
Inline error check.
Definition: definitions.h:601
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:412
continuous field
Definition: definitions.h:176
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 1015 of file ForcesAndSourcesCore.cpp.

1015  {
1017 
1018  const EntityType type = numeredEntFiniteElementPtr->getEntType();
1021 
1022  // Data on elements for proper spaces
1023  for (int space = H1; space != LASTSPACE; ++space) {
1024  dataOnElement[space]->setElementType(type);
1025  derivedDataOnElement[space]->setElementType(type);
1026  }
1027 
1029 
1031 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:180
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:412
continuous field
Definition: definitions.h:176
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 228 of file ForcesAndSourcesCore.cpp.

229  {
230  return getNodesIndices(field_name,
231  const_cast<FENumeredDofEntity_multiIndex &>(
232  numeredEntFiniteElementPtr->getColsDofs()),
233  data.dataOnEntities[MBVERTEX][0].getIndices(),
234  data.dataOnEntities[MBVERTEX][0].getLocalIndices());
235 }
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
MoFEMErrorCode getNodesIndices(const boost::string_ref field_name, FENumeredDofEntity_multiIndex &dofs, VectorInt &nodes_indices, VectorInt &local_nodes_indices) const
get node indices

◆ getElementPolynomialBase()

auto& MoFEM::ForcesAndSourcesCore::getElementPolynomialBase ( )

Get the Entity Polynomial Base object.

Returns
boost::shared_ptr<BaseFunction>&&

Definition at line 427 of file ForcesAndSourcesCore.hpp.

427 { 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 481 of file ForcesAndSourcesCore.hpp.

483  {
484  return dataOnElement[space]->dataOnEntities[type][side];
485  }
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 496 of file ForcesAndSourcesCore.hpp.

496  {
497  return dataOnElement[space]->dataOnEntities[type][side];
498  }
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 862 of file ForcesAndSourcesCore.hpp.

864  {
865  return getEntityIndices(data, field_name,
866  const_cast<FENumeredDofEntity_multiIndex &>(
867  numeredEntFiniteElementPtr->getColsDofs()),
868  type_lo, type_hi);
869 }
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
MoFEMErrorCode getEntityIndices(DataForcesAndSourcesCore &data, const std::string &field_name, FENumeredDofEntity_multiIndex &dofs, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) 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 119 of file ForcesAndSourcesCore.cpp.

121  {
123 
124  auto &side_table = numeredEntFiniteElementPtr->getSideNumberTable();
125 
126  for (unsigned int s = 0; s != data.size(); ++s)
127  data[s].getDataOrder() = 0;
128 
129  auto &fields_ents =
130  dataFieldEntsPtr->get<Composite_EntType_and_Space_mi_tag>();
131 
132  for (auto r = fields_ents.equal_range(boost::make_tuple(type, space));
133  r.first != r.second; ++r.first) {
134 
135  auto &e = **r.first;
136 
137  auto sit = side_table.find(e.getEnt());
138  if (sit != side_table.end()) {
139 
140  auto &side = *sit;
141  const int side_number = side->side_number;
142  if (side_number >= 0) {
143  ApproximationOrder ent_order = e.getMaxOrder();
144  auto &dat = data[side_number];
145  dat.getDataOrder() =
146  dat.getDataOrder() > ent_order ? dat.getDataOrder() : ent_order;
147  }
148  } else
149  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
150  "Entity on side of the element not found");
151  }
152 
153  for (auto r = side_table.get<2>().equal_range(type); r.first != r.second;
154  ++r.first) {
155  const int brother_side_number = (*r.first)->brother_side_number;
156  if (brother_side_number != -1) {
157  const int side_number = (*r.first)->side_number;
158  data[brother_side_number].getDataOrder() =
159  data[side_number].getDataOrder();
160  }
161  }
162 
164 }
int ApproximationOrder
Approximation on the entity.
Definition: Types.hpp:37
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
boost::shared_ptr< const FieldEntity_multiIndex_spaceType_view > dataFieldEntsPtr
Pointer to finite element field entities data view.
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

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

545  {
546  return getEntityDataOrder(type, space, data.dataOnEntities[type]);
547  }
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 542 of file ForcesAndSourcesCore.cpp.

544  {
546  for (EntityType t = type_lo; t != type_hi; ++t) {
547  for (auto &dat : data.dataOnEntities[t]) {
548  dat.getDataOrder() = 0;
549  dat.getBase() = NOBASE;
550  dat.getSpace() = NOSPACE;
551  dat.getFieldData().resize(0, false);
552  dat.getFieldDofs().resize(0, false);
553  }
554  }
555 
556  auto &dofs = const_cast<FEDofEntity_multiIndex &>(
557  numeredEntFiniteElementPtr->getDataDofs());
558  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
559  auto dit = dofs_by_type.lower_bound(boost::make_tuple(field_name, type_lo));
560  if (dit == dofs_by_type.end())
562  auto hi_dit =
563  dofs_by_type.lower_bound(boost::make_tuple(field_name, type_hi));
564  std::vector<boost::weak_ptr<FEDofEntity>> brother_dofs_vec;
565  for (; dit != hi_dit;) {
566 
567  auto &dof = **dit;
568  const int nb_dofs_on_ent = dof.getNbDofsOnEnt();
569  if (nb_dofs_on_ent) {
570 
571  const EntityType type = dof.getEntType();
572  const int side = dof.sideNumberPtr->side_number;
573  if (side >= 0) {
574 
575  auto &dat = data.dataOnEntities[type][side];
576  auto &ent_field_dofs = dat.getFieldDofs();
577  auto &ent_field_data = dat.getFieldData();
578  const int brother_side = dof.sideNumberPtr->brother_side_number;
579  if (brother_side != -1)
580  brother_dofs_vec.emplace_back(*dit);
581 
582  if (ent_field_data.empty()) {
583 
584  dat.getBase() = dof.getApproxBase();
585  dat.getSpace() = dof.getSpace();
586  const int ent_order = dof.getMaxOrder();
587  dat.getDataOrder() =
588  dat.getDataOrder() > ent_order ? dat.getDataOrder() : ent_order;
589  ent_field_data.resize(nb_dofs_on_ent, false);
590  noalias(ent_field_data) = dof.getEntFieldData();
591  ent_field_dofs.resize(nb_dofs_on_ent, false);
592  for (int ii = 0; ii != nb_dofs_on_ent; ++ii) {
593  ent_field_dofs[ii] = *dit;
594  ++dit;
595  }
596  }
597 
598  } else {
599 
600  for (int ii = 0; ii != nb_dofs_on_ent; ++ii)
601  ++dit;
602 
603  }
604 
605  }
606  }
607 
608  for (auto &dof_ptr : brother_dofs_vec) {
609  if (auto d = dof_ptr.lock()) {
610  const EntityType type = d->getEntType();
611  const int side = d->sideNumberPtr->side_number;
612  const int brother_side = d->sideNumberPtr->brother_side_number;
613  auto &dat = data.dataOnEntities[type][side];
614  auto &dat_brother = data.dataOnEntities[type][brother_side];
615  dat_brother.getBase() = dat.getBase();
616  dat_brother.getSpace() = dat.getSpace();
617  dat_brother.getDataOrder() = dat.getDataOrder();
618  dat_brother.getFieldData() = dat.getFieldData();
619  dat_brother.getFieldDofs() = dat.getFieldDofs();
620  }
621  }
622 
624 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ getEntityIndices()

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

Definition at line 237 of file ForcesAndSourcesCore.cpp.

240  {
242 
243  for (EntityType t = type_lo; t != type_hi; ++t) {
244  for (auto &dat : data.dataOnEntities[t]) {
245  dat.getIndices().resize(0, false);
246  dat.getLocalIndices().resize(0, false);
247  }
248  }
249 
250  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
251  auto dit = dofs_by_type.lower_bound(boost::make_tuple(field_name, type_lo));
252  if (dit == dofs_by_type.end())
254  auto hi_dit =
255  dofs_by_type.lower_bound(boost::make_tuple(field_name, type_hi));
256  for (; dit != hi_dit; ++dit) {
257  auto &dof = **dit;
258  const EntityType type = dof.getEntType();
259  const int side = dof.sideNumberPtr->side_number;
260 
261  if (side >= 0) {
262 
263  auto &dat = data.dataOnEntities[type][side];
264  const int nb_dofs_on_ent = dof.getNbDofsOnEnt();
265  if (nb_dofs_on_ent) {
266  const int brother_side = dof.sideNumberPtr->brother_side_number;
267  auto &ent_field_indices = dat.getIndices();
268  auto &ent_field_local_indices = dat.getLocalIndices();
269  if (ent_field_indices.empty()) {
270  ent_field_indices.resize(nb_dofs_on_ent, false);
271  ent_field_local_indices.resize(nb_dofs_on_ent, false);
272  std::fill(ent_field_indices.data().begin(),
273  ent_field_indices.data().end(), -1);
274  std::fill(ent_field_local_indices.data().begin(),
275  ent_field_local_indices.data().end(), -1);
276  }
277  const int idx = dof.getEntDofIdx();
278  ent_field_indices[idx] = dof.getPetscGlobalDofIdx();
279  ent_field_local_indices[idx] = dof.getPetscLocalDofIdx();
280  if (brother_side != -1) {
281  auto &dat_brother = data.dataOnEntities[type][brother_side];
282  dat_brother.getIndices().resize(nb_dofs_on_ent, false);
283  dat_brother.getLocalIndices().resize(nb_dofs_on_ent, false);
284  dat_brother.getIndices()[idx] = dat.getIndices()[idx];
285  dat_brother.getLocalIndices()[idx] = dat.getLocalIndices()[idx];
286  }
287  }
288 
289  }
290  }
291 
293 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

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

855  {
856  return getEntityIndices(data, field_name,
857  const_cast<FENumeredDofEntity_multiIndex &>(
858  numeredEntFiniteElementPtr->getRowsDofs()),
859  type_lo, type_hi);
860 }
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
MoFEMErrorCode getEntityIndices(DataForcesAndSourcesCore &data, const std::string &field_name, FENumeredDofEntity_multiIndex &dofs, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) 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 69 of file ForcesAndSourcesCore.cpp.

71  {
73 
74  auto &side_table = numeredEntFiniteElementPtr->getSideNumberTable().get<2>();
75  for (auto r = side_table.equal_range(type); r.first != r.second; ++r.first) {
76  const int side_number = (*r.first)->side_number;
77  if (side_number >= 0) {
78  const int brother_side_number = (*r.first)->brother_side_number;
79  const int sense = (*r.first)->sense;
80 
81  data[side_number].getSense() = sense;
82  if (brother_side_number != -1)
83  data[brother_side_number].getSense() = sense;
84  }
85  }
87 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
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:412

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

531  {
532  return getEntitySense(type, data.dataOnEntities[type]);
533  }
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 661 of file ForcesAndSourcesCore.cpp.

661  {
663  // PetscAttachDebugger();
664  data.facesNodes.resize(4, 3, false);
665  auto &side_table = const_cast<SideNumber_multiIndex &>(
666  numeredEntFiniteElementPtr->getSideNumberTable());
667  auto siit = side_table.get<1>().lower_bound(boost::make_tuple(MBTRI, 0));
668  auto hi_siit = side_table.get<1>().upper_bound(boost::make_tuple(MBTRI, 4));
669  if (std::distance(siit, hi_siit) != 4) {
670  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
671  "Should be 4 triangles on tet, side_table not initialized");
672  }
673  const int canonical_face_sense_p1[4][3] = {
674  {0, 1, 3},
675  {1, 2, 3},
676  {0, 3, 2} /**/,
677  {0, 2, 1} /**/}; // second index is offset (positive sense)
678  const int canonical_face_sense_m1[4][3] = {
679  {0, 3, 1},
680  {1, 3, 2},
681  {0, 2, 3},
682  {0, 1, 2}}; // second index is offset (negative sense
683  for (; siit != hi_siit; siit++) {
684  const boost::shared_ptr<SideNumber> side = *siit;
685  int face_conn[3] = {-1, -1, -1};
686  if (side->offset == 0) {
687  face_conn[0] = side->sense == 1
688  ? canonical_face_sense_p1[(int)side->side_number][0]
689  : canonical_face_sense_m1[(int)side->side_number][0];
690  face_conn[1] = side->sense == 1
691  ? canonical_face_sense_p1[(int)side->side_number][1]
692  : canonical_face_sense_m1[(int)side->side_number][1];
693  face_conn[2] = side->sense == 1
694  ? canonical_face_sense_p1[(int)side->side_number][2]
695  : canonical_face_sense_m1[(int)side->side_number][2];
696  }
697  if (side->offset == 1) {
698  face_conn[0] =
699  side->sense == 1
700  ? canonical_face_sense_p1[(int)side->side_number][1]
701  : canonical_face_sense_m1[(int)side->side_number][2] /**/;
702  face_conn[1] = side->sense == 1
703  ? canonical_face_sense_p1[(int)side->side_number][2]
704  : canonical_face_sense_m1[(int)side->side_number][0];
705  face_conn[2] = side->sense == 1
706  ? canonical_face_sense_p1[(int)side->side_number][0]
707  : canonical_face_sense_m1[(int)side->side_number][1];
708  }
709  if (side->offset == 2) {
710  face_conn[0] =
711  side->sense == 1
712  ? canonical_face_sense_p1[(int)side->side_number][2]
713  : canonical_face_sense_m1[(int)side->side_number][1] /**/;
714  face_conn[1] = side->sense == 1
715  ? canonical_face_sense_p1[(int)side->side_number][0]
716  : canonical_face_sense_m1[(int)side->side_number][2];
717  face_conn[2] = side->sense == 1
718  ? canonical_face_sense_p1[(int)side->side_number][1]
719  : canonical_face_sense_m1[(int)side->side_number][0];
720  }
721  for (int nn = 0; nn < 3; nn++)
722  data.facesNodes(side->side_number, nn) = face_conn[nn];
723  {
724  const EntityHandle *conn_tet;
725  int num_nodes_tet;
727  CHKERR mField.get_moab().get_connectivity(ent, conn_tet, num_nodes_tet,
728  true);
729  if (num_nodes_tet != 4)
730  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
731  "data inconsistency");
732  int num_nodes_face;
733  const EntityHandle *conn_face;
734  CHKERR mField.get_moab().get_connectivity(side->ent, conn_face,
735  num_nodes_face, true);
736  if (num_nodes_face != 3)
737  SETERRQ(PETSC_COMM_SELF, 1, "data inconsistency");
738  if (conn_face[0] != conn_tet[data.facesNodes(side->side_number, 0)])
739  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
740  "data inconsistency");
741  if (conn_face[1] != conn_tet[data.facesNodes(side->side_number, 1)])
742  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
743  "data inconsistency");
744  if (conn_face[2] != conn_tet[data.facesNodes(side->side_number, 2)])
745  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
746  "data inconsistency");
747  }
748  }
750 }
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:482
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ getMaxColOrder()

int MoFEM::ForcesAndSourcesCore::getMaxColOrder ( ) const

Get max order of approximation for field in columns.

Definition at line 115 of file ForcesAndSourcesCore.cpp.

115  {
116  return getMaxOrder(*colFieldEntsPtr);
117 }
boost::shared_ptr< const FieldEntity_vector_view > colFieldEntsPtr
Pointer to finite element field entities column view.
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 102 of file ForcesAndSourcesCore.cpp.

102  {
103  int max_order = 0;
104  for (auto e : *dataFieldEntsPtr) {
105  const int order = e->getMaxOrder();
106  max_order = (max_order < order) ? order : max_order;
107  }
108  return max_order;
109 }
constexpr int order
boost::shared_ptr< const FieldEntity_multiIndex_spaceType_view > dataFieldEntsPtr
Pointer to finite element field entities data view.

◆ getMaxRowOrder()

int MoFEM::ForcesAndSourcesCore::getMaxRowOrder ( ) const

Get max order of approximation for field in rows.

Definition at line 111 of file ForcesAndSourcesCore.cpp.

111  {
112  return getMaxOrder(*rowFieldEntsPtr);
113 }
boost::shared_ptr< const FieldEntity_vector_view > rowFieldEntsPtr
Pointer to finite element field entities row view.
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 454 of file ForcesAndSourcesCore.cpp.

455  {
456 
457  auto get_nodes_field_data = [&](FEDofEntity_multiIndex &dofs,
458  VectorDouble &nodes_data,
459  VectorDofs &nodes_dofs, FieldSpace &space,
461  VectorInt &bb_node_order) {
463 
464  auto &dofs_by_name_and_type = dofs.get<Composite_Name_And_Type_mi_tag>();
465  auto tuple = boost::make_tuple(field_name, MBVERTEX);
466  auto dit = dofs_by_name_and_type.lower_bound(tuple);
467  if (dit == dofs_by_name_and_type.end())
468  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
469  "No nodal dofs on element");
470  auto hi_dit = dofs.get<Composite_Name_And_Type_mi_tag>().upper_bound(tuple);
471 
472  if (dit != hi_dit) {
473  auto &first_dof = **dit;
474  space = first_dof.getSpace();
475  base = first_dof.getApproxBase();
476  int num_nodes;
477  CHKERR getNumberOfNodes(num_nodes);
478  bb_node_order.resize(num_nodes, false);
479  bb_node_order.clear();
480  const int nb_dof_idx = first_dof.getNbOfCoeffs();
481  const int max_nb_dofs = nb_dof_idx * num_nodes;
482  nodes_data.resize(max_nb_dofs, false);
483  nodes_dofs.resize(max_nb_dofs, false);
484  nodes_data.clear();
485 
486  std::vector<boost::weak_ptr<FEDofEntity>> brother_dofs_vec;
487  for (; dit != hi_dit;) {
488  const auto &dof_ptr = *dit;
489  const auto &dof = *dof_ptr;
490  const auto &sn = *dof.sideNumberPtr;
491  const int side_number = sn.side_number;
492  const int brother_side_number = sn.brother_side_number;
493  if (brother_side_number != -1)
494  brother_dofs_vec.emplace_back(dof_ptr);
495 
496  bb_node_order[side_number] = dof.getMaxOrder();
497  int pos = side_number * nb_dof_idx;
498  auto ent_filed_data_vec = dof.getEntFieldData();
499  for (int ii = 0; ii != nb_dof_idx; ++ii) {
500  nodes_data[pos] = ent_filed_data_vec[ii];
501  nodes_dofs[pos] = *dit;
502  ++pos;
503  ++dit;
504  }
505  }
506 
507  for (auto &dof_ptr : brother_dofs_vec) {
508  if (const auto d = dof_ptr.lock()) {
509  const auto &sn = d->sideNumberPtr;
510  const int side_number = sn->side_number;
511  const int brother_side_number = sn->brother_side_number;
512  bb_node_order[brother_side_number] = bb_node_order[side_number];
513  int pos = side_number * nb_dof_idx;
514  int brother_pos = brother_side_number * nb_dof_idx;
515  for (int ii = 0; ii != nb_dof_idx; ++ii) {
516  nodes_data[brother_pos] = nodes_data[pos];
517  nodes_dofs[brother_pos] = nodes_dofs[pos];
518  ++pos;
519  ++brother_pos;
520  }
521  }
522  }
523 
524  } else {
525  nodes_data.resize(0, false);
526  nodes_dofs.resize(0, false);
527  }
528 
530  };
531 
532  return get_nodes_field_data(
533  const_cast<FEDofEntity_multiIndex &>(
534  numeredEntFiniteElementPtr->getDataDofs()),
535  data.dataOnEntities[MBVERTEX][0].getFieldData(),
536  data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
537  data.dataOnEntities[MBVERTEX][0].getSpace(),
538  data.dataOnEntities[MBVERTEX][0].getBase(),
539  data.dataOnEntities[MBVERTEX][0].getBBNodeOrder());
540 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
const Tensor1_Expr< const dTensor0< T, Dim, i >, typename promote< T, double >::V, Dim, i > d(const Tensor0< T * > &a, const Index< i, Dim > index, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: dTensor0.hpp:27
FieldApproximationBase
approximation base
Definition: definitions.h:149
FieldSpace
approximation spaces
Definition: definitions.h:173
#define CHKERR
Inline error check.
Definition: definitions.h:601
multi_index_container< boost::shared_ptr< FEDofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< FEDofEntity::interface_type_DofEntity, const UId &, &FEDofEntity::getGlobalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FEDofEntity::interface_type_DofEntity, EntityHandle, &FEDofEntity::getEnt > >, ordered_non_unique< tag< FieldName_mi_tag >, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef > >, ordered_non_unique< tag< EntType_mi_tag >, const_mem_fun< FEDofEntity::interface_type_RefEntity, EntityType, &FEDofEntity::getEntType > >, ordered_non_unique< tag< Composite_Name_And_Type_mi_tag >, composite_key< FEDofEntity, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef >, const_mem_fun< FEDofEntity::interface_type_RefEntity, EntityType, &FEDofEntity::getEntType > > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< FEDofEntity, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef >, const_mem_fun< FEDofEntity::interface_type_DofEntity, EntityHandle, &FEDofEntity::getEnt > > >, ordered_non_unique< tag< Composite_Name_Type_And_Side_Number_mi_tag >, composite_key< FEDofEntity, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef >, const_mem_fun< FEDofEntity::interface_type_RefEntity, EntityType, &FEDofEntity::getEntType >, KeyFromKey< member< SideNumber, char, &SideNumber::side_number >, member< FEDofEntity::BaseFEEntity, boost::shared_ptr< SideNumber >, &FEDofEntity::sideNumberPtr > > > > > > FEDofEntity_multiIndex
MultiIndex container keeps FEDofEntity.
ublas::vector< int, IntAllocator > VectorInt
Definition: Types.hpp:72
ublas::vector< boost::shared_ptr< const FEDofEntity >, DofsAllocator > VectorDofs
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:73
MoFEMErrorCode getNumberOfNodes(int &num_nodes) const
Get number of DOFs on element.

◆ getNodesIndices()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNodesIndices ( const boost::string_ref  field_name,
FENumeredDofEntity_multiIndex dofs,
VectorInt nodes_indices,
VectorInt local_nodes_indices 
) const
protected

get node indices

Definition at line 168 of file ForcesAndSourcesCore.cpp.

170  {
172  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
173  auto tuple = boost::make_tuple(field_name, MBVERTEX);
174  auto dit = dofs_by_type.lower_bound(tuple);
175  auto hi_dit = dofs_by_type.upper_bound(tuple);
176 
177  if (dit != hi_dit) {
178 
179  int num_nodes;
180  CHKERR getNumberOfNodes(num_nodes);
181  int max_nb_dofs = 0;
182  const int nb_dofs_on_vert = (*dit)->getNbOfCoeffs();
183  max_nb_dofs = nb_dofs_on_vert * num_nodes;
184  nodes_indices.resize(max_nb_dofs, false);
185  local_nodes_indices.resize(max_nb_dofs, false);
186  if (std::distance(dit, hi_dit) != max_nb_dofs) {
187  std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
188  std::fill(local_nodes_indices.begin(), local_nodes_indices.end(), -1);
189  }
190 
191  for (; dit != hi_dit; dit++) {
192  auto &dof = **dit;
193  const int idx = dof.getPetscGlobalDofIdx();
194  const int local_idx = dof.getPetscLocalDofIdx();
195  const int side_number = dof.sideNumberPtr->side_number;
196  const int pos = side_number * nb_dofs_on_vert + dof.getDofCoeffIdx();
197  nodes_indices[pos] = idx;
198  local_nodes_indices[pos] = local_idx;
199  const int brother_side_number =
200  (*dit)->sideNumberPtr->brother_side_number;
201  if (brother_side_number != -1) {
202  const int elem_idx =
203  brother_side_number * nb_dofs_on_vert + (*dit)->getDofCoeffIdx();
204  nodes_indices[elem_idx] = idx;
205  local_nodes_indices[elem_idx] = local_idx;
206  }
207  }
208 
209  } else {
210  nodes_indices.resize(0, false);
211  local_nodes_indices.resize(0, false);
212  }
213 
215 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412
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 323 of file ForcesAndSourcesCore.cpp.

324  {
326  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
327  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
328  }
329  CHKERR getNoFieldIndices(field_name,
330  const_cast<FENumeredDofEntity_multiIndex &>(
331  numeredEntFiniteElementPtr->getColsDofs()),
332  data.dataOnEntities[MBENTITYSET][0].getIndices());
334 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
#define CHKERR
Inline error check.
Definition: definitions.h:601
MoFEMErrorCode getNoFieldIndices(const std::string &field_name, FENumeredDofEntity_multiIndex &dofs, VectorInt &nodes_indices) const
get NoField indices
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ getNoFieldFieldData() [1/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldFieldData ( const boost::string_ref  field_name,
FEDofEntity_multiIndex dofs,
VectorDouble ent_field_data,
VectorDofs ent_field_dofs 
) const
protected

Get field data on nodes.

Definition at line 626 of file ForcesAndSourcesCore.cpp.

628  {
630  auto dit = dofs.get<FieldName_mi_tag>().lower_bound(field_name);
631  auto hi_dit = dofs.get<FieldName_mi_tag>().upper_bound(field_name);
632  int size = std::distance(dit, hi_dit);
633  ent_field_data.resize(size, false);
634  ent_field_dofs.resize(size, false);
635  for (; dit != hi_dit; dit++) {
636  int idx = (*dit)->getDofCoeffIdx();
637  ent_field_data[idx] = (*dit)->getFieldData();
638  ent_field_dofs[idx] = *dit;
639  }
641 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513

◆ getNoFieldFieldData() [2/2]

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::getNoFieldFieldData ( DataForcesAndSourcesCore data,
const boost::string_ref  field_name 
) const
protected

Definition at line 643 of file ForcesAndSourcesCore.cpp.

644  {
646  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
647  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
648  }
650  field_name,
651  const_cast<FEDofEntity_multiIndex &>(
652  numeredEntFiniteElementPtr->getDataDofs()),
653  data.dataOnEntities[MBENTITYSET][0].getFieldData(),
654  data.dataOnEntities[MBENTITYSET][0].getFieldDofs());
656 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
MoFEMErrorCode getNoFieldFieldData(const boost::string_ref field_name, FEDofEntity_multiIndex &dofs, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs) const
Get field data on nodes.
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
#define CHKERR
Inline error check.
Definition: definitions.h:601
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ getNoFieldIndices()

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

get NoField indices

Definition at line 296 of file ForcesAndSourcesCore.cpp.

298  {
300  auto dit = dofs.get<FieldName_mi_tag>().lower_bound(field_name);
301  auto hi_dit = dofs.get<FieldName_mi_tag>().upper_bound(field_name);
302  indices.resize(std::distance(dit, hi_dit));
303  for (; dit != hi_dit; dit++) {
304  int idx = (*dit)->getPetscGlobalDofIdx();
305  indices[(*dit)->getDofCoeffIdx()] = idx;
306  }
308 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513

◆ getNoFieldRowIndices()

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

get col NoField indices

Definition at line 310 of file ForcesAndSourcesCore.cpp.

311  {
313  if (data.dataOnEntities[MBENTITYSET].size() == 0) {
314  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "data inconsistency");
315  }
316  CHKERR getNoFieldIndices(field_name,
317  const_cast<FENumeredDofEntity_multiIndex &>(
318  numeredEntFiniteElementPtr->getRowsDofs()),
319  data.dataOnEntities[MBENTITYSET][0].getIndices());
321 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
#define CHKERR
Inline error check.
Definition: definitions.h:601
MoFEMErrorCode getNoFieldIndices(const std::string &field_name, FENumeredDofEntity_multiIndex &dofs, VectorInt &nodes_indices) const
get NoField indices
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:412

◆ 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_skeleton_3d.cpp, continuity_check_on_skeleton_with_simple_2d.cpp, elasticity_mixed_formulation.cpp, EshelbianPlasticity.cpp, forces_and_sources_testing_edge_element.cpp, forces_and_sources_testing_flat_prism_element.cpp, MagneticElement.hpp, and simple_contact.cpp.

Definition at line 420 of file ForcesAndSourcesCore.hpp.

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

438  {
439  return getProblemNodesIndices(field_name, *(problemPtr->numeredDofsCols),
440  nodes_indices);
441 }
const Problem * problemPtr
raw pointer to problem
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsCols
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 338 of file ForcesAndSourcesCore.cpp.

340  {
342 
343  const Field *field_struture = mField.get_field_structure(field_name);
344  if (field_struture->getSpace() == H1) {
345 
346  int num_nodes;
347  CHKERR getNumberOfNodes(num_nodes);
348  nodes_indices.resize(field_struture->getNbOfCoeffs() * num_nodes, false);
349  std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
350 
351  auto &side_table = const_cast<SideNumber_multiIndex &>(
352  numeredEntFiniteElementPtr->getSideNumberTable());
353  auto siit = side_table.get<1>().lower_bound(boost::make_tuple(MBVERTEX, 0));
354  auto hi_siit =
355  side_table.get<1>().lower_bound(boost::make_tuple(MBVERTEX, 10000));
356 
357  int nn = 0;
358  for (; siit != hi_siit; siit++, nn++) {
359 
360  if (siit->get()->side_number == -1)
361  continue;
362 
363  const EntityHandle ent = siit->get()->ent;
364  auto dit =
365  dofs.get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().lower_bound(
366  boost::make_tuple(field_name, ent, 0));
367  auto hi_dit =
368  dofs.get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().upper_bound(
369  boost::make_tuple(field_name, ent, 10000)); /// very large number
370  for (; dit != hi_dit; dit++) {
371  nodes_indices[siit->get()->side_number * (*dit)->getNbOfCoeffs() +
372  (*dit)->getDofCoeffIdx()] =
373  (*dit)->getPetscGlobalDofIdx();
374  }
375  }
376  } else {
377  nodes_indices.resize(0, false);
378  }
379 
380 
382 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
#define CHKERR
Inline error check.
Definition: definitions.h:601
continuous field
Definition: definitions.h:176
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 423 of file ForcesAndSourcesCore.cpp.

424  {
425  return getProblemNodesIndices(field_name, *(problemPtr->numeredDofsRows),
426  nodes_indices);
427 }
const Problem * problemPtr
raw pointer to problem
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsRows
store DOFs on rows 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

◆ getProblemTypeColIndices()

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

Definition at line 444 of file ForcesAndSourcesCore.cpp.

446  {
447  return getProblemTypeIndices(field_name, *(problemPtr->numeredDofsCols), type,
448  side_number, indices);
449 }
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 > numeredDofsCols
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 384 of file ForcesAndSourcesCore.cpp.

386  {
388 
389  indices.resize(0);
390 
391  auto &side_table = const_cast<SideNumber_multiIndex &>(
392  numeredEntFiniteElementPtr->getSideNumberTable());
393  auto siit =
394  side_table.get<1>().lower_bound(boost::make_tuple(type, side_number));
395  auto hi_siit =
396  side_table.get<1>().upper_bound(boost::make_tuple(type, side_number));
397 
398  for (; siit != hi_siit; siit++) {
399 
400  if (siit->get()->side_number == -1)
401  continue;
402 
403  const EntityHandle ent = siit->get()->ent;
404  NumeredDofEntity_multiIndex::index<
405  Composite_Name_And_Ent_And_EntDofIdx_mi_tag>::type::iterator dit,
406  hi_dit;
407  dit = dofs.get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().lower_bound(
408  boost::make_tuple(field_name, ent, 0));
409  hi_dit =
410  dofs.get<Composite_Name_And_Ent_And_EntDofIdx_mi_tag>().upper_bound(
411  boost::make_tuple(field_name, ent, 10000)); /// very large number
412 
413  indices.resize(std::distance(dit, hi_dit));
414  for (; dit != hi_dit; dit++) {
415 
416  indices[(*dit)->getEntDofIdx()] = (*dit)->getPetscGlobalDofIdx();
417  }
418  }
419 
421 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr

◆ getProblemTypeRowIndices()

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

Definition at line 430 of file ForcesAndSourcesCore.cpp.

432  {
433  return getProblemTypeIndices(field_name, *(problemPtr->numeredDofsRows), type,
434  side_number, indices);
435 }
const Problem * problemPtr
raw pointer to problem
boost::shared_ptr< NumeredDofEntity_multiIndex > numeredDofsRows
store DOFs on rows for this 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 218 of file ForcesAndSourcesCore.cpp.

219  {
220  return getNodesIndices(field_name,
221  const_cast<FENumeredDofEntity_multiIndex &>(
222  numeredEntFiniteElementPtr->getRowsDofs()),
223  data.dataOnEntities[MBVERTEX][0].getIndices(),
224  data.dataOnEntities[MBVERTEX][0].getLocalIndices());
225 }
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
MoFEMErrorCode getNodesIndices(const boost::string_ref field_name, FENumeredDofEntity_multiIndex &dofs, VectorInt &nodes_indices, VectorInt &local_nodes_indices) const
get node indices

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

1354  {
1355  return getRuleHook ? getRuleHook(order_row, order_col, order_data)
1356  : getRule(order_data);
1357 }
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 754 of file ForcesAndSourcesCore.cpp.

755  {
757 
758  if (nInTheLoop == 0) {
759  data.sPace.reset();
760  data.bAse.reset();
761  for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
762  data.spacesOnEntities[t].reset();
763  data.basesOnEntities[t].reset();
764  }
765  for (int s = 0; s != LASTSPACE; ++s) {
766  data.basesOnSpaces[s].reset();
767  }
768  }
769 
770  if (dataFieldEntsPtr)
771  for (auto e : *dataFieldEntsPtr) {
772  // get data from entity
773  const EntityType type = e->getEntType();
774  const FieldSpace space = e->getSpace();
775  const FieldApproximationBase approx = e->getApproxBase();
776 
777  // set data
778  data.sPace.set(space);
779  data.bAse.set(approx);
780  data.spacesOnEntities[type].set(space);
781  data.basesOnEntities[type].set(approx);
782  data.basesOnSpaces[space].set(approx);
783  }
784  else
785  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
786  "data fields ents not allocated on element");
787 
789 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:180
FieldApproximationBase
approximation base
Definition: definitions.h:149
FieldSpace
approximation spaces
Definition: definitions.h:173
int nInTheLoop
number currently of processed method
boost::shared_ptr< const FieldEntity_multiIndex_spaceType_view > dataFieldEntsPtr
Pointer to finite element field entities data view.

◆ getUserPolynomialBase()

auto& MoFEM::ForcesAndSourcesCore::getUserPolynomialBase ( )

Get the User Polynomial Base object.

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

Definition at line 434 of file ForcesAndSourcesCore.hpp.

434 { 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 1063 of file ForcesAndSourcesCore.cpp.

1063  {
1065 
1066  const EntityType type = numeredEntFiniteElementPtr->getEntType();
1069  std::vector<std::string> last_eval_field_name(2);
1070 
1071  boost::ptr_vector<UserDataOperator>::iterator oit, hi_oit;
1072  oit = opPtrVector.begin();
1073  hi_oit = opPtrVector.end();
1074 
1075  for (; oit != hi_oit; oit++) {
1076 
1077  try {
1078 
1079  CHKERR oit->setPtrFE(this);
1080 
1081  if (oit->opType == UserDataOperator::OPLAST) {
1082 
1083  // Set field
1084  switch (oit->sPace) {
1085  case NOSPACE:
1086  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Unknown space");
1087  case NOFIELD:
1088  case H1:
1089  case HCURL:
1090  case HDIV:
1091  case L2:
1092  break;
1093  default:
1094  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1095  "Not implemented for this space", oit->sPace);
1096  }
1097 
1098  // Reseat all data which all field dependent
1099  dataOnElement[oit->sPace]->resetFieldDependentData();
1100  last_eval_field_name[0] = "";
1101 
1102  // Run operator
1103  try {
1104  CHKERR oit->opRhs(*dataOnElement[oit->sPace], false);
1105  }
1106  CATCH_OP_ERRORS(*oit);
1107 
1108  } else {
1109 
1110  boost::shared_ptr<DataForcesAndSourcesCore> op_data[2];
1111  std::array<bool, 2> base_swap;
1112  std::array<std::pair<std::string, FieldApproximationBase>, 2>
1113  base_swap_data;
1114  auto swap_bases = [&]() {
1116  for (size_t ss = 0; ss != 2; ++ss)
1117  if (base_swap[ss])
1118  CHKERR op_data[ss]->baseSwap(base_swap_data[ss].first,
1119  base_swap_data[ss].second);
1121  };
1122 
1123  for (size_t ss = 0; ss != 2; ss++) {
1124 
1125  const std::string field_name =
1126  !ss ? oit->rowFieldName : oit->colFieldName;
1127  if (field_name.empty()) {
1128  SETERRQ2(
1129  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1130  "No field name in operator %d (0-row, 1-column) in operator %s",
1131  ss,
1132  (boost::typeindex::type_id_runtime(*oit).pretty_name())
1133  .c_str());
1134  }
1135  const Field *field_struture = mField.get_field_structure(field_name);
1136  const BitFieldId data_id = field_struture->getId();
1137  const FieldSpace space = field_struture->getSpace();
1138  const FieldApproximationBase base = field_struture->getApproxBase();
1139  op_data[ss] =
1140  !ss ? dataOnElement[space] : derivedDataOnElement[space];
1141 
1142  switch (base) {
1144  base_swap_data[ss] = std::pair<std::string, FieldApproximationBase>(
1145  field_name, AINSWORTH_BERNSTEIN_BEZIER_BASE);
1146  base_swap[ss] = true;
1147  break;
1148  default:
1149  base_swap[ss] = false;
1150  };
1151 
1152  if ((oit->getNumeredEntFiniteElementPtr()->getBitFieldIdData() &
1153  data_id)
1154  .none()) {
1155  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1156  "no data field < %s > on finite element < %s >",
1157  field_name.c_str(), feName.c_str());
1158  }
1159 
1160  if (oit->getOpType() & types[ss] ||
1161  oit->getOpType() & UserDataOperator::OPROWCOL) {
1162 
1163  switch (space) {
1164  case NOSPACE:
1165  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1166  "unknown space");
1167  break;
1168  case NOFIELD:
1169  case H1:
1170  case HCURL:
1171  case HDIV:
1172  case L2:
1173  break;
1174  default:
1175  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1176  "Not implemented for this space", space);
1177  }
1178 
1179  if (last_eval_field_name[ss] != field_name) {
1180 
1181  CHKERR getEntityFieldData(*op_data[ss], field_name, MBEDGE);
1182  if (!ss)
1183  CHKERR getEntityRowIndices(*op_data[ss], field_name, MBEDGE);
1184  else
1185  CHKERR getEntityColIndices(*op_data[ss], field_name, MBEDGE);
1186 
1187  switch (space) {
1188  case NOSPACE:
1189  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
1190  "unknown space");
1191  break;
1192  case H1:
1193  if (!ss)
1194  CHKERR getRowNodesIndices(*op_data[ss], field_name);
1195  else
1196  CHKERR getColNodesIndices(*op_data[ss], field_name);
1197  CHKERR getNodesFieldData(*op_data[ss], field_name);
1198  break;
1199  case HCURL:
1200  case HDIV:
1201  break;
1202  case L2:
1203  switch (type) {
1204  case MBVERTEX:
1205  CHKERR getNodesFieldData(*op_data[ss], field_name);
1206  break;
1207  default:
1208  break;
1209  }
1210  break;
1211  case NOFIELD:
1212  if (!getNinTheLoop()) {
1213  // NOFIELD data are the same for each element, can be
1214  // retrieved only once
1215  if (!ss) {
1216  CHKERR getNoFieldRowIndices(*op_data[ss], field_name);
1217  } else {
1218  CHKERR getNoFieldColIndices(*op_data[ss], field_name);
1219  }
1220  CHKERR getNoFieldFieldData(*op_data[ss], field_name);
1221  }
1222  break;
1223  default:
1224  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
1225  "Not implemented for this space", space);
1226  }
1227  last_eval_field_name[ss] = field_name;
1228  }
1229  }
1230  }
1231 
1232  CHKERR swap_bases();
1233 
1234  if (oit->getOpType() & UserDataOperator::OPROW) {
1235  try {
1236  CHKERR oit->opRhs(*op_data[0], false);
1237  }
1238  CATCH_OP_ERRORS(*oit);
1239  }
1240 
1241  if (oit->getOpType() & UserDataOperator::OPCOL) {
1242  try {
1243  CHKERR oit->opRhs(*op_data[1], false);
1244  }
1245  CATCH_OP_ERRORS(*oit);
1246  }
1247 
1248  if (oit->getOpType() & UserDataOperator::OPROWCOL) {
1249  try {
1250  CHKERR oit->opLhs(*op_data[0], *op_data[1]);
1251  }
1252  CATCH_OP_ERRORS(*oit);
1253  }
1254 
1255  CHKERR swap_bases();
1256  }
1257  }
1258  CATCH_OP_ERRORS(*oit);
1259  }
1261 }
boost::ptr_vector< UserDataOperator > opPtrVector
Vector of finite element users data operators.
field with continuous normal traction
Definition: definitions.h:178
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:506
scalar or vector of scalars describe (no true field)
Definition: definitions.h:175
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:482
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:513
MoFEMErrorCode getNoFieldRowIndices(DataForcesAndSourcesCore &data, const std::string &field_name) const
get col NoField indices
OpType
Controls loop over entities on element.
MoFEMErrorCode getNoFieldFieldData(const boost::string_ref field_name, FEDofEntity_multiIndex &dofs, VectorDouble &ent_field_data, VectorDofs &ent_field_dofs) const
Get field data on nodes.
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)
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnElement
Entity data on element entity columns fields.
FieldApproximationBase
approximation base
Definition: definitions.h:149
field with continuous tangents
Definition: definitions.h:177
FieldSpace
approximation spaces
Definition: definitions.h:173
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
Definition: Types.hpp:52
#define CHKERR
Inline error check.
Definition: definitions.h:601
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:412
continuous field
Definition: definitions.h:176
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:179

◆ 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::VolumeElementForcesAndSourcesCoreOnSideSwitch< SWITCH >, MoFEM::FatPrismElementForcesAndSourcesCore, MoFEM::ContactPrismElementForcesAndSourcesCore, MoFEM::EdgeElementForcesAndSourcesCoreSwitch< SWITCH >, MoFEM::FaceElementForcesAndSourcesCoreOnSideSwitch< SWITCH >, MoFEM::FaceElementForcesAndSourcesCoreOnSideSwitch< FaceElementForcesAndSourcesCore::NO_HO_GEOMETRY >, MoFEM::FlatPrismElementForcesAndSourcesCore, NeumannForcesSurfaceComplexForLazy::MyTriangleSpatialFE, and MoFEM::VertexElementForcesAndSourcesCore.

Definition at line 1399 of file ForcesAndSourcesCore.cpp.

1399  {
1401  if (operatorHook) {
1402  ierr = operatorHook();
1403  CHKERRG(ierr);
1404  }
1406 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:549
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
boost::function< MoFEMErrorCode()> operatorHook
Hook function for operator.

◆ postProcess()

MoFEMErrorCode MoFEM::ForcesAndSourcesCore::postProcess ( )
virtual

function is run at the end of loop

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

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

Reimplemented from MoFEM::BasicMethod.

Reimplemented in 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 1407 of file ForcesAndSourcesCore.cpp.

1407  {
1409  if (postProcessHook) {
1410  ierr = postProcessHook();
1411  CHKERRG(ierr);
1412  }
1414 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:549
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
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 EdgeSlidingConstrains::MyEdgeFE, PostProcEdgeOnRefinedMesh, PostProcFatPrismOnRefinedMesh, PostProcTemplateVolumeOnRefinedMesh< MoFEM::VolumeElementForcesAndSourcesCore >, GelModule::Gel::GelFE, SurfaceSlidingConstrains::MyTriangleFE, KelvinVoigtDamper::DamperFE, SolidShellModule::SolidShellPrismElement::SolidShellError, SolidShellModule::SolidShellPrismElement::SolidShell, NeumannForcesSurfaceComplexForLazy::MyTriangleSpatialFE, Smoother::MyVolumeFE, and NonlinearElasticElement::MyVolumeFE.

Definition at line 1391 of file ForcesAndSourcesCore.cpp.

1391  {
1393  if (preProcessHook) {
1394  ierr = preProcessHook();
1395  CHKERRG(ierr);
1396  }
1398 }
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:506
#define CHKERRG(n)
Check error code of MoFEM/MOAB/PETSc function.
Definition: definitions.h:549
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513
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 1359 of file ForcesAndSourcesCore.cpp.

1360  {
1361  return setRuleHook ? setRuleHook(order_row, order_col, order_data)
1362  : setGaussPts(order_data);
1363 }
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, and MoFEM::FaceElementForcesAndSourcesCoreOnSideBase.

Definition at line 1369 of file ForcesAndSourcesCore.cpp.

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

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

1298  {
1300  sidePtrFE = const_cast<ForcesAndSourcesCore *>(side_fe_ptr);
1302 }
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:506
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:513

Friends And Related Function Documentation

◆ FaceElementForcesAndSourcesCoreOnSideBase

Definition at line 847 of file ForcesAndSourcesCore.hpp.

◆ UserDataOperator

friend class UserDataOperator
friend

◆ VolumeElementForcesAndSourcesCoreOnSideBase

Definition at line 846 of file ForcesAndSourcesCore.hpp.

Member Data Documentation

◆ dataH1

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataH1
protected

Definition at line 797 of file ForcesAndSourcesCore.hpp.

◆ dataHcurl

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataHcurl
protected

Definition at line 798 of file ForcesAndSourcesCore.hpp.

◆ dataHdiv

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataHdiv
protected

Definition at line 799 of file ForcesAndSourcesCore.hpp.

◆ dataL2

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataL2
protected

Definition at line 800 of file ForcesAndSourcesCore.hpp.

◆ dataNoField

DataForcesAndSourcesCore& MoFEM::ForcesAndSourcesCore::dataNoField
protected

Definition at line 796 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 787 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 794 of file ForcesAndSourcesCore.hpp.

◆ elementPolynomialBasePtr

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

Pointer to entity polynomial base.

Definition at line 822 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 446 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.

Definition at line 50 of file ForcesAndSourcesCore.hpp.

◆ lastEvaluatedElementEntityType

EntityType MoFEM::ForcesAndSourcesCore::lastEvaluatedElementEntityType
protected

Last evaluated type of element entity.

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

◆ userPolynomialBasePtr

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

Pointer to user polynomail base.

Definition at line 827 of file ForcesAndSourcesCore.hpp.


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