v0.9.1
Classes | Public Member Functions | Protected Member Functions | Protected 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...
 
MoFEMErrorCode getNumberOfNodes (int &num_nodes) const
 Get number of DOFs on element. 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 ()
 
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 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...
 

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:178
std::array< boost::ptr_vector< EntData >, MBMAXTYPE > dataOnEntities
scalar or vector of scalars describe (no true field)
Definition: definitions.h:175
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:177
continuous field
Definition: definitions.h:176
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnMaster
Entity data on element entity rows fields.
field with C-1 continuity
Definition: definitions.h:179

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 663 of file ContactPrismElementForcesAndSourcesCore.cpp.

666  {
668 
669  auto reset_data = [type_lo, type_hi](auto &data) {
670  for (EntityType t = type_lo; t != type_hi; ++t) {
671  for (auto &dat : data.dataOnEntities[t]) {
672  dat.getDataOrder() = 0;
673  dat.getBase() = NOBASE;
674  dat.getSpace() = NOSPACE;
675  dat.getFieldData().resize(0, false);
676  dat.getFieldDofs().resize(0, false);
677  }
678  }
679  };
680  reset_data(master_data);
681  reset_data(slave_data);
682 
683  auto &dofs = const_cast<FEDofEntity_multiIndex &>(
684  numeredEntFiniteElementPtr->getDataDofs());
685  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
686  auto dit = dofs_by_type.lower_bound(boost::make_tuple(field_name, type_lo));
687  if (dit == dofs_by_type.end())
689  auto hi_dit =
690  dofs_by_type.lower_bound(boost::make_tuple(field_name, type_hi));
691 
692  auto get_data = [&](auto &data, auto &dof, auto type, auto side) {
693  auto &dat = data.dataOnEntities[type][side];
694  auto &ent_field_dofs = dat.getFieldDofs();
695  auto &ent_field_data = dat.getFieldData();
696  if (ent_field_data.empty()) {
697  dat.getBase() = dof.getApproxBase();
698  dat.getSpace() = dof.getSpace();
699  const int ent_order = dof.getMaxOrder();
700  dat.getDataOrder() =
701  dat.getDataOrder() > ent_order ? dat.getDataOrder() : ent_order;
702  const auto dof_ent_field_data = dof.getEntFieldData();
703  const int nb_dofs_on_ent = dof.getNbDofsOnEnt();
704  ent_field_data.resize(nb_dofs_on_ent, false);
705  noalias(ent_field_data) = dof.getEntFieldData();
706  ent_field_dofs.resize(nb_dofs_on_ent, false);
707  for (int ii = 0; ii != nb_dofs_on_ent; ++ii) {
708  ent_field_dofs[ii] = *dit;
709  ++dit;
710  }
711  }
712  };
713 
714  for (; dit != hi_dit;) {
715 
716  auto &dof = **dit;
717  if (dof.getNbDofsOnEnt()) {
718 
719  const EntityType type = dof.getEntType();
720  const int side = dof.sideNumberPtr->side_number;
721 
722  switch (type) {
723  case MBEDGE:
724 
725  if (side < 3)
726  get_data(master_data, dof, type, side);
727  else if (side > 5)
728  get_data(slave_data, dof, type, side - 6);
729 
730  break;
731  case MBTRI:
732 
733  if (side == 3)
734  get_data(master_data, dof, type, 0);
735  if (side == 4)
736  get_data(slave_data, dof, type, 0);
737 
738  break;
739  default:
740  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
741  "Entity type not implemented");
742  };
743 
744  const int brother_side = dof.sideNumberPtr->brother_side_number;
745  if (brother_side != -1)
746  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
747  "Case with brother side not implemented");
748  }
749  }
750 
752 }
#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
#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::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 827 of file ContactPrismElementForcesAndSourcesCore.cpp.

