v0.9.1
Classes | Public Member Functions | Protected Member Functions | Protected Attributes | Private Attributes | List of all members
MoFEM::ContactPrismElementForcesAndSourcesCore Struct Reference

ContactPrism finite elementUser is implementing own operator at Gauss points level, by own class derived from ContactPrismElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary number of operator added pushing instances to rowOpPtrVector and rowColOpPtrVector. More...

#include <src/finite_elements/ContactPrismElementForcesAndSourcesCore.hpp>

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

Classes

struct  UserDataOperator
 default operator for Contact Prism element More...
 

Public Member Functions

 ContactPrismElementForcesAndSourcesCore (Interface &m_field)
 
MoFEMErrorCode operator() ()
 function is run for every finite element More...
 
- Public Member Functions inherited from MoFEM::ForcesAndSourcesCore
 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 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)
 

Protected Member Functions

MoFEMErrorCode setDefaultGaussPts (const int rule)
 
MoFEMErrorCode loopOverOperators ()
 Iterate user data operators. More...
 
MoFEMErrorCode getEntityFieldData (DataForcesAndSourcesCore &master_data, DataForcesAndSourcesCore &slave_data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
 function that gets entity field data. More...
 
MoFEMErrorCode getEntityIndices (DataForcesAndSourcesCore &master_data, DataForcesAndSourcesCore &slave_data, const std::string &field_name, FENumeredDofEntity_multiIndex &dofs, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
 function that gets entity indices. More...
 
MoFEMErrorCode getNodesIndices (const boost::string_ref field_name, FENumeredDofEntity_multiIndex &dofs, VectorInt &master_nodes_indices, VectorInt &master_local_nodes_indices, VectorInt &slave_nodes_indices, VectorInt &slave_local_nodes_indices) const
 function that gets nodes indices. More...
 
MoFEMErrorCode getNodesFieldData (const boost::string_ref field_name, FEDofEntity_multiIndex &dofs, VectorDouble &master_nodes_data, VectorDouble &slave_nodes_data, VectorDofs &master_nodes_dofs, VectorDofs &slave_nodes_dofs, FieldSpace &master_space, FieldSpace &slave_space, FieldApproximationBase &master_base, FieldApproximationBase &slave_base) const
 function that gets nodes field data. More...
 
- Protected Member Functions inherited from MoFEM::ForcesAndSourcesCore
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...
 
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...
 
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
 
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
 
virtual int getRule (int order)
 
virtual MoFEMErrorCode setGaussPts (int order)
 
- 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...
 

Protected Attributes

std::array< double, 2 > aRea
 Array storing master and slave faces areas. More...
 
VectorDouble normal
 vector storing vector normal to master or slave element More...
 
VectorDouble coords
 
MatrixDouble coordsAtGaussPtsMaster
 
MatrixDouble coordsAtGaussPtsSlave
 
MatrixDouble gaussPtsMaster
 
MatrixDouble gaussPtsSlave
 
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACEdataOnMaster
 Entity data on element entity rows fields. More...
 
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACEdataOnSlave
 
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACEderivedDataOnMaster
 Entity data on element entity columns fields. More...
 
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACEderivedDataOnSlave
 
DataForcesAndSourcesCoredataH1Master
 
DataForcesAndSourcesCoredataH1Slave
 
DataForcesAndSourcesCoredataNoFieldMaster
 
DataForcesAndSourcesCoredataNoFieldSlave
 
DataForcesAndSourcesCoredataHcurlMaster
 
DataForcesAndSourcesCoredataHcurlSlave
 
DataForcesAndSourcesCoredataHdivMaster
 
DataForcesAndSourcesCoredataHdivSlave
 
DataForcesAndSourcesCoredataL2Master
 
DataForcesAndSourcesCoredataL2Slave
 
- Protected Attributes inherited from MoFEM::ForcesAndSourcesCore
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...
 

Private Attributes

int nbGaussPts
 

Additional Inherited Members

- Public Types inherited from MoFEM::ForcesAndSourcesCore
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 Attributes inherited from MoFEM::ForcesAndSourcesCore
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...
 
- 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

ContactPrism finite element

User is implementing own operator at Gauss points level, by own class derived from ContactPrismElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary number of operator added pushing instances to rowOpPtrVector and rowColOpPtrVector.

Definition at line 37 of file ContactPrismElementForcesAndSourcesCore.hpp.

Constructor & Destructor Documentation

◆ ContactPrismElementForcesAndSourcesCore()

MoFEM::ContactPrismElementForcesAndSourcesCore::ContactPrismElementForcesAndSourcesCore ( Interface m_field)

Definition at line 24 of file ContactPrismElementForcesAndSourcesCore.cpp.

25  : ForcesAndSourcesCore(m_field),
27 
28  nullptr,
29  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // NOFIELD
30  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // H1
31  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // HCURL
32  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // HDIV
33  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET) // L2
34 
35  },
37 
38  nullptr,
39  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // NOFIELD
40  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // H1
41  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // HCURL
42  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET), // HDIV
43  boost::make_shared<DataForcesAndSourcesCore>(MBENTITYSET) // L2
44 
45  },
47 
48  nullptr,
49  boost::make_shared<DerivedDataForcesAndSourcesCore>(
51  boost::make_shared<DerivedDataForcesAndSourcesCore>(dataOnMaster[H1]),
52  boost::make_shared<DerivedDataForcesAndSourcesCore>(
54  boost::make_shared<DerivedDataForcesAndSourcesCore>(
56  boost::make_shared<DerivedDataForcesAndSourcesCore>(dataOnMaster[L2])
57 
58  },
60 
61  nullptr,
62  boost::make_shared<DerivedDataForcesAndSourcesCore>(
64  boost::make_shared<DerivedDataForcesAndSourcesCore>(dataOnSlave[H1]),
65  boost::make_shared<DerivedDataForcesAndSourcesCore>(
67  boost::make_shared<DerivedDataForcesAndSourcesCore>(
68  dataOnSlave[HDIV]),
69  boost::make_shared<DerivedDataForcesAndSourcesCore>(dataOnSlave[L2])
70 
71  },
72  dataH1Master(*dataOnMaster[H1].get()),
73  dataH1Slave(*dataOnSlave[H1].get()),
79  dataL2Master(*dataOnMaster[L2].get()),
81  dataL2Slave(*dataOnSlave[L2].get()) {
82 
84  boost::shared_ptr<BaseFunction>(new TriPolynomialBase());
85 
86  dataH1Master.dataOnEntities[MBVERTEX].push_back(
88  dataH1Slave.dataOnEntities[MBVERTEX].push_back(
90 
91  for (int ee = 0; ee != 3; ++ee) {
92  dataH1Master.dataOnEntities[MBEDGE].push_back(
94  dataH1Slave.dataOnEntities[MBEDGE].push_back(
96  }
97 
98  dataH1Master.dataOnEntities[MBTRI].push_back(
100  dataH1Slave.dataOnEntities[MBTRI].push_back(
102 
103  // Data on elements for proper spaces
104  dataOnMaster[H1]->setElementType(MBTRI);
105  derivedDataOnMaster[H1]->setElementType(MBTRI);
106  dataOnSlave[H1]->setElementType(MBTRI);
107  derivedDataOnSlave[H1]->setElementType(MBTRI);
108 }
ForcesAndSourcesCore(Interface &m_field)
field with continuous normal traction
Definition: definitions.h:179
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
auto & getUserPolynomialBase()
Get the User Polynomial Base object.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnMaster
Entity data on element entity columns fields.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnSlave
DataForcesAndSourcesCore::EntData EntData
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnSlave
field with continuous tangents
Definition: definitions.h:178
continuous field
Definition: definitions.h:177
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnMaster
Entity data on element entity rows fields.
field with C-1 continuity
Definition: definitions.h:180

Member Function Documentation

◆ getEntityFieldData()

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

function that gets entity field data.

Parameters
master_datadata fot master face
slave_datadata fot master face
field_namefield name of interest
type_lolowest dimension entity type to be searched
type_hihighest dimension entity type to be searched

Definition at line 676 of file ContactPrismElementForcesAndSourcesCore.cpp.

679  {
681 
682  auto reset_data = [type_lo, type_hi](auto &data) {
683  for (EntityType t = type_lo; t != type_hi; ++t) {
684  for (auto &dat : data.dataOnEntities[t]) {
685  dat.getDataOrder() = 0;
686  dat.getBase() = NOBASE;
687  dat.getSpace() = NOSPACE;
688  dat.getFieldData().resize(0, false);
689  dat.getFieldDofs().resize(0, false);
690  }
691  }
692  };
693  reset_data(master_data);
694  reset_data(slave_data);
695 
696  auto &dofs = const_cast<FEDofEntity_multiIndex &>(
697  numeredEntFiniteElementPtr->getDataDofs());
698  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
699  auto dit = dofs_by_type.lower_bound(boost::make_tuple(field_name, type_lo));
700  if (dit == dofs_by_type.end())
702  auto hi_dit =
703  dofs_by_type.lower_bound(boost::make_tuple(field_name, type_hi));
704 
705  auto get_data = [&](auto &data, auto &dof, auto type, auto side) {
706  auto &dat = data.dataOnEntities[type][side];
707  auto &ent_field_dofs = dat.getFieldDofs();
708  auto &ent_field_data = dat.getFieldData();
709  if (ent_field_data.empty()) {
710  dat.getBase() = dof.getApproxBase();
711  dat.getSpace() = dof.getSpace();
712  const int ent_order = dof.getMaxOrder();
713  dat.getDataOrder() =
714  dat.getDataOrder() > ent_order ? dat.getDataOrder() : ent_order;
715  const auto dof_ent_field_data = dof.getEntFieldData();
716  const int nb_dofs_on_ent = dof.getNbDofsOnEnt();
717  ent_field_data.resize(nb_dofs_on_ent, false);
718  noalias(ent_field_data) = dof.getEntFieldData();
719  ent_field_dofs.resize(nb_dofs_on_ent, false);
720  for (int ii = 0; ii != nb_dofs_on_ent; ++ii) {
721  const_cast<FEDofEntity *&>(ent_field_dofs[ii]) = (*dit).get();
722  ++dit;
723  }
724  }
725  };
726 
727  for (; dit != hi_dit;) {
728 
729  auto &dof = **dit;
730  if (dof.getNbDofsOnEnt()) {
731 
732  const EntityType type = dof.getEntType();
733  const int side = dof.sideNumberPtr->side_number;
734 
735  switch (type) {
736  case MBEDGE:
737 
738  if (side < 3)
739  get_data(master_data, dof, type, side);
740  else if (side > 5)
741  get_data(slave_data, dof, type, side - 6);
742 
743  break;
744  case MBTRI:
745 
746  if (side == 3)
747  get_data(master_data, dof, type, 0);
748  if (side == 4)
749  get_data(slave_data, dof, type, 0);
750 
751  break;
752  default:
753  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
754  "Entity type not implemented");
755  };
756 
757  const int brother_side = dof.sideNumberPtr->brother_side_number;
758  if (brother_side != -1)
759  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
760  "Case with brother side not implemented");
761  }
762  }
763 
765 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
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:413

◆ getEntityIndices()

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

function that gets entity indices.

Parameters
master_datadata fot master face
slave_datadata fot master face
field_namefield name of interest
dofsMultiIndex container keeping FENumeredDofEntity.
type_lolowest dimension entity type to be searched
type_hihighest dimension entity type to be searched

Definition at line 840 of file ContactPrismElementForcesAndSourcesCore.cpp.

843  {
845 
846  auto clear_data = [type_lo, type_hi](auto &data) {
847  for (EntityType t = type_lo; t != type_hi; ++t) {
848  for (auto &d : data.dataOnEntities[t]) {
849  d.getIndices().resize(0, false);
850  d.getLocalIndices().resize(0, false);
851  }
852  }
853  };
854 
855  clear_data(master_data);
856  clear_data(slave_data);
857 
858  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
859  auto dit = dofs_by_type.lower_bound(boost::make_tuple(field_name, type_lo));
860  if (dit == dofs_by_type.end())
862  auto hi_dit =
863  dofs_by_type.lower_bound(boost::make_tuple(field_name, type_hi));
864 
865  auto get_indices = [&](auto &data, auto &dof, const auto type,
866  const auto side) {
867  auto &dat = data.dataOnEntities[type][side];
868  auto &ent_field_indices = dat.getIndices();
869  auto &ent_field_local_indices = dat.getLocalIndices();
870  if (ent_field_indices.empty()) {
871  const int nb_dofs_on_ent = dof.getNbDofsOnEnt();
872  ent_field_indices.resize(nb_dofs_on_ent, false);
873  ent_field_local_indices.resize(nb_dofs_on_ent, false);
874  std::fill(ent_field_indices.data().begin(),
875  ent_field_indices.data().end(), -1);
876  std::fill(ent_field_local_indices.data().begin(),
877  ent_field_local_indices.data().end(), -1);
878  }
879  const int idx = dof.getEntDofIdx();
880  ent_field_indices[idx] = dof.getPetscGlobalDofIdx();
881  ent_field_local_indices[idx] = dof.getPetscLocalDofIdx();
882  };
883 
884  for (; dit != hi_dit; ++dit) {
885 
886  auto &dof = **dit;
887 
888  if (dof.getNbDofsOnEnt()) {
889 
890  const EntityType type = dof.getEntType();
891  const int side = dof.sideNumberPtr->side_number;
892 
893  switch (type) {
894  case MBEDGE:
895 
896  if (side < 3)
897  get_indices(master_data, dof, type, side);
898  else if (side > 5)
899  get_indices(slave_data, dof, type, side - 6);
900 
901  break;
902  case MBTRI:
903 
904  if (side == 3)
905  get_indices(master_data, dof, type, 0);
906  if (side == 4)
907  get_indices(slave_data, dof, type, 0);
908 
909  break;
910  default:
911  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
912  "Entity type not implemented");
913  };
914 
915  const int brother_side = dof.sideNumberPtr->brother_side_number;
916  if (brother_side != -1)
917  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented case");
918  }
919  }
920 
922 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
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:413

◆ getNodesFieldData()

MoFEMErrorCode MoFEM::ContactPrismElementForcesAndSourcesCore::getNodesFieldData ( const boost::string_ref  field_name,
FEDofEntity_multiIndex dofs,
VectorDouble master_nodes_data,
VectorDouble slave_nodes_data,
VectorDofs master_nodes_dofs,
VectorDofs slave_nodes_dofs,
FieldSpace master_space,
FieldSpace slave_space,
FieldApproximationBase master_base,
FieldApproximationBase slave_base 
) const
protected

function that gets nodes field data.

Parameters
field_namefield name of interest
dofsMultiIndex container keeping FENumeredDofEntity.
master_nodes_datavector containing master nodes data
slave_nodes_datavector containing master nodes data
master_nodes_dofsvector containing master nodes dofs
slave_nodes_dofsvector containing slave nodes dofs
master_spaceapproximation energy space at master
slave_spaceapproximation energy space at slave
master_basebase for master face
slave_basebase for slave face

Definition at line 767 of file ContactPrismElementForcesAndSourcesCore.cpp.

773  {
775 
776  auto set_zero = [](auto &nodes_data, auto &nodes_dofs) {
777  nodes_data.resize(0, false);
778  nodes_dofs.resize(0, false);
779  };
780  set_zero(master_nodes_data, master_nodes_dofs);
781  set_zero(slave_nodes_data, slave_nodes_dofs);
782 
783  auto &dofs_by_name_and_type = dofs.get<Composite_Name_And_Type_mi_tag>();
784  auto tuple = boost::make_tuple(field_name, MBVERTEX);
785  auto dit = dofs_by_name_and_type.lower_bound(tuple);
786  if (dit == dofs_by_name_and_type.end())
787  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
788  "No nodal dofs on element");
789  auto hi_dit = dofs.get<Composite_Name_And_Type_mi_tag>().upper_bound(tuple);
790 
791  if (dit != hi_dit) {
792  auto &first_dof = **dit;
793  const int nb_dof_idx = first_dof.getNbOfCoeffs();
794 
795  auto init_set = [&](auto &nodes_data, auto &nodes_dofs, auto &space,
796  auto &base) {
797  constexpr int num_nodes = 3;
798  const int max_nb_dofs = nb_dof_idx * num_nodes;
799  space = first_dof.getSpace();
800  base = first_dof.getApproxBase();
801  nodes_data.resize(max_nb_dofs, false);
802  nodes_dofs.resize(max_nb_dofs, false);
803  nodes_data.clear();
804  };
805 
806  init_set(master_nodes_data, master_nodes_dofs, master_space, master_base);
807  init_set(slave_nodes_data, slave_nodes_dofs, slave_space, slave_base);
808 
809  std::vector<boost::weak_ptr<FEDofEntity>> brother_dofs_vec;
810  for (; dit != hi_dit;) {
811  const auto &dof_ptr = *dit;
812  const auto &dof = *dof_ptr;
813  const auto &sn = *dof.sideNumberPtr;
814  int side = sn.side_number;
815 
816  auto set_data = [&](auto &nodes_data, auto &nodes_dofs, int pos) {
817  auto ent_filed_data_vec = dof.getEntFieldData();
818  for (int ii = 0; ii != nb_dof_idx; ++ii) {
819  nodes_data[pos] = ent_filed_data_vec[ii];
820  const_cast<FEDofEntity *&>(nodes_dofs[pos]) = (*dit).get();
821  ++pos;
822  ++dit;
823  }
824  };
825 
826  if (side < 3)
827  set_data(master_nodes_data, master_nodes_dofs, side * nb_dof_idx);
828  else
829  set_data(slave_nodes_data, slave_nodes_dofs, (side - 3) * nb_dof_idx);
830 
831  const int brother_side = sn.brother_side_number;
832  if (brother_side != -1)
833  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented");
834  }
835  }
836 
838 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ getNodesIndices()

MoFEMErrorCode MoFEM::ContactPrismElementForcesAndSourcesCore::getNodesIndices ( const boost::string_ref  field_name,
FENumeredDofEntity_multiIndex dofs,
VectorInt master_nodes_indices,
VectorInt master_local_nodes_indices,
VectorInt slave_nodes_indices,
VectorInt slave_local_nodes_indices 
) const
protected

function that gets nodes indices.

Parameters
field_namefield name of interest
dofsMultiIndex container keeping FENumeredDofEntity.
master_nodes_indicesvector containing global master nodes indices
master_local_nodes_indicesvector containing local master nodes indices
slave_nodes_indicesvector containing global master nodes indices
slave_local_nodes_indicesvector containing local master nodes indices

Definition at line 924 of file ContactPrismElementForcesAndSourcesCore.cpp.

928  {
930 
931  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
932  auto tuple = boost::make_tuple(field_name, MBVERTEX);
933  auto dit = dofs_by_type.lower_bound(tuple);
934  auto hi_dit = dofs_by_type.upper_bound(tuple);
935 
936  master_nodes_indices.resize(0, false);
937  master_local_nodes_indices.resize(0, false);
938  slave_nodes_indices.resize(0, false);
939  slave_local_nodes_indices.resize(0, false);
940 
941  if (dit != hi_dit) {
942 
943  const int num_nodes = 3;
944  int max_nb_dofs = 0;
945  const int nb_dofs_on_vert = (*dit)->getNbOfCoeffs();
946  max_nb_dofs = nb_dofs_on_vert * num_nodes;
947 
948  auto set_vec_size = [&](auto &nodes_indices, auto &local_nodes_indices) {
949  nodes_indices.resize(max_nb_dofs, false);
950  local_nodes_indices.resize(max_nb_dofs, false);
951  if (std::distance(dit, hi_dit) != max_nb_dofs) {
952  std::fill(nodes_indices.begin(), nodes_indices.end(), -1);
953  std::fill(local_nodes_indices.begin(), local_nodes_indices.end(), -1);
954  }
955  };
956 
957  set_vec_size(master_nodes_indices, master_local_nodes_indices);
958  set_vec_size(slave_nodes_indices, slave_local_nodes_indices);
959 
960  auto get_indices = [&](auto &nodes_indices, auto &local_nodes_indices,
961  auto &dof, const auto side) {
962  const int idx = dof.getPetscGlobalDofIdx();
963  const int local_idx = dof.getPetscLocalDofIdx();
964  const int pos = side * nb_dofs_on_vert + dof.getDofCoeffIdx();
965  nodes_indices[pos] = idx;
966  local_nodes_indices[pos] = local_idx;
967  };
968 
969  for (; dit != hi_dit; dit++) {
970  auto &dof = **dit;
971  const int side = dof.sideNumberPtr->side_number;
972 
973  if (side < 3)
974  get_indices(master_nodes_indices, master_local_nodes_indices, dof,
975  side);
976  else if (side > 2)
977  get_indices(slave_nodes_indices, slave_local_nodes_indices, dof,
978  side - 3);
979  else
980  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Impossible case");
981 
982  const int brother_side = (*dit)->sideNumberPtr->brother_side_number;
983  if (brother_side != -1)
984  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented case");
985  }
986  }
987 
989 }
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

◆ loopOverOperators()

MoFEMErrorCode MoFEM::ContactPrismElementForcesAndSourcesCore::loopOverOperators ( )
protected

Iterate user data operators.

Returns
MoFEMErrorCode

Definition at line 397 of file ContactPrismElementForcesAndSourcesCore.cpp.

397  {
399 
400  const EntityType type = numeredEntFiniteElementPtr->getEntType();
401  constexpr std::array<UserDataOperator::OpType, 2> types{
403  std::array<std::string, 2> last_eval_field_name{std::string(), std::string()};
404 
405  auto oit = opPtrVector.begin();
406  auto hi_oit = opPtrVector.end();
407 
408  for (; oit != hi_oit; oit++) {
409 
410  oit->setPtrFE(this);
411 
412  if (oit->opType == UserDataOperator::OPLAST) {
413 
414  // Set field
415  switch (oit->sPace) {
416  case NOSPACE:
417  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Unknown space");
418  case NOFIELD:
419  case H1:
420  case HCURL:
421  case HDIV:
422  case L2:
423  break;
424  default:
425  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
426  "Not implemented for this space", oit->sPace);
427  }
428 
429  // Reseat all data which all field dependent
430  dataOnMaster[oit->sPace]->resetFieldDependentData();
431  dataOnSlave[oit->sPace]->resetFieldDependentData();
432 
433  // Run operator
434  try {
435  CHKERR oit->opRhs(*dataOnMaster[oit->sPace], false);
436  CHKERR oit->opRhs(*dataOnSlave[oit->sPace], false);
437  }
438  CATCH_OP_ERRORS(*oit);
439 
440  } else {
441  boost::shared_ptr<DataForcesAndSourcesCore> op_master_data[2];
442  boost::shared_ptr<DataForcesAndSourcesCore> op_slave_data[2];
443 
444  for (int ss = 0; ss != 2; ss++) {
445 
446  const std::string field_name =
447  !ss ? oit->rowFieldName : oit->colFieldName;
448  const Field *field_struture = mField.get_field_structure(field_name);
449  const BitFieldId data_id = field_struture->getId();
450  const FieldSpace space = field_struture->getSpace();
451  op_master_data[ss] =
452  !ss ? dataOnMaster[space] : derivedDataOnMaster[space];
453  op_slave_data[ss] =
454  !ss ? dataOnSlave[space] : derivedDataOnSlave[space];
455 
456  if ((oit->getNumeredEntFiniteElementPtr()->getBitFieldIdData() &
457  data_id)
458  .none())
459  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
460  "no data field < %s > on finite element < %s >",
461  field_name.c_str(), feName.c_str());
462 
463  if (oit->getOpType() & types[ss] ||
464  oit->getOpType() & UserDataOperator::OPROWCOL) {
465 
466  switch (space) {
467  case NOSPACE:
468  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown space");
469  break;
470  case NOFIELD:
471  case H1:
472  case HCURL:
473  case HDIV:
474  case L2:
475  break;
476  default:
477  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
478  "Not implemented for this space", space);
479  }
480 
481  if (last_eval_field_name[ss] != field_name) {
482  CHKERR getEntityFieldData(*op_master_data[ss], *op_slave_data[ss],
483  field_name, MBEDGE);
484 
485  if (!ss)
487  *op_master_data[ss], *op_slave_data[ss], field_name,
488  const_cast<FENumeredDofEntity_multiIndex &>(
489  numeredEntFiniteElementPtr->getRowsDofs()),
490  MBEDGE);
491  else
493  *op_master_data[ss], *op_slave_data[ss], field_name,
494  const_cast<FENumeredDofEntity_multiIndex &>(
495  numeredEntFiniteElementPtr->getColsDofs()),
496  MBEDGE);
497 
498  switch (space) {
499  case NOSPACE:
500  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
501  "unknown space");
502  break;
503  case H1: {
504 
505  auto get_indices = [&](auto &master, auto &slave, auto &dofs) {
506  return getNodesIndices(
507  field_name, dofs,
508  master.dataOnEntities[MBVERTEX][0].getIndices(),
509  master.dataOnEntities[MBVERTEX][0].getLocalIndices(),
510  slave.dataOnEntities[MBVERTEX][0].getIndices(),
511  slave.dataOnEntities[MBVERTEX][0].getLocalIndices());
512  };
513 
514  auto get_data = [&](DataForcesAndSourcesCore &master_data,
515  DataForcesAndSourcesCore &slave_data) {
516  return getNodesFieldData(
517  field_name,
518  const_cast<FEDofEntity_multiIndex &>(
519  numeredEntFiniteElementPtr->getDataDofs()),
520  master_data.dataOnEntities[MBVERTEX][0].getFieldData(),
521  slave_data.dataOnEntities[MBVERTEX][0].getFieldData(),
522  master_data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
523  slave_data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
524  master_data.dataOnEntities[MBVERTEX][0].getSpace(),
525  slave_data.dataOnEntities[MBVERTEX][0].getSpace(),
526  master_data.dataOnEntities[MBVERTEX][0].getBase(),
527  slave_data.dataOnEntities[MBVERTEX][0].getBase());
528  };
529 
530  if (!ss)
531  CHKERR get_indices(
532  *op_master_data[ss], *op_slave_data[ss],
533  const_cast<FENumeredDofEntity_multiIndex &>(
534  numeredEntFiniteElementPtr->getRowsDofs()));
535  else
536  CHKERR get_indices(
537  *op_master_data[ss], *op_slave_data[ss],
538  const_cast<FENumeredDofEntity_multiIndex &>(
539  numeredEntFiniteElementPtr->getColsDofs()));
540 
541  CHKERR get_data(*op_master_data[ss], *op_slave_data[ss]);
542 
543  } break;
544  case HCURL:
545  case HDIV:
546  case L2:
547  break;
548  case NOFIELD:
549  if (!getNinTheLoop()) {
550  // NOFIELD data are the same for each element, can be
551  // retrieved only once
552  if (!ss) {
553  CHKERR getNoFieldRowIndices(*op_master_data[ss], field_name);
554  CHKERR getNoFieldRowIndices(*op_slave_data[ss], field_name);
555  } else {
556  CHKERR getNoFieldColIndices(*op_master_data[ss], field_name);
557  CHKERR getNoFieldColIndices(*op_slave_data[ss], field_name);
558  }
559  CHKERR getNoFieldFieldData(*op_master_data[ss], field_name);
560  CHKERR getNoFieldFieldData(*op_slave_data[ss], field_name);
561  }
562  break;
563  default:
564  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
565  "Not implemented for this space", space);
566  }
567  last_eval_field_name[ss] = field_name;
568  }
569  }
570  }
571 
572  int type;
573 
574  if (UserDataOperator *cast_oit =
575  dynamic_cast<UserDataOperator *>(&*oit)) {
576  type = cast_oit->getFaceType();
577  if (((oit->getOpType() & UserDataOperator::OPROW) ||
578  (oit->getOpType() & UserDataOperator::OPCOL)) &&
583  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
584  "Wrong combination of FaceType and OpType, OPROW or OPCOL "
585  "combined with face-face OpType");
586  }
587 
588  if ((oit->getOpType() & UserDataOperator::OPROWCOL) &&
591  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
592  "Wrong combination of FaceType and OpType, OPROWCOL "
593  "combined with face-face OpType");
594  }
595 
596  if (!type) {
597  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
598  "Face type is not set");
599  }
600  } else {
606  }
607 
608  if (oit->getOpType() & UserDataOperator::OPROW &&
610  try {
611  CHKERR oit->opRhs(*op_master_data[0], false);
612  }
613  CATCH_OP_ERRORS(*oit);
614  }
615 
616  if (oit->getOpType() & UserDataOperator::OPROW &&
618  try {
619  CHKERR oit->opRhs(*op_slave_data[0], false);
620  }
621  CATCH_OP_ERRORS(*oit);
622  }
623 
624  if (oit->getOpType() & UserDataOperator::OPCOL &&
626  try {
627  CHKERR oit->opRhs(*op_master_data[1], false);
628  }
629  CATCH_OP_ERRORS(*oit);
630  }
631 
632  if (oit->getOpType() & UserDataOperator::OPCOL &&
634  try {
635  CHKERR oit->opRhs(*op_slave_data[1], false);
636  }
637  CATCH_OP_ERRORS(*oit);
638  }
639 
640  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
642  try {
643  CHKERR oit->opLhs(*op_master_data[0], *op_master_data[1]);
644  }
645  CATCH_OP_ERRORS(*oit);
646  }
647 
648  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
650  try {
651  CHKERR oit->opLhs(*op_master_data[0], *op_slave_data[1]);
652  }
653  CATCH_OP_ERRORS(*oit);
654  }
655 
656  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
658  try {
659  CHKERR oit->opLhs(*op_slave_data[0], *op_master_data[1]);
660  }
661  CATCH_OP_ERRORS(*oit);
662  }
663 
664  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
666  try {
667  CHKERR oit->opLhs(*op_slave_data[0], *op_slave_data[1]);
668  }
669  CATCH_OP_ERRORS(*oit);
670  }
671  }
672  }
674 }
boost::ptr_vector< UserDataOperator > opPtrVector
Vector of finite element users data operators.
field with continuous normal traction
Definition: definitions.h:179
MoFEMErrorCode getNoFieldColIndices(DataForcesAndSourcesCore &data, const std::string &field_name) const
get col NoField indices
scalar or vector of scalars describe (no true field)
Definition: definitions.h:176
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
std::string feName
Name of finite element.
MoFEMErrorCode getNodesFieldData(const boost::string_ref field_name, FEDofEntity_multiIndex &dofs, VectorDouble &master_nodes_data, VectorDouble &slave_nodes_data, VectorDofs &master_nodes_dofs, VectorDofs &slave_nodes_dofs, FieldSpace &master_space, FieldSpace &slave_space, FieldApproximationBase &master_base, FieldApproximationBase &slave_base) const
function that gets nodes field data.
int getNinTheLoop() const
get number of evaluated element in the loop
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnMaster
Entity data on element entity columns fields.
MoFEMErrorCode getNoFieldRowIndices(DataForcesAndSourcesCore &data, const std::string &field_name) const
get col NoField indices
MoFEMErrorCode getEntityFieldData(DataForcesAndSourcesCore &master_data, DataForcesAndSourcesCore &slave_data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
function that gets entity field 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.
MoFEMErrorCode getNodesIndices(const boost::string_ref field_name, FENumeredDofEntity_multiIndex &dofs, VectorInt &master_nodes_indices, VectorInt &master_local_nodes_indices, VectorInt &slave_nodes_indices, VectorInt &slave_local_nodes_indices) const
function that gets nodes indices.
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
#define CATCH_OP_ERRORS(OP)
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnSlave
MoFEMErrorCode getEntityIndices(DataForcesAndSourcesCore &master_data, DataForcesAndSourcesCore &slave_data, const std::string &field_name, FENumeredDofEntity_multiIndex &dofs, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
function that gets entity indices.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnSlave
field with continuous tangents
Definition: definitions.h:178
FieldSpace
approximation spaces
Definition: definitions.h:174
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
Definition: Types.hpp:52
#define CHKERR
Inline error check.
Definition: definitions.h:602
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
continuous field
Definition: definitions.h:177
virtual const Field * get_field_structure(const std::string &name)=0
get field structure
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnMaster
Entity data on element entity rows fields.
field with C-1 continuity
Definition: definitions.h:180

◆ operator()()

MoFEMErrorCode MoFEM::ContactPrismElementForcesAndSourcesCore::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::ForcesAndSourcesCore.

Definition at line 159 of file ContactPrismElementForcesAndSourcesCore.cpp.

159  {
161 
162  if (numeredEntFiniteElementPtr->getEntType() != MBPRISM)
165 
166  DataForcesAndSourcesCore &data_div = *dataOnElement[HDIV];
167  DataForcesAndSourcesCore &data_curl = *dataOnElement[HCURL];
168  DataForcesAndSourcesCore &data_l2 = *dataOnElement[HCURL];
169 
171  auto get_coord_and_normal = [&]() {
173  int num_nodes;
174  const EntityHandle *conn;
175  CHKERR mField.get_moab().get_connectivity(ent, conn, num_nodes, true);
176  coords.resize(num_nodes * 3, false);
177  CHKERR mField.get_moab().get_coords(conn, num_nodes,
178  &*coords.data().begin());
179  normal.resize(6, false);
182  aRea[0] = cblas_dnrm2(3, &normal[0], 1) * 0.5;
183  aRea[1] = cblas_dnrm2(3, &normal[3], 1) * 0.5;
185  };
186  CHKERR get_coord_and_normal();
187 
189 
190  // H1
191  if ((dataH1.spacesOnEntities[MBVERTEX]).test(H1)) {
192  CHKERR getEntitySense<MBEDGE>(dataH1);
193  CHKERR getEntityDataOrder<MBEDGE>(dataH1, H1);
194  CHKERR getEntitySense<MBTRI>(dataH1);
195  CHKERR getEntityDataOrder<MBTRI>(dataH1, H1);
196  }
197 
198  // H1
199  if ((dataH1.spacesOnEntities[MBEDGE]).test(HCURL)) {
200  CHKERR getEntitySense<MBEDGE>(data_curl);
201  CHKERR getEntityDataOrder<MBEDGE>(data_curl, HCURL);
202  CHKERR getEntitySense<MBTRI>(data_curl);
203  CHKERR getEntityDataOrder<MBTRI>(data_curl, HCURL);
204  }
205 
206  // Hdiv
207  if ((dataH1.spacesOnEntities[MBTRI]).test(HDIV)) {
208  CHKERR getEntitySense<MBTRI>(data_div);
209  CHKERR getEntityDataOrder<MBTRI>(data_div, HDIV);
210  }
211 
212  // L2
213  if ((dataH1.spacesOnEntities[MBTRI]).test(L2)) {
214  CHKERR getEntitySense<MBTRI>(data_l2);
215  CHKERR getEntityDataOrder<MBTRI>(data_l2, L2);
216  }
217 
218  auto clean_data = [](DataForcesAndSourcesCore &data) {
220  data.sPace.reset();
221  data.bAse.reset();
222  for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
223  data.spacesOnEntities[t].reset();
224  data.basesOnEntities[t].reset();
225  }
226  for (int s = 0; s != LASTSPACE; ++s)
227  data.basesOnSpaces[s].reset();
228 
230  };
231 
232  auto copy_data = [](DataForcesAndSourcesCore &data,
233  DataForcesAndSourcesCore &copy_data, const int shift) {
235 
236  if (shift != 0 && shift != 6) {
237  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
238  "Wrong shift for contact prism element");
239  }
240 
241  data.sPace = copy_data.sPace;
242  data.bAse = copy_data.bAse;
243  data.spacesOnEntities[MBVERTEX] = copy_data.spacesOnEntities[MBVERTEX];
244  data.spacesOnEntities[MBEDGE] = copy_data.spacesOnEntities[MBEDGE];
245  data.spacesOnEntities[MBTRI] = copy_data.spacesOnEntities[MBTRI];
246 
247  data.basesOnEntities[MBVERTEX] = copy_data.basesOnEntities[MBVERTEX];
248  data.basesOnEntities[MBEDGE] = copy_data.basesOnEntities[MBEDGE];
249  data.basesOnEntities[MBTRI] = copy_data.basesOnEntities[MBTRI];
250 
251  data.basesOnSpaces[MBVERTEX] = copy_data.basesOnSpaces[MBVERTEX];
252  data.basesOnSpaces[MBEDGE] = copy_data.basesOnSpaces[MBEDGE];
253  data.basesOnSpaces[MBTRI] = copy_data.basesOnSpaces[MBTRI];
254 
255  for (int ii = 0; ii != 3; ++ii) {
256  data.dataOnEntities[MBEDGE][ii].getSense() =
257  copy_data.dataOnEntities[MBEDGE][ii + shift].getSense();
258  data.dataOnEntities[MBEDGE][ii].getDataOrder() =
259  copy_data.dataOnEntities[MBEDGE][ii + shift].getDataOrder();
260  }
261 
262  if (shift == 0) {
263  data.dataOnEntities[MBTRI][0].getSense() =
264  copy_data.dataOnEntities[MBTRI][3].getSense();
265  data.dataOnEntities[MBTRI][0].getDataOrder() =
266  copy_data.dataOnEntities[MBTRI][3].getDataOrder();
267  } else {
268  data.dataOnEntities[MBTRI][0].getSense() =
269  copy_data.dataOnEntities[MBTRI][4].getSense();
270  data.dataOnEntities[MBTRI][0].getDataOrder() =
271  copy_data.dataOnEntities[MBTRI][4].getDataOrder();
272  }
273 
275  };
276 
277  CHKERR clean_data(dataH1Slave);
278  CHKERR copy_data(dataH1Slave, dataH1, 6);
279  CHKERR clean_data(dataH1Master);
280  CHKERR copy_data(dataH1Master, dataH1, 0);
281 
282  int order_data = getMaxDataOrder();
283  int order_row = getMaxRowOrder();
284  int order_col = getMaxColOrder(); // maybe two different rules?
285  int rule = getRule(order_row, order_col, order_data);
286 
287  if (rule >= 0) {
288 
290 
291  } else {
292 
293  // Master-Slave
294  if (gaussPtsMaster.size2() != gaussPtsSlave.size2())
295  SETERRQ(
296  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
297  "Number of Gauss Points at Master triangle is different than slave");
298 
299  CHKERR setGaussPts(order_row, order_col, order_data);
300  nbGaussPts = gaussPtsMaster.size2();
301  dataH1Master.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nbGaussPts, 3,
302  false);
303  dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nbGaussPts, 3,
304  false);
305 
306  if (nbGaussPts) {
307  CHKERR Tools::shapeFunMBTRI<1>(&*dataH1Master.dataOnEntities[MBVERTEX][0]
308  .getN(NOBASE)
309  .data()
310  .begin(),
311  &gaussPtsMaster(0, 0),
312  &gaussPtsMaster(1, 0), nbGaussPts);
313 
314  CHKERR Tools::shapeFunMBTRI<1>(
315  &*dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).data().begin(),
316  &gaussPtsSlave(0, 0), &gaussPtsSlave(1, 0), nbGaussPts);
317  }
318  }
319 
320  if (nbGaussPts == 0)
322 
323  // Get coordinates on slave and master
324  {
325  coordsAtGaussPtsMaster.resize(nbGaussPts, 3, false);
326  coordsAtGaussPtsSlave.resize(nbGaussPts, 3, false);
327  for (int gg = 0; gg < nbGaussPts; gg++) {
328  for (int dd = 0; dd < 3; dd++) {
330  3, &dataH1Master.dataOnEntities[MBVERTEX][0].getN(NOBASE)(gg, 0), 1,
331  &coords[dd], 3);
333  3, &dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE)(gg, 0), 1,
334  &coords[9 + dd], 3);
335  }
336  }
337  }
338 
339  for (int space = HCURL; space != LASTSPACE; ++space)
340  if (dataOnElement[space]) {
341  dataH1Master.dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
342  dataOnMaster[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
343  dataH1Slave.dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
344  dataOnSlave[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
345  }
346 
347  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
348  if (dataH1.bAse.test(b)) {
349  switch (static_cast<FieldApproximationBase>(b)) {
352  if (dataH1.spacesOnEntities[MBVERTEX].test(H1)) {
353 
354  CHKERR getUserPolynomialBase()->getValue(
356  boost::shared_ptr<BaseFunctionCtx>(new EntPolynomialBaseCtx(
357  dataH1Master, H1, static_cast<FieldApproximationBase>(b),
358  NOBASE)));
359 
360  CHKERR getUserPolynomialBase()->getValue(
362  boost::shared_ptr<BaseFunctionCtx>(new EntPolynomialBaseCtx(
363  dataH1Slave, H1, static_cast<FieldApproximationBase>(b),
364  NOBASE)));
365  }
366 
367  if (dataH1.spacesOnEntities[MBTRI].test(HDIV)) {
368  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
369  "Not yet implemented");
370  }
371  if (dataH1.spacesOnEntities[MBEDGE].test(HCURL)) {
372  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
373  "Not yet implemented");
374  }
375  if (dataH1.spacesOnEntities[MBTET].test(L2)) {
376  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
377  "Not yet implemented");
378  }
379  break;
380  default:
381  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
382  "Not yet implemented");
383  }
384  }
385  }
386 
387  if (dataH1.spacesOnEntities[MBTRI].test(HDIV)) {
388  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented yet");
389  }
390 
391  // Iterate over operators
393 
395 }
int getMaxRowOrder() const
Get max order of approximation for field in rows.
field with continuous normal traction
Definition: definitions.h:179
virtual moab::Interface & get_moab()=0
virtual int getRule(int order_row, int order_col, int order_data)
another variant of getRule
MoFEMErrorCode createDataOnElement()
Create a entity data on element object.
DataForcesAndSourcesCore & dataH1
VectorDouble normal
vector storing vector normal to master or slave element
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
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:483
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:514
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:181
double cblas_dnrm2(const int N, const double *X, const int incX)
Definition: cblas_dnrm2.c:12
double cblas_ddot(const int N, const double *X, const int incX, const double *Y, const int incY)
Definition: cblas_ddot.c:12
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
std::array< double, 2 > aRea
Array storing master and slave faces areas.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnSlave
const Tensor2_symmetric_Expr< const ddTensor0< T, Dim, i, j >, typename promote< T, double >::V, Dim, i, j > dd(const Tensor0< T * > &a, const Index< i, Dim > index1, const Index< j, Dim > index2, const Tensor1< int, Dim > &d_ijk, const Tensor1< double, Dim > &d_xyz)
Definition: ddTensor0.hpp:33
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:152
static MoFEMErrorCode getTriNormal(const double *coords, double *normal)
Get the Tri Normal objectGet triangle normal.
Definition: Tools.cpp:258
field with continuous tangents
Definition: definitions.h:178
virtual MoFEMErrorCode setGaussPts(int order_row, int order_col, int order_data)
set user specific integration rule
#define CHKERR
Inline error check.
Definition: definitions.h:602
std::array< std::bitset< LASTSPACE >, MBMAXTYPE > spacesOnEntities
spaces on entity types
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413
continuous field
Definition: definitions.h:177
MoFEMErrorCode getSpacesAndBaseOnEntities(DataForcesAndSourcesCore &data) const
Get field approximation space and base on entities.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnElement
Entity data on element entity rows fields.
std::bitset< LASTBASE > bAse
bases on element
int getMaxColOrder() const
Get max order of approximation for field in columns.
int getMaxDataOrder() const
Get max order of approximation for data fields.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnMaster
Entity data on element entity rows fields.
field with C-1 continuity
Definition: definitions.h:180

