v0.9.1
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 
28 /** \brief ContactPrism finite element
29  \ingroup mofem_forces_and_sources_prism_element
30 
31  User is implementing own operator at Gauss points level, by own class
32  derived from ContactPrismElementForcesAndSourcesCoreL::UserDataOperator.
33  Arbitrary number of operator added pushing instances to rowOpPtrVector and
34  rowColOpPtrVector.
35 
36  */
38 
40 
41  /** \brief default operator for Contact Prism element
42  * \ingroup mofem_forces_and_sources_prism_element
43  */
45 
48 
49  UserDataOperator(const std::string &field_name, const char type)
50  : ForcesAndSourcesCore::UserDataOperator(field_name, type) {}
51 
52  UserDataOperator(const std::string &row_field_name,
53  const std::string &col_field_name, const char type)
54  : ForcesAndSourcesCore::UserDataOperator(row_field_name, col_field_name,
55  type) {}
56 
57  UserDataOperator(const std::string &row_field_name,
58  const std::string &col_field_name, const char type,
59  const char face_type)
60  : ForcesAndSourcesCore::UserDataOperator(row_field_name, col_field_name,
61  type),
62  faceType(face_type) {}
63 
64  UserDataOperator(const std::string &field_name, const char type,
65  const char face_type)
67  faceType(face_type) {}
68 
69  /**
70  * \brief Controls loop over faces and face combination on element
71  *
72  * FACEMASTER is used when column or row data needs to be accessed
73  located at master face
74  * FACESLAVE is used when column or row data needs to be accessed
75  located at slave face
76  * FACEMASTERMASTER is used for accessing simultaneously row and col
77  data located at master face.
78  * FACEMASTERSLAVE is used for accessing simultaneously row data that
79  is located on master face and col data located at slave face.
80  * FACESLAVEMASTER is used for accessing simultaneously row data that
81  is located on slave face and col data located at master face.
82  * FACESLAVESLAVE is used for accessing simultaneously row and col
83  data located at slave face.
84  *
85  */
86  enum FaceType {
87  FACEMASTER = 1 << 0,
88  FACESLAVE = 1 << 1,
89  FACEMASTERMASTER = 1 << 2,
90  FACEMASTERSLAVE = 1 << 3,
91  FACESLAVEMASTER = 1 << 4,
92  FACESLAVESLAVE = 1 << 5,
93  FACELAST = 1 << 6
94  };
95 
96  char faceType;
97 
98  /**
99  * \brief Get operator types
100  * @return Return operator type
101  */
102  inline int getFaceType() const;
103 
104  /** \brief get face aRea Master
105  */
106  inline double getAreaMaster();
107 
108  /** \brief get face aRea Slave
109  */
110  inline double getAreaSlave();
111 
112  /** \brief get face normal vector to Master face
113  */
115 
116  /** \brief get face normal vector to Slave face
117  */
118  inline VectorAdaptor getNormalSlave();
119 
120  /** \brief get Gauss point at Master face
121  */
123 
124  /** \brief get Gauss point at Slave face
125  */
126  inline MatrixDouble &getGaussPtsSlave();
127 
128  /** \brief get triangle coordinates
129 
130  Vector has 9 elements, i.e. coordinates on Master face
131 
132  */
133  inline VectorDouble getCoordsMaster();
134 
135  /** \brief get triangle coordinates
136 
137  Vector has 9 elements, i.e. coordinates on Slave face
138 
139  */
140  inline VectorDouble getCoordsSlave();
141 
142  /** \brief get coordinates at Gauss pts on full prism.
143 
144  Matrix has size (nb integration points on master)x(3),
145  i.e. coordinates on face Master
146 
147  */
149 
150  /** \brief get coordinates at Gauss pts on full prism.
151 
152  Matrix has size (nb integration points on slave)x(3),
153  i.e. coordinates on face Slave
154 
155  */
157 
158  /** \brief return pointer to triangle finite element object
159  */
162  };
163 
165 
166 protected:
167  std::array<double, 2> aRea; ///< Array storing master and slave faces areas
168 
170  normal; ///< vector storing vector normal to master or slave element
172  MatrixDouble coordsAtGaussPtsMaster; ///< matrix storing master Gauss points
173  ///< global coordinates
174  MatrixDouble coordsAtGaussPtsSlave; ///< matrix storing slave Gauss points
175  ///< global coordinates
176 
177  MatrixDouble gaussPtsMaster; ///< matrix storing master Gauss points local
178  ///< coordinates and weights
179  MatrixDouble gaussPtsSlave; ///< matrix storing slave Gauss points local
180  ///< coordinates and weights
181 
182  /**
183  * @brief Entity data on element entity rows fields
184  *
185  *
186  * FIXME: that should be moved to private class data and acessed only by
187  * member function
188  */
189  const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE>
191  const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE>
193 
194  /**
195  * @brief Entity data on element entity columns fields
196  *
197  * FIXME: that should be moved to private class data and acessed only by
198  * member function
199  */
200  const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE>
202  const std::array<boost::shared_ptr<DataForcesAndSourcesCore>, LASTSPACE>
204 
207 
216 
217  /**
218  * @brief Iterate user data operators
219  *
220  * @return MoFEMErrorCode
221  */
223 
224  /** \brief function that gets entity field data.
225  *
226  * \param master_data data fot master face
227  * \param slave_data data fot master face
228  * \param field_name field name of interest
229  * \param type_lo lowest dimension entity type to be searched
230  * \param type_hi highest dimension entity type to be searched
231  */
234  DataForcesAndSourcesCore &slave_data,
235  const std::string &field_name,
236  const EntityType type_lo = MBVERTEX,
237  const EntityType type_hi = MBPOLYHEDRON) const;
238 
239  /** \brief function that gets entity indices.
240  *
241  * \param master_data data fot master face
242  * \param slave_data data fot master face
243  * \param field_name field name of interest
244  * \param dofs MultiIndex container keeping FENumeredDofEntity.
245  * \param type_lo lowest dimension entity type to be searched
246  * \param type_hi highest dimension entity type to be searched
247  */
249  DataForcesAndSourcesCore &master_data,
250  DataForcesAndSourcesCore &slave_data, const std::string &field_name,
251  FENumeredDofEntity_multiIndex &dofs, const EntityType type_lo = MBVERTEX,
252  const EntityType type_hi = MBPOLYHEDRON) const;
253 
254  /** \brief function that gets nodes indices.
255  *
256  * \param field_name field name of interest
257  * \param dofs MultiIndex container keeping FENumeredDofEntity.
258  * \param master_nodes_indices vector containing global master nodes indices
259  * \param master_local_nodes_indices vector containing local master nodes
260  * indices
261  * \param slave_nodes_indices vector containing global master nodes indices
262  * \param slave_local_nodes_indices vector containing local master nodes
263  * indices
264  */
265  MoFEMErrorCode getNodesIndices(const boost::string_ref field_name,
267  VectorInt &master_nodes_indices,
268  VectorInt &master_local_nodes_indices,
269  VectorInt &slave_nodes_indices,
270  VectorInt &slave_local_nodes_indices) const;
271 
272  /** \brief function that gets nodes field data.
273  *
274  * \param field_name field name of interest
275  * \param dofs MultiIndex container keeping FENumeredDofEntity.
276  * \param master_nodes_data vector containing master nodes data
277  * \param slave_nodes_data vector containing master nodes data
278  * \param master_nodes_dofs vector containing master nodes dofs
279  * \param slave_nodes_dofs vector containing slave nodes dofs
280  * \param master_space approximation energy space at master
281  * \param slave_space approximation energy space at slave
282  * \param master_base base for master face
283  * \param slave_base base for slave face
284  */
286  const boost::string_ref field_name, FEDofEntity_multiIndex &dofs,
287  VectorDouble &master_nodes_data, VectorDouble &slave_nodes_data,
288  VectorDofs &master_nodes_dofs, VectorDofs &slave_nodes_dofs,
289  FieldSpace &master_space, FieldSpace &slave_space,
290  FieldApproximationBase &master_base,
291  FieldApproximationBase &slave_base) const;
292 };
293 
294 inline int
296  return faceType;
297 }
298 
299 inline double
301  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[0];
302 }
303 
304 inline double
306  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->aRea[1];
307 }
308 
309 inline VectorAdaptor
311  double *data = &(
312  static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->normal[0]);
313  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
314 }
315 
316 inline VectorAdaptor
318  double *data = &(
319  static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->normal[3]);
320  return VectorAdaptor(3, ublas::shallow_array_adaptor<double>(3, data));
321 }
322 
323 inline MatrixDouble &
325  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
326  ->gaussPtsMaster;
327 }
328 
329 inline MatrixDouble &
331  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
332  ->gaussPtsSlave;
333 }
334 
335 inline VectorDouble
337  double *data = &(
338  static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->coords[0]);
339  return VectorAdaptor(9, ublas::shallow_array_adaptor<double>(9, data));
340 }
341 
342 inline VectorDouble
344  double *data = &(
345  static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)->coords[9]);
346  return VectorAdaptor(9, ublas::shallow_array_adaptor<double>(9, data));
347 }
348 
351  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
352  ->coordsAtGaussPtsMaster;
353 }
354 
357  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE)
358  ->coordsAtGaussPtsSlave;
359 }
360 
364  return static_cast<ContactPrismElementForcesAndSourcesCore *>(ptrFE);
365 }
366 
367 } // namespace MoFEM
368 
369 #endif //__CONTACTPRISMELEMENTFORCESANDSURCESCORE_HPP__
Deprecated interface functions.
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
const ContactPrismElementForcesAndSourcesCore * getContactPrismElementForcesAndSourcesCore()
return pointer to triangle finite element object
ublas::matrix< double, ublas::row_major, DoubleAllocator > MatrixDouble
Definition: Types.hpp:74
multi_index_container< boost::shared_ptr< FENumeredDofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< FENumeredDofEntity::interface_type_DofEntity, const UId &, &FENumeredDofEntity::getGlobalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FENumeredDofEntity::interface_type_DofEntity, EntityHandle, &FENumeredDofEntity::getEnt > >, ordered_non_unique< tag< FieldName_mi_tag >, const_mem_fun< FENumeredDofEntity::interface_type_Field, boost::string_ref, &FENumeredDofEntity::getNameRef > >, ordered_non_unique< tag< Composite_Name_Type_And_Side_Number_mi_tag >, composite_key< FENumeredDofEntity, const_mem_fun< FENumeredDofEntity::interface_type_Field, boost::string_ref, &FENumeredDofEntity::getNameRef >, const_mem_fun< FENumeredDofEntity::interface_type_RefEntity, EntityType, &FENumeredDofEntity::getEntType >, KeyFromKey< member< SideNumber, char, &SideNumber::side_number >, member< FENumeredDofEntity::BaseFEEntity, boost::shared_ptr< SideNumber >, &FENumeredDofEntity::sideNumberPtr > > > >, ordered_non_unique< tag< Composite_Name_And_Type_mi_tag >, composite_key< FENumeredDofEntity, const_mem_fun< FENumeredDofEntity::interface_type_Field, boost::string_ref, &FENumeredDofEntity::getNameRef >, const_mem_fun< FENumeredDofEntity::interface_type_RefEntity, EntityType, &FENumeredDofEntity::getEntType > > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< FENumeredDofEntity, const_mem_fun< FENumeredDofEntity::interface_type_Field, boost::string_ref, &FENumeredDofEntity::getNameRef >, const_mem_fun< FENumeredDofEntity::interface_type_DofEntity, EntityHandle, &FENumeredDofEntity::getEnt > > > > > FENumeredDofEntity_multiIndex
MultiIndex container keeps FENumeredDofEntity.
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.
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:180
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.
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.
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > dataOnSlave
FieldApproximationBase
approximation base
Definition: definitions.h:149
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:66
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.
VectorShallowArrayAdaptor< double > VectorAdaptor
Definition: Types.hpp:104
const std::array< boost::shared_ptr< DataForcesAndSourcesCore >, LASTSPACE > derivedDataOnSlave
FieldSpace
approximation spaces
Definition: definitions.h:173
UserDataOperator(const std::string &row_field_name, const std::string &col_field_name, const char type, const char face_type)
multi_index_container< boost::shared_ptr< FEDofEntity >, indexed_by< ordered_unique< tag< Unique_mi_tag >, const_mem_fun< FEDofEntity::interface_type_DofEntity, const UId &, &FEDofEntity::getGlobalUniqueId > >, ordered_non_unique< tag< Ent_mi_tag >, const_mem_fun< FEDofEntity::interface_type_DofEntity, EntityHandle, &FEDofEntity::getEnt > >, ordered_non_unique< tag< FieldName_mi_tag >, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef > >, ordered_non_unique< tag< EntType_mi_tag >, const_mem_fun< FEDofEntity::interface_type_RefEntity, EntityType, &FEDofEntity::getEntType > >, ordered_non_unique< tag< Composite_Name_And_Type_mi_tag >, composite_key< FEDofEntity, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef >, const_mem_fun< FEDofEntity::interface_type_RefEntity, EntityType, &FEDofEntity::getEntType > > >, ordered_non_unique< tag< Composite_Name_And_Ent_mi_tag >, composite_key< FEDofEntity, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef >, const_mem_fun< FEDofEntity::interface_type_DofEntity, EntityHandle, &FEDofEntity::getEnt > > >, ordered_non_unique< tag< Composite_Name_Type_And_Side_Number_mi_tag >, composite_key< FEDofEntity, const_mem_fun< FEDofEntity::interface_type_Field, boost::string_ref, &FEDofEntity::getNameRef >, const_mem_fun< FEDofEntity::interface_type_RefEntity, EntityType, &FEDofEntity::getEntType >, KeyFromKey< member< SideNumber, char, &SideNumber::side_number >, member< FEDofEntity::BaseFEEntity, boost::shared_ptr< SideNumber >, &FEDofEntity::sideNumberPtr > > > > > > FEDofEntity_multiIndex
MultiIndex container keeps FEDofEntity.
ublas::vector< int, IntAllocator > VectorInt
Definition: Types.hpp:71
ublas::vector< boost::shared_ptr< const FEDofEntity >, DofsAllocator > VectorDofs
ContactPrism finite elementUser is implementing own operator at Gauss points level,...
MatrixDouble & getCoordsAtGaussPtsMaster()
get coordinates at Gauss pts on full prism.
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
MatrixDouble & getCoordsAtGaussPtsSlave()
get coordinates at Gauss pts on full prism.
ublas::vector< double, DoubleAllocator > VectorDouble
Definition: Types.hpp:72
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.