830  {
832 
833  auto clear_data = [type_lo, type_hi](auto &data) {
834  for (EntityType t = type_lo; t != type_hi; ++t) {
835  for (auto &d : data.dataOnEntities[t]) {
836  d.getIndices().resize(0, false);
837  d.getLocalIndices().resize(0, false);
838  }
839  }
840  };
841 
842  clear_data(master_data);
843  clear_data(slave_data);
844 
845  auto &dofs_by_type = dofs.get<Composite_Name_And_Type_mi_tag>();
846  auto dit = dofs_by_type.lower_bound(boost::make_tuple(field_name, type_lo));
847  if (dit == dofs_by_type.end())
849  auto hi_dit =
850  dofs_by_type.lower_bound(boost::make_tuple(field_name, type_hi));
851 
852  auto get_indices = [&](auto &data, auto &dof, const auto type,
853  const auto side) {
854  auto &dat = data.dataOnEntities[type][side];
855  auto &ent_field_indices = dat.getIndices();
856  auto &ent_field_local_indices = dat.getLocalIndices();
857  if (ent_field_indices.empty()) {
858  const int nb_dofs_on_ent = dof.getNbDofsOnEnt();
859  ent_field_indices.resize(nb_dofs_on_ent, false);
860  ent_field_local_indices.resize(nb_dofs_on_ent, false);
861  std::fill(ent_field_indices.data().begin(),
862  ent_field_indices.data().end(), -1);
863  std::fill(ent_field_local_indices.data().begin(),
864  ent_field_local_indices.data().end(), -1);
865  }
866  const int idx = dof.getEntDofIdx();
867  ent_field_indices[idx] = dof.getPetscGlobalDofIdx();
868  ent_field_local_indices[idx] = dof.getPetscLocalDofIdx();
869  };
870 
871  for (; dit != hi_dit; ++dit) {
872 
873  auto &dof = **dit;
874 
875  if (dof.getNbDofsOnEnt()) {
876 
877  const EntityType type = dof.getEntType();
878  const int side = dof.sideNumberPtr->side_number;
879 
880  switch (type) {
881  case MBEDGE:
882 
883  if (side < 3)
884  get_indices(master_data, dof, type, side);
885  else if (side > 5)
886  get_indices(slave_data, dof, type, side - 6);
887 
888  break;
889  case MBTRI:
890 
891  if (side == 3)
892  get_indices(master_data, dof, type, 0);
893  if (side == 4)
894  get_indices(slave_data, dof, type, 0);
895 
896  break;
897  default:
898  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
899  "Entity type not implemented");
900  };
901 
902  const int brother_side = dof.sideNumberPtr->brother_side_number;
903  if (brother_side != -1)
904  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "Not implemented case");
905  }
906  }
907 
909 }
#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
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

◆ 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 754 of file ContactPrismElementForcesAndSourcesCore.cpp.

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

◆ 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 911 of file ContactPrismElementForcesAndSourcesCore.cpp.

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

◆ loopOverOperators()

MoFEMErrorCode MoFEM::ContactPrismElementForcesAndSourcesCore::loopOverOperators ( )
protected

Iterate user data operators.

Returns
MoFEMErrorCode

Definition at line 388 of file ContactPrismElementForcesAndSourcesCore.cpp.

