v0.10.0
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 /* This file is part of MoFEM.
10  * MoFEM is free software: you can redistribute it and/or modify it under
11  * the terms of the GNU Lesser General Public License as published by the
12  * Free Software Foundation, either version 3 of the License, or (at your
13  * option) any later version.
14  *
15  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
16  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
17  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
18  * License for more details.
19  *
20  * You should have received a copy of the GNU Lesser General Public
21  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
22 
23 #ifndef __CONTACTPRISMELEMENTFORCESANDSURCESCORE_HPP__
24 #define __CONTACTPRISMELEMENTFORCESANDSURCESCORE_HPP__
25 
26 namespace MoFEM {
27 template <int SWITCH>
29 
30 /** \brief ContactPrism finite element
31  \ingroup mofem_forces_and_sources_prism_element
32 
33  User is implementing own operator at Gauss points level, by own class
34  derived from ContactPrismElementForcesAndSourcesCoreL::UserDataOperator.
35  Arbitrary number of operator added pushing instances to rowOpPtrVector and
36  rowColOpPtrVector.
37 
38  */
40 
42 
43  /** \brief default operator for Contact Prism element
44  * \ingroup mofem_forces_and_sources_prism_element
45  */
47 
50 
51  UserDataOperator(const std::string &field_name, const char type)
52  : ForcesAndSourcesCore::UserDataOperator(field_name, type) {}
53 
54  UserDataOperator(const std::string &row_field_name,
55  const std::string &col_field_name, const char type)
56  : ForcesAndSourcesCore::UserDataOperator(row_field_name, col_field_name,
57  type) {}
58 
59  UserDataOperator(const std::string &row_field_name,
60  const std::string &col_field_name, const char type,
61  const char face_type)
62  : ForcesAndSourcesCore::UserDataOperator(row_field_name, col_field_name,
63  type),
64  faceType(face_type) {}
65 
66  UserDataOperator(const std::string &field_name, const char type,
67  const char face_type)
69  faceType(face_type) {}
70 
71  /**
72  * \brief Controls loop over faces and face combination on element
73  *
74  * FACEMASTER is used when column or row data needs to be accessed
75  located at master face
76  * FACESLAVE is used when column or row data needs to be accessed
77  located at slave face
78  * FACEMASTERMASTER is used for accessing simultaneously row and col
79  data located at master face.
80  * FACEMASTERSLAVE is used for accessing simultaneously row data that
81  is located on master face and col data located at slave face.
82  * FACESLAVEMASTER is used for accessing simultaneously row data that
83  is located on slave face and col data located at master face.
84  * FACESLAVESLAVE is used for accessing simultaneously row and col
85  data located at slave face.
86  *
87  */
88  enum FaceType {
89  FACEMASTER = 1 << 0,
90  FACESLAVE = 1 << 1,
91  FACEMASTERMASTER = 1 << 2,
92  FACEMASTERSLAVE = 1 << 3,
93  FACESLAVEMASTER = 1 << 4,
94  FACESLAVESLAVE = 1 << 5,
95  FACELAST = 1 << 6
96  };
97 
98  char faceType;
99 
100  /**
101  * \brief Get operator types
102  * @return Return operator type
103  */
104  inline int getFaceType() const;
105 
106  inline boost::shared_ptr<const NumeredEntFiniteElement>
108 
109  /** \brief get face aRea Master
110  */
111  inline double getAreaMaster();
112 
113  /** \brief get face aRea Slave
114  */
115  inline double getAreaSlave();
116 
117  /** \brief get face normal vector to Master face
118  */
120 
121  /** \brief get first face tangent vector to Master face
122  */
124 
125  /** \brief get second face tangent vector to Master face
126  */
128 
129  /** \brief get face normal vector to Slave face
130  */
131  inline VectorAdaptor getNormalSlave();
132 
133  /** \brief get first face tangent vector to Slave face
134  */
136 
137  /** \brief get second face tangent vector to Slave face
138  */
140 
141  /** \brief get Gauss point at Master face
142  */
144 
145  /** \brief get Gauss point at Slave face
146  */
147  inline MatrixDouble &getGaussPtsSlave();
148 
149  /**
150  * @brief Get integration weights for slave side
151  *
152  * \code
153  * auto t_w = getFTensor0IntegrationWeight();
154  * for(int gg = 0; gg!=getGaussPts.size2(); ++gg) {
155  * // integrate something
156  * ++t_w;
157  * }
158  * \endcode
159  *
160  * @return FTensor::Tensor0<FTensor::PackPtr<double *, 1>>
161  */
162  inline auto getFTensor0IntegrationWeightSlave();
163 
164  /**
165  * @brief Get integration weights for master side
166  *
167  * \code
168  * auto t_w = getFTensor0IntegrationWeight();
169  * for(int gg = 0; gg!=getGaussPts.size2(); ++gg) {
170  * // integrate something
171  * ++t_w;
172  * }
173  * \endcode
174  *
175  * @return FTensor::Tensor0<FTensor::PackPtr<double *, 1>>
176  */
178 
179  /** \brief get triangle coordinates
180 
181  Vector has 9 elements, i.e. coordinates on Master face
182 
183  */
184  inline VectorDouble getCoordsMaster();
185 
186  /** \brief get triangle coordinates
187 
188  Vector has 9 elements, i.e. coordinates on Slave face
189 
190  */
191  inline VectorDouble getCoordsSlave();
192 
193  /** \brief get coordinates at Gauss pts on full prism.
194 
195  Matrix has size (nb integration points on master)x(3),
196  i.e. coordinates on face Master
197 
198  */
200 
201  /** \brief get coordinates at Gauss pts on full prism.
202 
203  Matrix has size (nb integration points on slave)x(3),
204  i.e. coordinates on face Slave
205 
206  */
208 
209  /** \brief return pointer to triangle finite element object
210  */
213 
214  /**
215  *
216  * User call this function to loop over elements on the side of face. This
217  * function calls MoFEM::VolumeElementForcesAndSourcesCoreOnContactPrismSide with
218  * is operator to do calculations.
219  *
220  * @param fe_name Name of the element
221  * @param method Finite element object
222  * @param side_type states the side from which side element will work (0
223  * for master 1 for slave)
224  * @return error code
225  */
226  template <int SWITCH>
228  const string &fe_name,
230  const int side_type, const EntityHandle ent_for_side);
231 
232  protected:
233  inline ForcesAndSourcesCore *getSidePtrFE() const;
234  };
235 
237 
238  inline const std::array<boost::shared_ptr<DataForcesAndSourcesCore>,
239  LASTSPACE>
241  return dataOnMaster;
242  }
243 
244  inline const std::array<boost::shared_ptr<DataForcesAndSourcesCore>,
245  LASTSPACE>
247  return dataOnSlave;
248  }
249 
251 
253 
254 protected:
255  std::array<double, 2> aRea; ///< Array storing master and slave faces areas
256 
258  normal; ///< vector storing vector normal to master or slave element
260  MatrixDouble coordsAtGaussPtsMaster; ///< matrix storing master Gauss points
261  ///< global coordinates
262  MatrixDouble coordsAtGaussPtsSlave; ///< matrix storing slave Gauss points
263  ///< global coordinates
264 
265  MatrixDouble gaussPtsMaster; ///< matrix storing master Gauss points local
266  ///< coordinates and weights
267  MatrixDouble gaussPtsSlave; ///< matrix storing slave Gauss points local
268  ///< coordinates and weights
269 
272 
276 
278 
279  /**
280  * @brief Entity data on element entity rows fields
281  *
282  *
283  * FIXME: that should be moved to private class data and acessed only by
284  * member function
285  */
286  const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE>
288  const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE>
290 
291  /**
292  * @brief Entity data on element entity columns fields
293  *
294  * FIXME: that should be moved to private class data and acessed only by
295  * member function
296  */
297  const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE>
299  const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE>
301 
304 
313 
314  MoFEMErrorCode setDefaultGaussPts(const int rule);
315 
316  /**
317  * @brief Iterate user data operators
318  *
319  * @return MoFEMErrorCode
320  */
322 
323  /**
324  * @brief Iterate user data operators
325  *
326  * @return MoFEMErrorCode
327  */
329  FieldApproximationBase m_s_base,
330  DataForcesAndSourcesCore &m_s_data);
331 
332  /** \brief function that gets entity field data.
333  *
334  * \param master_data data fot master face
335  * \param slave_data data fot master face
336  * \param field_name field name of interest
337  * \param type_lo lowest dimension entity type to be searched
338  * \param type_hi highest dimension entity type to be searched
339  */
342  DataForcesAndSourcesCore &slave_data,
343  const std::string &field_name,
344  const EntityType type_lo = MBVERTEX,
345  const EntityType type_hi = MBPOLYHEDRON) const;
346 
347  /** \brief function that gets entity indices.
348  *
349  * \param master_data data fot master face
350  * \param slave_data data fot master face
351  * \param field_name field name of interest
352  * \param dofs MultiIndex container keeping FENumeredDofEntity.
353  * \param type_lo lowest dimension entity type to be searched
354  * \param type_hi highest dimension entity type to be searched
355  */
356  template <typename EXTRACTOR>
358  DataForcesAndSourcesCore &slave_data,
359  const std::string &field_name,
360  FieldEntity_vector_view &ents_field,
361  const EntityType type_lo,
362  const EntityType type_hi,
363  EXTRACTOR &&extractor) const;
364 
365  /** \brief function that gets nodes indices.
366  *
367  * \param field_name field name of interest
368  * \param dofs MultiIndex container keeping FENumeredDofEntity.
369  * \param master_nodes_indices vector containing global master nodes indices
370  * \param master_local_nodes_indices vector containing local master nodes
371  * indices
372  * \param slave_nodes_indices vector containing global master nodes indices
373  * \param slave_local_nodes_indices vector containing local master nodes
374  * indices
375  */
376  template <typename EXTRACTOR>
377  MoFEMErrorCode getNodesIndices(const std::string field_name,
378  FieldEntity_vector_view &ents_field,
379  VectorInt &master_nodes_indices,
380  VectorInt &master_local_nodes_indices,
381  VectorInt &slave_nodes_indices,
382  VectorInt &slave_local_nodes_indices,
383  EXTRACTOR &&extractor) const;
384 
385  /** \brief function that gets nodes field data.
386  *
387  * \param field_name field name of interest
388  * \param dofs MultiIndex container keeping FENumeredDofEntity.
389  * \param master_nodes_data vector containing master nodes data
390  * \param slave_nodes_data vector containing master nodes data
391  * \param master_nodes_dofs vector containing master nodes dofs
392  * \param slave_nodes_dofs vector containing slave nodes dofs
393  * \param master_space approximation energy space at master
394  * \param slave_space approximation energy space at slave
395  * \param master_base base for master face
396  * \param slave_base base for slave face
397  */
399  const std::string field_name, VectorDouble &master_nodes_data,
400  VectorDouble &slave_nodes_data, VectorDofs &master_nodes_dofs,
401  VectorDofs &slave_nodes_dofs,
402 
403  VectorFieldEntities &master_field_entities, VectorFieldEntities &slave_field_entities,
404 
405 
406  FieldSpace &master_space,
407  FieldSpace &slave_space, FieldApproximationBase &master_base,
408  FieldApproximationBase &slave_base) const;
409 
410 private:
412 
413 };
414 
415 boost::shared_ptr<const NumeredEntFiniteElement>
418  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
419  ->numeredEntFiniteElementPtr;
420 };
421 
422 template <int SWITCH>
425  const string &fe_name,
427  const int side_type, const EntityHandle ent_for_side) {
428  return loopSide(fe_name, &fe_method, side_type, ent_for_side);
429 }
430 
431 inline int
433  return faceType;
434 }
435 
436 inline double
438  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[0];
439 }
440 
441 inline double
443  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[1];
444 }
445 
446 inline VectorAdaptor
448  double *data = &(
449  static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->normal[0]);
450  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
451 }
452 
455  double *data = &(static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
456  ->tangentMasterOne[0]);
457  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
458 }
459 
462  double *data = &(static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
463  ->tangentMasterTwo[0]);
464  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
465 }
466 
467 inline VectorAdaptor
469  double *data = &(
470  static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->normal[3]);
471  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
472 }
473 
476  double *data = &(static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
477  ->tangentSlaveOne[0]);
478  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
479 }
480 
483  double *data = &(static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
484  ->tangentSlaveTwo[0]);
485  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
486 }
487 
488 inline MatrixDouble &
490  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
491  ->gaussPtsMaster;
492 }
493 
494 inline MatrixDouble &
496  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
497  ->gaussPtsSlave;
498 }
499 
503  &(getGaussPtsSlave()(getGaussPtsSlave().size1() - 1, 0)));
504 }
505 
509  &(getGaussPtsMaster()(getGaussPtsMaster().size1() - 1, 0)));
510 }
511 
512 inline VectorDouble
514  double *data = &(
515  static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->coords[0]);
516  return VectorAdaptor(9, ublas::shallow_array_adaptor<double>(9, data));
517 }
518 
519 inline VectorDouble
521  double *data = &(
522  static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->coords[9]);
523  return VectorAdaptor(9, ublas::shallow_array_adaptor<double>(9, data));
524 }
525 
528  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
529  ->coordsAtGaussPtsMaster;
530 }
531 
534  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
535  ->coordsAtGaussPtsSlave;
536 }
537 
541  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE);
542 }
543 
544 } // namespace MoFEM
545 
546 #endif //__CONTACTPRISMELEMENTFORCESANDSURCESCORE_HPP__
VectorAdaptor getTangentMasterTwo()
get second face tangent vector to Master face
Deprecated interface functions.
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.
data structure for finite element entityIt keeps that about indices of degrees of freedom,...
VectorDouble normal
vector storing vector normal to master or slave element
std::vector< boost::weak_ptr< FieldEntity > > FieldEntity_vector_view
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:76
MoFEMErrorCode loopSideVolumes(const string &fe_name, VolumeElementForcesAndSourcesCoreOnContactPrismSideSwitch< SWITCH > &fe_method, const int side_type, const EntityHandle ent_for_side)
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnMaster
Entity data on element entity columns fields.
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.
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:21
FieldSpace in [ 0, LASTSPACE )
Definition: definitions.h:181
const ContactPrismElementForcesAndSourcesCore * getContactPrismElementForcesAndSourcesCore()
return pointer to triangle finite element object
UserDataOperator(const std::string &field_name, const char type, const char face_type)
std::array< double, 2 > aRea
Array storing master and slave faces areas.
ublas::vector< FieldEntity *, FieldEntAllocator > VectorFieldEntities
VectorAdaptor getTangentMasterOne()
get first face tangent vector to Master face
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnSlave
MatrixDouble & getCoordsAtGaussPtsMaster()
get coordinates at Gauss pts on full prism.
FieldApproximationBase
approximation base
Definition: definitions.h:150
MoFEMErrorCode getValueHdivDemkowiczBase(MatrixDouble &pts, FieldApproximationBase m_s_base, DataForcesAndSourcesCore &m_s_data)
Iterate user data operators.
ublas::vector< FEDofEntity *, DofsAllocator > VectorDofs
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
transform Hdiv base fluxes from reference element to physical triangle
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:109
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnSlave
boost::shared_ptr< const NumeredEntFiniteElement > getNumeredEntFiniteElementPtr() const
FieldSpace
approximation spaces
Definition: definitions.h:174
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.
VectorAdaptor getTangentSlaveOne()
get first face tangent vector to Slave face
UserDataOperator(const std::string &row_field_name, const std::string &col_field_name, const char type, const char face_type)
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > getDataOnSlaveFromEleSide()
MoFEMErrorCode getEntityIndices(DataForcesAndSourcesCore &master_data, DataForcesAndSourcesCore &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.
MatrixDouble & getCoordsAtGaussPtsSlave()
get coordinates at Gauss pts on full prism.
ublas::vector< int, IntAllocator > VectorInt
Definition: Types.hpp:73
ContactPrism finite elementUser is implementing own operator at Gauss points level,...
MoFEMErrorCode operator()()
function is run for every finite element
UserDataOperator(const std::string &row_field_name, const std::string &col_field_name, const char type)
structure to get information form mofem into DataForcesAndSourcesCore
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:74
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > getDataOnMasterFromEleSide()
VectorAdaptor getTangentSlaveTwo()
get second face tangent vector to Slave face
Data operator to do calculations at integration points.Is inherited and implemented by user to do cal...
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnMaster
Entity data on element entity rows fields.