v0.14.0
Loading...
Searching...
No Matches
ContactPrismElementForcesAndSourcesCore.hpp
Go to the documentation of this file.
1/** \file ContactPrismElementForcesAndSourcesCore.hpp
2 \brief Implementation of the contact prism element.
3
4 These elements are used to enforce contact constraints in the interface
5 between two solids.
6
7*/
8
9
10
11#ifndef __CONTACTPRISMELEMENTFORCESANDSURCESCORE_HPP__
12#define __CONTACTPRISMELEMENTFORCESANDSURCESCORE_HPP__
13
14namespace MoFEM {
15
16struct VolumeElementForcesAndSourcesCoreOnContactPrismSide;
17
18/** \brief ContactPrism finite element
19 \ingroup mofem_forces_and_sources_prism_element
20
21 User is implementing own operator at Gauss points level, by own class
22 derived from ContactPrismElementForcesAndSourcesCoreL::UserDataOperator.
23 Arbitrary number of operator added pushing instances to rowOpPtrVector and
24 rowColOpPtrVector.
25
26 */
28
30
31 /** \brief default operator for Contact Prism element
32 * \ingroup mofem_forces_and_sources_prism_element
33 */
34 struct UserDataOperator;
35
37
38 inline const std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE>
42
43 inline const std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE>
47
49
51
52protected:
53 std::array<double, 2> aRea; ///< Array storing master and slave faces areas
54
56 normal; ///< vector storing vector normal to master or slave element
58 MatrixDouble coordsAtGaussPtsMaster; ///< matrix storing master Gauss points
59 ///< global coordinates
60 MatrixDouble coordsAtGaussPtsSlave; ///< matrix storing slave Gauss points
61 ///< global coordinates
62
63 MatrixDouble gaussPtsMaster; ///< matrix storing master Gauss points local
64 ///< coordinates and weights
65 MatrixDouble gaussPtsSlave; ///< matrix storing slave Gauss points local
66 ///< coordinates and weights
67
71
72 /**
73 * @brief Entity data on element entity rows fields
74 *
75 *
76 * FIXME: that should be moved to private class data and acessed only by
77 * member function
78 */
79 const std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE>
81 const std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE> dataOnSlave;
82
83 /**
84 * @brief Entity data on element entity columns fields
85 *
86 * FIXME: that should be moved to private class data and acessed only by
87 * member function
88 */
89 const std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE>
91 const std::array<boost::shared_ptr<EntitiesFieldData>, LASTSPACE>
93
96
105
106 MoFEMErrorCode setDefaultGaussPts(const int rule);
107
108 /**
109 * @brief Iterate user data operators
110 *
111 * @return MoFEMErrorCode
112 */
114
115 /**
116 * @brief Iterate user data operators
117 *
118 * @return MoFEMErrorCode
119 */
121 FieldApproximationBase m_s_base,
122 EntitiesFieldData &m_s_data);
123
124 /** \brief function that gets entity field data.
125 *
126 * \param master_data data fot master face
127 * \param slave_data data fot master face
128 * \param field_name field name of interest
129 * \param type_lo lowest dimension entity type to be searched
130 * \param type_hi highest dimension entity type to be searched
131 */
133 EntitiesFieldData &master_data, EntitiesFieldData &slave_data,
134 const std::string &field_name, const EntityType type_lo = MBVERTEX,
135 const EntityType type_hi = MBPOLYHEDRON) const;
136
137 /** \brief function that gets entity indices.
138 *
139 * \param master_data data fot master face
140 * \param slave_data data fot master face
141 * \param field_name field name of interest
142 * \param dofs MultiIndex container keeping FENumeredDofEntity.
143 * \param type_lo lowest dimension entity type to be searched
144 * \param type_hi highest dimension entity type to be searched
145 */
146 template <typename EXTRACTOR>
149 EntitiesFieldData &slave_data, const std::string &field_name,
150 FieldEntity_vector_view &ents_field,
151 const EntityType type_lo, const EntityType type_hi,
152 EXTRACTOR &&extractor) const;
153
154 /** \brief function that gets nodes indices.
155 *
156 * \param field_name field name of interest
157 * \param dofs MultiIndex container keeping FENumeredDofEntity.
158 * \param master_nodes_indices vector containing global master nodes indices
159 * \param master_local_nodes_indices vector containing local master nodes
160 * indices
161 * \param slave_nodes_indices vector containing global master nodes indices
162 * \param slave_local_nodes_indices vector containing local master nodes
163 * indices
164 */
165 template <typename EXTRACTOR>
166 MoFEMErrorCode getNodesIndices(const std::string field_name,
167 FieldEntity_vector_view &ents_field,
168 VectorInt &master_nodes_indices,
169 VectorInt &master_local_nodes_indices,
170 VectorInt &slave_nodes_indices,
171 VectorInt &slave_local_nodes_indices,
172 EXTRACTOR &&extractor) const;
173
174 /** \brief function that gets nodes field data.
175 *
176 * \param field_name field name of interest
177 * \param dofs MultiIndex container keeping FENumeredDofEntity.
178 * \param master_nodes_data vector containing master nodes data
179 * \param slave_nodes_data vector containing master nodes data
180 * \param master_nodes_dofs vector containing master nodes dofs
181 * \param slave_nodes_dofs vector containing slave nodes dofs
182 * \param master_space approximation energy space at master
183 * \param slave_space approximation energy space at slave
184 * \param master_base base for master face
185 * \param slave_base base for slave face
186 */
188 VectorDouble &master_nodes_data,
189 VectorDouble &slave_nodes_data,
190 VectorDofs &master_nodes_dofs,
191 VectorDofs &slave_nodes_dofs,
192
193 VectorFieldEntities &master_field_entities,
194 VectorFieldEntities &slave_field_entities,
195
196 FieldSpace &master_space,
197 FieldSpace &slave_space,
198 FieldApproximationBase &master_base,
199 FieldApproximationBase &slave_base) const;
200
201private:
203};
204
205/** \brief default operator for Contact Prism element
206 * \ingroup mofem_forces_and_sources_prism_element
207 */
210
213
214 UserDataOperator(const std::string &field_name, const char type)
216
217 UserDataOperator(const std::string &row_field_name,
218 const std::string &col_field_name, const char type)
219 : ForcesAndSourcesCore::UserDataOperator(row_field_name, col_field_name,
220 type) {}
221
222 UserDataOperator(const std::string &row_field_name,
223 const std::string &col_field_name, const char type,
224 const char face_type)
225 : ForcesAndSourcesCore::UserDataOperator(row_field_name, col_field_name,
226 type),
227 faceType(face_type) {}
228
229 UserDataOperator(const std::string &field_name, const char type,
230 const char face_type)
232 faceType(face_type) {}
233
234 /**
235 * \brief Controls loop over faces and face combination on element
236 *
237 * FACEMASTER is used when column or row data needs to be accessed
238 located at master face
239 * FACESLAVE is used when column or row data needs to be accessed
240 located at slave face
241 * FACEMASTERMASTER is used for accessing simultaneously row and col
242 data located at master face.
243 * FACEMASTERSLAVE is used for accessing simultaneously row data that
244 is located on master face and col data located at slave face.
245 * FACESLAVEMASTER is used for accessing simultaneously row data that
246 is located on slave face and col data located at master face.
247 * FACESLAVESLAVE is used for accessing simultaneously row and col
248 data located at slave face.
249 *
250 */
251 enum FaceType {
252 FACEMASTER = 1 << 0,
253 FACESLAVE = 1 << 1,
258 FACELAST = 1 << 6
259 };
260
262
263 /**
264 * \brief Get operator types
265 * @return Return operator type
266 */
267 inline int getFaceType() const;
268
269 inline boost::shared_ptr<const NumeredEntFiniteElement>
271
272 /** \brief get face aRea Master
273 */
274 inline double getAreaMaster();
275
276 /** \brief get face aRea Slave
277 */
278 inline double getAreaSlave();
279
280 /** \brief get face normal vector to Master face
281 */
283
284 /** \brief get first face tangent vector to Master face
285 */
287
288 /** \brief get second face tangent vector to Master face
289 */
291
292 /** \brief get face normal vector to Slave face
293 */
295
296 /** \brief get first face tangent vector to Slave face
297 */
299
300 /** \brief get second face tangent vector to Slave face
301 */
303
304 /** \brief get Gauss point at Master face
305 */
307
308 /** \brief get Gauss point at Slave face
309 */
311
312 /**
313 * @brief Get integration weights for slave side
314 *
315 * \code
316 * auto t_w = getFTensor0IntegrationWeight();
317 * for(int gg = 0; gg!=getGaussPts.size2(); ++gg) {
318 * // integrate something
319 * ++t_w;
320 * }
321 * \endcode
322 *
323 * @return FTensor::Tensor0<FTensor::PackPtr<double *, 1>>
324 */
326
327 /**
328 * @brief Get integration weights for master side
329 *
330 * \code
331 * auto t_w = getFTensor0IntegrationWeight();
332 * for(int gg = 0; gg!=getGaussPts.size2(); ++gg) {
333 * // integrate something
334 * ++t_w;
335 * }
336 * \endcode
337 *
338 * @return FTensor::Tensor0<FTensor::PackPtr<double *, 1>>
339 */
341
342 /** \brief get triangle coordinates
343
344 Vector has 9 elements, i.e. coordinates on Master face
345
346 */
348
349 /** \brief get triangle coordinates
350
351 Vector has 9 elements, i.e. coordinates on Slave face
352
353 */
355
356 /** \brief get coordinates at Gauss pts on full prism.
357
358 Matrix has size (nb integration points on master)x(3),
359 i.e. coordinates on face Master
360
361 */
363
364 /** \brief get coordinates at Gauss pts on full prism.
365
366 Matrix has size (nb integration points on slave)x(3),
367 i.e. coordinates on face Slave
368
369 */
371
372 /** \brief return pointer to triangle finite element object
373 */
376
377 /**
378 *
379 * User call this function to loop over elements on the side of face. This
380 * function calls MoFEM::VolumeElementForcesAndSourcesCoreOnContactPrismSide
381 * with is operator to do calculations.
382 *
383 * @param fe_name Name of the element
384 * @param method Finite element object
385 * @param side_type states the side from which side element will work (0
386 * for master 1 for slave)
387 * @return error code
388 */
390 const string fe_name,
392 const int side_type, const EntityHandle ent_for_side);
393
394protected:
396};
397
398boost::shared_ptr<const NumeredEntFiniteElement>
399ContactPrismElementForcesAndSourcesCore::UserDataOperator::
400 getNumeredEntFiniteElementPtr() const {
403};
404
405inline int
406ContactPrismElementForcesAndSourcesCore::UserDataOperator::getFaceType() const {
407 return faceType;
408}
409
410inline double
411ContactPrismElementForcesAndSourcesCore::UserDataOperator::getAreaMaster() {
412 return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[0];
413}
414
415inline double
416ContactPrismElementForcesAndSourcesCore::UserDataOperator::getAreaSlave() {
417 return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[1];
418}
419
420inline VectorAdaptor
421ContactPrismElementForcesAndSourcesCore::UserDataOperator::getNormalMaster() {
422 double *data = &(
423 static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->normal[0]);
424 return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
425}
426
427inline VectorAdaptor ContactPrismElementForcesAndSourcesCore::UserDataOperator::
428 getTangentMasterOne() {
429 double *data = &(static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
430 ->tangentMasterOne[0]);
431 return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
432}
433
434inline VectorAdaptor ContactPrismElementForcesAndSourcesCore::UserDataOperator::
435 getTangentMasterTwo() {
436 double *data = &(static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
437 ->tangentMasterTwo[0]);
438 return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
439}
440
441inline VectorAdaptor
442ContactPrismElementForcesAndSourcesCore::UserDataOperator::getNormalSlave() {
443 double *data = &(
444 static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->normal[3]);
445 return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
446}
447
448inline VectorAdaptor ContactPrismElementForcesAndSourcesCore::UserDataOperator::
449 getTangentSlaveOne() {
450 double *data = &(static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
451 ->tangentSlaveOne[0]);
452 return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
453}
454
455inline VectorAdaptor ContactPrismElementForcesAndSourcesCore::UserDataOperator::
456 getTangentSlaveTwo() {
457 double *data = &(static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
458 ->tangentSlaveTwo[0]);
459 return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
460}
461
462inline MatrixDouble &
463ContactPrismElementForcesAndSourcesCore::UserDataOperator::getGaussPtsMaster() {
464 return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
466}
467
468inline MatrixDouble &
469ContactPrismElementForcesAndSourcesCore::UserDataOperator::getGaussPtsSlave() {
470 return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
472}
473
474auto ContactPrismElementForcesAndSourcesCore::UserDataOperator::
475 getFTensor0IntegrationWeightSlave() {
477 &(getGaussPtsSlave()(getGaussPtsSlave().size1() - 1, 0)));
478}
479
480auto ContactPrismElementForcesAndSourcesCore::UserDataOperator::
481 getFTensor0IntegrationWeightMaster() {
483 &(getGaussPtsMaster()(getGaussPtsMaster().size1() - 1, 0)));
484}
485
486inline VectorDouble
487ContactPrismElementForcesAndSourcesCore::UserDataOperator::getCoordsMaster() {
488 double *data = &(
489 static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->coords[0]);
490 return VectorAdaptor(9, ublas::shallow_array_adaptor<double>(9, data));
491}
492
493inline VectorDouble
494ContactPrismElementForcesAndSourcesCore::UserDataOperator::getCoordsSlave() {
495 double *data = &(
496 static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->coords[9]);
497 return VectorAdaptor(9, ublas::shallow_array_adaptor<double>(9, data));
498}
499
500inline MatrixDouble &ContactPrismElementForcesAndSourcesCore::UserDataOperator::
501 getCoordsAtGaussPtsMaster() {
502 return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
504}
505
506inline MatrixDouble &ContactPrismElementForcesAndSourcesCore::UserDataOperator::
507 getCoordsAtGaussPtsSlave() {
508 return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
510}
511
513ContactPrismElementForcesAndSourcesCore::UserDataOperator::
514 getContactPrismElementForcesAndSourcesCore() {
515 return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE);
516}
517
518} // namespace MoFEM
519
520#endif //__CONTACTPRISMELEMENTFORCESANDSURCESCORE_HPP__
FieldApproximationBase
approximation base
Definition definitions.h:58
FieldSpace
approximation spaces
Definition definitions.h:82
@ LASTSPACE
FieldSpace in [ 0, LASTSPACE )
Definition definitions.h:89
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition Types.hpp:115
UBlasVector< int > VectorInt
Definition Types.hpp:67
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
std::vector< boost::weak_ptr< FieldEntity > > FieldEntity_vector_view
ublas::vector< FieldEntity *, FieldEntAllocator > VectorFieldEntities
ublas::vector< FEDofEntity *, DofsAllocator > VectorDofs
constexpr auto field_name
UserDataOperator(const std::string &field_name, const char type, const char face_type)
boost::shared_ptr< const NumeredEntFiniteElement > getNumeredEntFiniteElementPtr() const
UserDataOperator(const std::string &row_field_name, const std::string &col_field_name, const char type, const char face_type)
VectorAdaptor getTangentMasterOne()
get first face tangent vector to Master face
VectorAdaptor getTangentSlaveTwo()
get second face tangent vector to Slave face
const ContactPrismElementForcesAndSourcesCore * getContactPrismElementForcesAndSourcesCore()
return pointer to triangle finite element object
UserDataOperator(const std::string &row_field_name, const std::string &col_field_name, const char type)
MoFEMErrorCode loopSideVolumes(const string fe_name, VolumeElementForcesAndSourcesCoreOnContactPrismSide &fe_method, const int side_type, const EntityHandle ent_for_side)
VectorAdaptor getTangentSlaveOne()
get first face tangent vector to Slave face
VectorAdaptor getTangentMasterTwo()
get second face tangent vector to Master face
MatrixDouble & getCoordsAtGaussPtsMaster()
get coordinates at Gauss pts on full prism.
MatrixDouble & getCoordsAtGaussPtsSlave()
get coordinates at Gauss pts on full prism.
MoFEMErrorCode getValueHdivDemkowiczBase(MatrixDouble &pts, FieldApproximationBase m_s_base, EntitiesFieldData &m_s_data)
Iterate user data operators.
MoFEMErrorCode operator()()
function is run for every finite element
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > getDataOnMasterFromEleSide()
VectorDouble normal
vector storing vector normal to master or slave element
MoFEMErrorCode getEntityIndices(EntitiesFieldData &master_data, EntitiesFieldData &slave_data, const std::string &field_name, FieldEntity_vector_view &ents_field, const EntityType type_lo, const EntityType type_hi, EXTRACTOR &&extractor) const
function that gets entity indices.
MoFEMErrorCode getNodesIndices(const std::string field_name, FieldEntity_vector_view &ents_field, VectorInt &master_nodes_indices, VectorInt &master_local_nodes_indices, VectorInt &slave_nodes_indices, VectorInt &slave_local_nodes_indices, EXTRACTOR &&extractor) const
function that gets nodes indices.
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > derivedDataOnSlave
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > dataOnMaster
Entity data on element entity rows fields.
std::array< double, 2 > aRea
Array storing master and slave faces areas.
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > derivedDataOnMaster
Entity data on element entity columns fields.
MoFEMErrorCode getNodesFieldData(const std::string field_name, VectorDouble &master_nodes_data, VectorDouble &slave_nodes_data, VectorDofs &master_nodes_dofs, VectorDofs &slave_nodes_dofs, VectorFieldEntities &master_field_entities, VectorFieldEntities &slave_field_entities, FieldSpace &master_space, FieldSpace &slave_space, FieldApproximationBase &master_base, FieldApproximationBase &slave_base) const
function that gets nodes field data.
MoFEMErrorCode getEntityFieldData(EntitiesFieldData &master_data, EntitiesFieldData &slave_data, const std::string &field_name, const EntityType type_lo=MBVERTEX, const EntityType type_hi=MBPOLYHEDRON) const
function that gets entity field data.
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > dataOnSlave
const std::array< boost::shared_ptr< EntitiesFieldData >, LASTSPACE > getDataOnSlaveFromEleSide()
Deprecated interface functions.
data structure for finite element entity
boost::shared_ptr< const NumeredEntFiniteElement > numeredEntFiniteElementPtr
structure to get information form mofem into EntitiesFieldData
transform Hdiv base fluxes from reference element to physical triangle
Base volume element used to integrate on contact surface (could be extended to other volume elements)...