388  {
390 
391  const EntityType type = numeredEntFiniteElementPtr->getEntType();
392  constexpr std::array<UserDataOperator::OpType, 2> types{
394  std::array<std::string, 2> last_eval_field_name{std::string(), std::string()};
395 
396  auto oit = opPtrVector.begin();
397  auto hi_oit = opPtrVector.end();
398 
399  for (; oit != hi_oit; oit++) {
400 
401  oit->setPtrFE(this);
402 
403  if (oit->opType == UserDataOperator::OPLAST) {
404 
405  // Set field
406  switch (oit->sPace) {
407  case NOSPACE:
408  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "Unknown space");
409  case NOFIELD:
410  case H1:
411  case HCURL:
412  case HDIV:
413  case L2:
414  break;
415  default:
416  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
417  "Not implemented for this space", oit->sPace);
418  }
419 
420  // Reseat all data which all field dependent
421  dataOnMaster[oit->sPace]->resetFieldDependentData();
422  dataOnSlave[oit->sPace]->resetFieldDependentData();
423 
424  // Run operator
425  try {
426  CHKERR oit->opRhs(*dataOnMaster[oit->sPace], false);
427  CHKERR oit->opRhs(*dataOnSlave[oit->sPace], false);
428  }
429  CATCH_OP_ERRORS(*oit);
430 
431  } else {
432  boost::shared_ptr<DataForcesAndSourcesCore> op_master_data[2];
433  boost::shared_ptr<DataForcesAndSourcesCore> op_slave_data[2];
434 
435  for (int ss = 0; ss != 2; ss++) {
436 
437  const std::string field_name =
438  !ss ? oit->rowFieldName : oit->colFieldName;
439  const Field *field_struture = mField.get_field_structure(field_name);
440  const BitFieldId data_id = field_struture->getId();
441  const FieldSpace space = field_struture->getSpace();
442  op_master_data[ss] =
443  !ss ? dataOnMaster[space] : derivedDataOnMaster[space];
444  op_slave_data[ss] =
445  !ss ? dataOnSlave[space] : derivedDataOnSlave[space];
446 
447  if ((oit->getNumeredEntFiniteElementPtr()->getBitFieldIdData() &
448  data_id)
449  .none())
450  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
451  "no data field < %s > on finite element < %s >",
452  field_name.c_str(), feName.c_str());
453 
454  if (oit->getOpType() & types[ss] ||
455  oit->getOpType() & UserDataOperator::OPROWCOL) {
456 
457  if (UserDataOperator *cast_oit =
458  dynamic_cast<UserDataOperator *>(&*oit)) {
459  } else {
460  printf("Check\n");
461  }
462  switch (space) {
463  case NOSPACE:
464  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "unknown space");
465  break;
466  case NOFIELD:
467  case H1:
468  case HCURL:
469  case HDIV:
470  case L2:
471  break;
472  default:
473  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
474  "Not implemented for this space", space);
475  }
476 
477  if (last_eval_field_name[ss] != field_name) {
478  CHKERR getEntityFieldData(*op_master_data[ss], *op_slave_data[ss],
479  field_name, MBEDGE);
480 
481  if (!ss)
483  *op_master_data[ss], *op_slave_data[ss], field_name,
484  const_cast<FENumeredDofEntity_multiIndex &>(
485  numeredEntFiniteElementPtr->getRowsDofs()),
486  MBEDGE);
487  else
489  *op_master_data[ss], *op_slave_data[ss], field_name,
490  const_cast<FENumeredDofEntity_multiIndex &>(
491  numeredEntFiniteElementPtr->getColsDofs()),
492  MBEDGE);
493 
494  switch (space) {
495  case NOSPACE:
496  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
497  "unknown space");
498  break;
499  case H1: {
500 
501  auto get_indices = [&](auto &master, auto &slave, auto &dofs) {
502  return getNodesIndices(
503  field_name, dofs,
504  master.dataOnEntities[MBVERTEX][0].getIndices(),
505  master.dataOnEntities[MBVERTEX][0].getLocalIndices(),
506  slave.dataOnEntities[MBVERTEX][0].getIndices(),
507  slave.dataOnEntities[MBVERTEX][0].getLocalIndices());
508  };
509 
510  auto get_data = [&](DataForcesAndSourcesCore &master_data,
511  DataForcesAndSourcesCore &slave_data) {
512  return getNodesFieldData(
513  field_name,
514  const_cast<FEDofEntity_multiIndex &>(
515  numeredEntFiniteElementPtr->getDataDofs()),
516  master_data.dataOnEntities[MBVERTEX][0].getFieldData(),
517  slave_data.dataOnEntities[MBVERTEX][0].getFieldData(),
518  master_data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
519  slave_data.dataOnEntities[MBVERTEX][0].getFieldDofs(),
520  master_data.dataOnEntities[MBVERTEX][0].getSpace(),
521  slave_data.dataOnEntities[MBVERTEX][0].getSpace(),
522  master_data.dataOnEntities[MBVERTEX][0].getBase(),
523  slave_data.dataOnEntities[MBVERTEX][0].getBase());
524  };
525 
526  if (!ss)
527  CHKERR get_indices(
528  *op_master_data[ss], *op_slave_data[ss],
529  const_cast<FENumeredDofEntity_multiIndex &>(
530  numeredEntFiniteElementPtr->getRowsDofs()));
531  else
532  CHKERR get_indices(
533  *op_master_data[ss], *op_slave_data[ss],
534  const_cast<FENumeredDofEntity_multiIndex &>(
535  numeredEntFiniteElementPtr->getColsDofs()));
536 
537  CHKERR get_data(*op_master_data[ss], *op_slave_data[ss]);
538 
539  } break;
540  case HCURL:
541  case HDIV:
542  case L2:
543  break;
544  case NOFIELD:
545  if (!getNinTheLoop()) {
546  // NOFIELD data are the same for each element, can be
547  // retrieved only once
548  if (!ss) {
549  CHKERR getNoFieldRowIndices(*op_master_data[ss], field_name);
550  CHKERR getNoFieldRowIndices(*op_slave_data[ss], field_name);
551  } else {
552  CHKERR getNoFieldColIndices(*op_master_data[ss], field_name);
553  CHKERR getNoFieldColIndices(*op_slave_data[ss], field_name);
554  }
555  CHKERR getNoFieldFieldData(*op_master_data[ss], field_name);
556  CHKERR getNoFieldFieldData(*op_slave_data[ss], field_name);
557  }
558  break;
559  default:
560  SETERRQ1(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED,
561  "Not implemented for this space", space);
562  }
563  last_eval_field_name[ss] = field_name;
564  }
565  }
566  }
567 
568  int type;
569 
570  if (UserDataOperator *cast_oit =
571  dynamic_cast<UserDataOperator *>(&*oit)) {
572  type = cast_oit->getFaceType();
573  if (((oit->getOpType() & UserDataOperator::OPROW) ||
574  (oit->getOpType() & UserDataOperator::OPCOL)) &&
579  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
580  "Wrong combination of FaceType and OpType, OPROW or OPCOL "
581  "combined with face-face OpType");
582  }
583  if (!type) {
584  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
585  "Face type is not set");
586  }
587  } else {
593  }
594 
595  if (oit->getOpType() & UserDataOperator::OPROW &&
597  try {
598  CHKERR oit->opRhs(*op_master_data[0], false);
599  }
600  CATCH_OP_ERRORS(*oit);
601  }
602 
603  if (oit->getOpType() & UserDataOperator::OPROW &&
605  try {
606  CHKERR oit->opRhs(*op_slave_data[0], false);
607  }
608  CATCH_OP_ERRORS(*oit);
609  }
610 
611  if (oit->getOpType() & UserDataOperator::OPCOL &&
613  try {
614  CHKERR oit->opRhs(*op_master_data[1], false);
615  }
616  CATCH_OP_ERRORS(*oit);
617  }
618 
619  if (oit->getOpType() & UserDataOperator::OPCOL &&
621  try {
622  CHKERR oit->opRhs(*op_slave_data[1], false);
623  }
624  CATCH_OP_ERRORS(*oit);
625  }
626 
627  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
629  try {
630  CHKERR oit->opLhs(*op_master_data[0], *op_master_data[1]);
631  }
632  CATCH_OP_ERRORS(*oit);
633  }
634 
635  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
637  try {
638  CHKERR oit->opLhs(*op_master_data[0], *op_slave_data[1]);
639  }
640  CATCH_OP_ERRORS(*oit);
641  }
642 
643  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
645  try {
646  CHKERR oit->opLhs(*op_slave_data[0], *op_master_data[1]);
647  }
648  CATCH_OP_ERRORS(*oit);
649  }
650 
651  if (oit->getOpType() & UserDataOperator::OPROWCOL &&
653  try {
654  CHKERR oit->opLhs(*op_slave_data[0], *op_slave_data[1]);
655  }
656  CATCH_OP_ERRORS(*oit);
657  }
658  }
659  }
661 }
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
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.
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:177
FieldSpace
approximation spaces
Definition: definitions.h:173
std::bitset< BITFIELDID_SIZE > BitFieldId
Field Id.
Definition: Types.hpp:53
#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
continuous field
Definition: definitions.h:176
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:179