◆ setDefaultGaussPts()

MoFEMErrorCode MoFEM::ContactPrismElementForcesAndSourcesCore::setDefaultGaussPts ( const int  rule)
protected

Definition at line 111 of file ContactPrismElementForcesAndSourcesCore.cpp.

111  {
113 
114  if (rule < QUAD_2D_TABLE_SIZE) {
115  if (QUAD_2D_TABLE[rule]->dim != 2) {
116  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong dimension");
117  }
118  if (QUAD_2D_TABLE[rule]->order < rule) {
119  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
120  "wrong order %d != %d", QUAD_2D_TABLE[rule]->order, rule);
121  }
123  // For master and slave
124  gaussPtsMaster.resize(3, nbGaussPts, false);
125  gaussPtsSlave.resize(3, nbGaussPts, false);
126 
127  cblas_dcopy(nbGaussPts, &QUAD_2D_TABLE[rule]->points[1], 3,
128  &gaussPtsMaster(0, 0), 1);
129  cblas_dcopy(nbGaussPts, &QUAD_2D_TABLE[rule]->points[2], 3,
130  &gaussPtsMaster(1, 0), 1);
131  cblas_dcopy(nbGaussPts, QUAD_2D_TABLE[rule]->weights, 1,
132  &gaussPtsMaster(2, 0), 1);
133 
135 
136  dataH1Master.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nbGaussPts, 3,
137  false);
138 
139  dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nbGaussPts, 3,
140  false);
141 
142  double *shape_ptr_master =
143  &*dataH1Master.dataOnEntities[MBVERTEX][0].getN(NOBASE).data().begin();
144  cblas_dcopy(3 * nbGaussPts, QUAD_2D_TABLE[rule]->points, 1,
145  shape_ptr_master, 1);
146  double *shape_ptr_slave =
147  &*dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).data().begin();
148  cblas_dcopy(3 * nbGaussPts, QUAD_2D_TABLE[rule]->points, 1, shape_ptr_slave,
149  1);
150  } else {
151  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
152  "rule > quadrature order %d < %d", rule, QUAD_2D_TABLE_SIZE);
153  nbGaussPts = 0;
154  }
155 
157 }
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
int npoints
Definition: quad.h:29
void cblas_dcopy(const int N, const double *X, const int incX, double *Y, const int incY)
Definition: cblas_dcopy.c:11
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:483
#define QUAD_2D_TABLE_SIZE
Definition: quad.h:174
static QUAD *const QUAD_2D_TABLE[]
Definition: quad.h:175
const int dim
constexpr int order
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:413