◆ 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 110 of file ContactPrismElementForcesAndSourcesCore.cpp.

110  {
112 
113  if (numeredEntFiniteElementPtr->getEntType() != MBPRISM)
116 
117  DataForcesAndSourcesCore &data_div = *dataOnElement[HDIV];
118  DataForcesAndSourcesCore &data_curl = *dataOnElement[HCURL];
119  DataForcesAndSourcesCore &data_l2 = *dataOnElement[HCURL];
120 
122  auto get_coord_and_normal = [&]() {
124  int num_nodes;
125  const EntityHandle *conn;
126  CHKERR mField.get_moab().get_connectivity(ent, conn, num_nodes, true);
127  coords.resize(num_nodes * 3, false);
128  CHKERR mField.get_moab().get_coords(conn, num_nodes,
129  &*coords.data().begin());
130  normal.resize(6, false);
133  aRea[0] = cblas_dnrm2(3, &normal[0], 1) * 0.5;
134  aRea[1] = cblas_dnrm2(3, &normal[3], 1) * 0.5;
136  };
137  CHKERR get_coord_and_normal();
138 
140 
141  // H1
142  if ((dataH1.spacesOnEntities[MBVERTEX]).test(H1)) {
143  CHKERR getEntitySense<MBEDGE>(dataH1);
144  CHKERR getEntityDataOrder<MBEDGE>(dataH1, H1);
145  CHKERR getEntitySense<MBTRI>(dataH1);
146  CHKERR getEntityDataOrder<MBTRI>(dataH1, H1);
147  }
148 
149  // H1
150  if ((dataH1.spacesOnEntities[MBEDGE]).test(HCURL)) {
151  CHKERR getEntitySense<MBEDGE>(data_curl);
152  CHKERR getEntityDataOrder<MBEDGE>(data_curl, HCURL);
153  CHKERR getEntitySense<MBTRI>(data_curl);
154  CHKERR getEntityDataOrder<MBTRI>(data_curl, HCURL);
155  }
156 
157  // Hdiv
158  if ((dataH1.spacesOnEntities[MBTRI]).test(HDIV)) {
159  CHKERR getEntitySense<MBTRI>(data_div);
160  CHKERR getEntityDataOrder<MBTRI>(data_div, HDIV);
161  }
162 
163  // L2
164  if ((dataH1.spacesOnEntities[MBTRI]).test(L2)) {
165  CHKERR getEntitySense<MBTRI>(data_l2);
166  CHKERR getEntityDataOrder<MBTRI>(data_l2, L2);
167  }
168 
169  auto clean_data = [](DataForcesAndSourcesCore &data) {
171  data.sPace.reset();
172  data.bAse.reset();
173  for (EntityType t = MBVERTEX; t != MBMAXTYPE; ++t) {
174  data.spacesOnEntities[t].reset();
175  data.basesOnEntities[t].reset();
176  }
177  for (int s = 0; s != LASTSPACE; ++s)
178  data.basesOnSpaces[s].reset();
179 
181  };
182 
183  auto copy_data = [](DataForcesAndSourcesCore &data,
184  DataForcesAndSourcesCore &copy_data, const int shift) {
186 
187  if (shift != 0 && shift != 6) {
188  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
189  "Wrong shift for contact prism element");
190  }
191 
192  data.sPace = copy_data.sPace;
193  data.bAse = copy_data.bAse;
194  data.spacesOnEntities[MBVERTEX] = copy_data.spacesOnEntities[MBVERTEX];
195  data.spacesOnEntities[MBEDGE] = copy_data.spacesOnEntities[MBEDGE];
196  data.spacesOnEntities[MBTRI] = copy_data.spacesOnEntities[MBTRI];
197 
198  data.basesOnEntities[MBVERTEX] = copy_data.basesOnEntities[MBVERTEX];
199  data.basesOnEntities[MBEDGE] = copy_data.basesOnEntities[MBEDGE];
200  data.basesOnEntities[MBTRI] = copy_data.basesOnEntities[MBTRI];
201 
202  data.basesOnSpaces[MBVERTEX] = copy_data.basesOnSpaces[MBVERTEX];
203  data.basesOnSpaces[MBEDGE] = copy_data.basesOnSpaces[MBEDGE];
204  data.basesOnSpaces[MBTRI] = copy_data.basesOnSpaces[MBTRI];
205 
206  for (int ii = 0; ii != 3; ++ii) {
207  data.dataOnEntities[MBEDGE][ii].getSense() =
208  copy_data.dataOnEntities[MBEDGE][ii + shift].getSense();
209  data.dataOnEntities[MBEDGE][ii].getDataOrder() =
210  copy_data.dataOnEntities[MBEDGE][ii + shift].getDataOrder();
211  }
212 
213  if (shift == 0) {
214  data.dataOnEntities[MBTRI][0].getSense() =
215  copy_data.dataOnEntities[MBTRI][3].getSense();
216  data.dataOnEntities[MBTRI][0].getDataOrder() =
217  copy_data.dataOnEntities[MBTRI][3].getDataOrder();
218  } else {
219  data.dataOnEntities[MBTRI][0].getSense() =
220  copy_data.dataOnEntities[MBTRI][4].getSense();
221  data.dataOnEntities[MBTRI][0].getDataOrder() =
222  copy_data.dataOnEntities[MBTRI][4].getDataOrder();
223  }
224 
226  };
227 
228  CHKERR clean_data(dataH1Slave);
229  CHKERR copy_data(dataH1Slave, dataH1, 6);
230  CHKERR clean_data(dataH1Master);
231  CHKERR copy_data(dataH1Master, dataH1, 0);
232 
233  int order_data = getMaxDataOrder();
234  int order_row = getMaxRowOrder();
235  int order_col = getMaxColOrder(); // maybe two different rules?
236  int rule = getRule(order_row, order_col, order_data);
237 
238  int nb_gauss_pts;
239  if (rule >= 0) {
240  if (rule < QUAD_2D_TABLE_SIZE) {
241  if (QUAD_2D_TABLE[rule]->dim != 2) {
242  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY, "wrong dimension");
243  }
244  if (QUAD_2D_TABLE[rule]->order < rule) {
245  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
246  "wrong order %d != %d", QUAD_2D_TABLE[rule]->order, rule);
247  }
248  nb_gauss_pts = QUAD_2D_TABLE[rule]->npoints;
249  // For master and slave
250  gaussPtsMaster.resize(3, nb_gauss_pts, false);
251  gaussPtsSlave.resize(3, nb_gauss_pts, false);
252 
253  cblas_dcopy(nb_gauss_pts, &QUAD_2D_TABLE[rule]->points[1], 3,
254  &gaussPtsMaster(0, 0), 1);
255  cblas_dcopy(nb_gauss_pts, &QUAD_2D_TABLE[rule]->points[2], 3,
256  &gaussPtsMaster(1, 0), 1);
257  cblas_dcopy(nb_gauss_pts, QUAD_2D_TABLE[rule]->weights, 1,
258  &gaussPtsMaster(2, 0), 1);
259 
261 
262  dataH1Master.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nb_gauss_pts,
263  3, false);
264 
265  dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nb_gauss_pts,
266  3, false);
267 
268  double *shape_ptr_master = &*dataH1Master.dataOnEntities[MBVERTEX][0]
269  .getN(NOBASE)
270  .data()
271  .begin();
272  cblas_dcopy(3 * nb_gauss_pts, QUAD_2D_TABLE[rule]->points, 1,
273  shape_ptr_master, 1);
274  double *shape_ptr_slave =
275  &*dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).data().begin();
276  cblas_dcopy(3 * nb_gauss_pts, QUAD_2D_TABLE[rule]->points, 1,
277  shape_ptr_slave, 1);
278  } else {
279  SETERRQ2(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
280  "rule > quadrature order %d < %d", rule, QUAD_2D_TABLE_SIZE);
281  nb_gauss_pts = 0;
282  }
283 
284  } else {
285 
286  // Master-Slave
287  if (gaussPtsMaster.size2() != gaussPtsSlave.size2())
288  SETERRQ(
289  PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
290  "Number of Gauss Points at Master triangle is different than slave");
291 
292  CHKERR setGaussPts(order_row, order_col, order_data);
293  nb_gauss_pts = gaussPtsMaster.size2();
294  dataH1Master.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nb_gauss_pts,
295  3, false);
296  dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).resize(nb_gauss_pts, 3,
297  false);
298 
299  if (nb_gauss_pts) {
301  .getN(NOBASE)
302  .data()
303  .begin(),
304  &gaussPtsMaster(0, 0), &gaussPtsMaster(1, 0),
305  nb_gauss_pts);
306 
308  &*dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE).data().begin(),
309  &gaussPtsSlave(0, 0), &gaussPtsSlave(1, 0), nb_gauss_pts);
310  }
311  }
312  if (nb_gauss_pts == 0)
314 
315  {
316  coordsAtGaussPtsMaster.resize(nb_gauss_pts, 3, false);
317  coordsAtGaussPtsSlave.resize(nb_gauss_pts, 3, false);
318  for (int gg = 0; gg < nb_gauss_pts; gg++) {
319  for (int dd = 0; dd < 3; dd++) {
321  3, &dataH1Master.dataOnEntities[MBVERTEX][0].getN(NOBASE)(gg, 0), 1,
322  &coords[dd], 3);
324  3, &dataH1Slave.dataOnEntities[MBVERTEX][0].getN(NOBASE)(gg, 0), 1,
325  &coords[9 + dd], 3);
326  }
327  }
328  }
329 
330  for (int space = HCURL; space != LASTSPACE; ++space)
331  if (dataOnElement[space]) {
332  dataH1Master.dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
333  dataOnMaster[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
334  dataH1Slave.dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE) =
335  dataOnSlave[H1]->dataOnEntities[MBVERTEX][0].getNSharedPtr(NOBASE);
336  }
337 
338  for (int b = AINSWORTH_LEGENDRE_BASE; b != LASTBASE; b++) {
339  if (dataH1.bAse.test(b)) {
340  switch (static_cast<FieldApproximationBase>(b)) {
343  if (dataH1.spacesOnEntities[MBVERTEX].test(H1)) {
344 
345  CHKERR getUserPolynomialBase()->getValue(
347  boost::shared_ptr<BaseFunctionCtx>(new EntPolynomialBaseCtx(
348  dataH1Master, H1, static_cast<FieldApproximationBase>(b),
349  NOBASE)));
350 
351  CHKERR getUserPolynomialBase()->getValue(
353  boost::shared_ptr<BaseFunctionCtx>(new EntPolynomialBaseCtx(
354  dataH1Slave, H1, static_cast<FieldApproximationBase>(b),
355  NOBASE)));
356  }
357 
358  if (dataH1.spacesOnEntities[MBTRI].test(HDIV)) {
359  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
360  "Not yet implemented");
361  }
362  if (dataH1.spacesOnEntities[MBEDGE].test(HCURL)) {
363  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
364  "Not yet implemented");
365  }
366  if (dataH1.spacesOnEntities[MBTET].test(L2)) {
367  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
368  "Not yet implemented");
369  }
370  break;
371  default:
372  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
373  "Not yet implemented");
374  }
375  }
376  }
377 
378  if (dataH1.spacesOnEntities[MBTRI].test(HDIV)) {
379  SETERRQ(PETSC_COMM_SELF, MOFEM_NOT_IMPLEMENTED, "not implemented yet");
380  }
381 
382  // Iterate over operators
384 
386 }
int getMaxRowOrder() const
Get max order of approximation for field in rows.
field with continuous normal traction
Definition: definitions.h:178
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
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
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
#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
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
#define QUAD_2D_TABLE_SIZE
Definition: quad.h:174
static QUAD *const QUAD_2D_TABLE[]
Definition: quad.h:175
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
constexpr int order
Ainsworth Cole (Legendre) approx. base .
Definition: definitions.h:151
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:177
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:601
PetscErrorCode ShapeMBTRI(double *N, const double *X, const double *Y, const int G_DIM)
calculate shape functions on triangle
Definition: fem_tools.c:194
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:412
continuous field
Definition: definitions.h:176
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:179

Member Data Documentation

◆ aRea

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

Array storing master and slave faces areas.

Definition at line 167 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ coords

VectorDouble MoFEM::ContactPrismElementForcesAndSourcesCore::coords
protected

Definition at line 171 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ coordsAtGaussPtsMaster

MatrixDouble MoFEM::ContactPrismElementForcesAndSourcesCore::coordsAtGaussPtsMaster
protected

matrix storing master Gauss points global coordinates

Definition at line 172 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ coordsAtGaussPtsSlave

MatrixDouble MoFEM::ContactPrismElementForcesAndSourcesCore::coordsAtGaussPtsSlave
protected

matrix storing slave Gauss points global coordinates

Definition at line 174 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataH1Master

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataH1Master
protected

Definition at line 205 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataH1Slave

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataH1Slave
protected

Definition at line 206 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataHcurlMaster

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataHcurlMaster
protected

Definition at line 210 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataHcurlSlave

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataHcurlSlave
protected

Definition at line 211 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataHdivMaster

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataHdivMaster
protected

Definition at line 212 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataHdivSlave

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataHdivSlave
protected

Definition at line 213 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataL2Master

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataL2Master
protected

Definition at line 214 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataL2Slave

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataL2Slave
protected

Definition at line 215 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataNoFieldMaster

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataNoFieldMaster
protected

Definition at line 208 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataNoFieldSlave

DataForcesAndSourcesCore& MoFEM::ContactPrismElementForcesAndSourcesCore::dataNoFieldSlave
protected

Definition at line 209 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 190 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ dataOnSlave

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

Definition at line 192 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 201 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ derivedDataOnSlave

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

Definition at line 203 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ gaussPtsMaster

MatrixDouble MoFEM::ContactPrismElementForcesAndSourcesCore::gaussPtsMaster
protected

matrix storing master Gauss points local coordinates and weights

Definition at line 177 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ gaussPtsSlave

MatrixDouble MoFEM::ContactPrismElementForcesAndSourcesCore::gaussPtsSlave
protected

matrix storing slave Gauss points local coordinates and weights

Definition at line 179 of file ContactPrismElementForcesAndSourcesCore.hpp.

◆ normal

VectorDouble MoFEM::ContactPrismElementForcesAndSourcesCore::normal
protected

vector storing vector normal to master or slave element

Definition at line 170 of file ContactPrismElementForcesAndSourcesCore.hpp.


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