Member Data Documentation

◆ aRea

std::array<double, 2> MoFEM::ContactPrismElementForcesAndSourcesCore::aRea
protected

Array storing master and slave faces areas.

Definition at line 197 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ coords

VectorDouble MoFEM::ContactPrismElementForcesAndSourcesCore::coords
protected

Definition at line 201 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ coordsAtGaussPtsMaster

MatrixDouble MoFEM::ContactPrismElementForcesAndSourcesCore::coordsAtGaussPtsMaster
protected

matrix storing master Gauss points global coordinates

Definition at line 202 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ coordsAtGaussPtsSlave

MatrixDouble MoFEM::ContactPrismElementForcesAndSourcesCore::coordsAtGaussPtsSlave
protected

matrix storing slave Gauss points global coordinates

Definition at line 204 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataH1Master

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataH1Master
protected

Definition at line 235 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataH1Slave

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataH1Slave
protected

Definition at line 236 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataHcurlMaster

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataHcurlMaster
protected

Definition at line 240 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataHcurlSlave

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataHcurlSlave
protected

Definition at line 241 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataHdivMaster

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataHdivMaster
protected

Definition at line 242 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataHdivSlave

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataHdivSlave
protected

Definition at line 243 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataL2Master

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataL2Master
protected

Definition at line 244 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataL2Slave

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataL2Slave
protected

Definition at line 245 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataNoFieldMaster

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataNoFieldMaster
protected

Definition at line 238 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataNoFieldSlave

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataNoFieldSlave
protected

Definition at line 239 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataOnMaster

const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE> MoFEM::ContactPrismElementForcesAndSourcesCore::dataOnMaster
protected

Entity data on element entity rows fields.

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

Definition at line 220 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataOnSlave

const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE> MoFEM::ContactPrismElementForcesAndSourcesCore::dataOnSlave
protected

Definition at line 222 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ derivedDataOnMaster

const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE> MoFEM::ContactPrismElementForcesAndSourcesCore::derivedDataOnMaster
protected

Entity data on element entity columns fields.

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

Definition at line 231 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ derivedDataOnSlave

const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE> MoFEM::ContactPrismElementForcesAndSourcesCore::derivedDataOnSlave
protected

Definition at line 233 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ gaussPtsMaster

MatrixDouble MoFEM::ContactPrismElementForcesAndSourcesCore::gaussPtsMaster
protected

matrix storing master Gauss points local coordinates and weights

Definition at line 207 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ gaussPtsSlave

MatrixDouble MoFEM::ContactPrismElementForcesAndSourcesCore::gaussPtsSlave
protected

matrix storing slave Gauss points local coordinates and weights

Definition at line 209 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ nbGaussPts

int MoFEM::ContactPrismElementForcesAndSourcesCore::nbGaussPts
private

Definition at line 326 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ normal

VectorDouble MoFEM::ContactPrismElementForcesAndSourcesCore::normal
protected

vector storing vector normal to master or slave element

Definition at line 200 of file ContactPrismElementForcesAndSourcesCore.hpp.


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