v0.14.0
EdgeElementForcesAndSourcesCore.hpp
Go to the documentation of this file.
1 /** \file EdgeElementForcesAndSourcesCore.hpp
2 
3  \brief Implementation of elements on entities.
4 
5  Those element are inherited by user to implement specific implementation of
6  particular problem.
7 
8 */
9 
10 
11 
12 #ifndef __EDGEELEMENTFORCESANDSURCESCORE_HPP__
13 #define __EDGEELEMENTFORCESANDSURCESCORE_HPP__
14 
15 using namespace boost::numeric;
16 
17 namespace MoFEM {
18 
19 struct FaceElementForcesAndSourcesCoreOnSide;
20 
21 /** \brief Edge finite element
22  * \ingroup mofem_forces_and_sources_edge_element
23  *
24  * User is implementing own operator at Gauss points level, by own object
25  * derived from EdgeElementForcesAndSourcesCoreL::UserDataOperator. Arbitrary
26  * number of operator added pushing objects to rowOpPtrVector and
27  * rowColOpPtrVector.
28  *
29  */
32 
34 
35  /** \brief default operator for EDGE element
36  \ingroup mofem_forces_and_sources_edge_element
37  */
38  struct UserDataOperator;
39 
40  enum Switches {};
41 
42  MoFEMErrorCode operator()();
43 
45 
46 protected:
47 
50 
51  double &lEngth;
52 
53  int numNodes;
57 
58  MoFEMErrorCode calculateEdgeDirection();
59  MoFEMErrorCode setIntegrationPts();
60  MoFEMErrorCode calculateCoordsAtIntegrationPts();
61 
63 };
64 
65 /** \brief default operator for EDGE element
66  \ingroup mofem_forces_and_sources_edge_element
67  */
70 
72 
73  /** \brief get element connectivity
74  */
75  inline const EntityHandle *getConn();
76 
77  /**
78  * \brief get edge length
79  */
80  inline double getLength();
81 
82  /**
83  * \brief get edge direction
84  */
85  inline VectorDouble &getDirection();
86 
87  /**
88  * \brief get edge normal
89  * NOTE: it should be used only in 2D analysis
90  */
91  inline auto getFTensor1Normal();
92 
93  /**
94  * @brief get ftensor1 edge normal
95  *
96  * @param vec vector in third direction
97  * @return auto
98  */
99  inline auto getFTensor1Normal(const FTensor::Tensor1<double, 3> &vec);
100 
101  inline auto
102  getFTensor1NormalsAtGaussPts(const FTensor::Tensor1<double, 3> &vec);
103 
104  /** \brief get normal at integration points
105 
106  Example:
107  \code
108  double nrm2;
109  FTensor::Index<'i',3> i;
110  auto t_normal = getFTensor1NormalsAtGaussPts();
111  for(int gg = gg!=data.getN().size1();gg++) {
112  nrm2 = sqrt(t_normal(i)*t_normal(i));
113  ++t_normal;
114  }
115  \endcode
116 
117  */
118  inline auto getFTensor1NormalsAtGaussPts();
119 
120  /**
121  * \brief get edge node coordinates
122  */
123  inline VectorDouble &getCoords();
124 
125  /**
126  * \brief get tangent vector to edge curve at integration points
127  */
128  inline MatrixDouble &getTangentAtGaussPts();
129 
131  return getTangentAtGaussPts();
132  }
133 
134  /**
135  * \brief get pointer to this finite element
136  */
137  inline const EdgeElementForcesAndSourcesCore *getEdgeFE();
138 
139  inline FTensor::Tensor1<double, 3> getFTensor1Direction();
140 
141  /**
142  * \brief get get coords at gauss points
143 
144  \code
145  FTensor::Index<'i',3> i;
146  auto t_center;
147  auto t_coords = getFTensor1Coords();
148  t_center(i) = 0;
149  for(int nn = 0;nn!=2;nn++) {
150  t_center(i) += t_coords(i);
151  ++t_coords;
152  }
153  t_center(i) /= 2;
154  \endcode
155 
156  */
157  inline FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3> getFTensor1Coords();
158 
159  /**
160  * @deprecated use getFTensor1Coords
161  */
164  return getFTensor1Coords();
165  }
166 
167  /**
168  * @brief Get tangent vector at integration points
169  *
170  * @return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, DIM>
171  */
172  template <int DIM = 3>
174  getFTensor1TangentAtGaussPts();
175 
177  loopSideFaces(const string fe_name,
179 
180 protected:
181  MoFEMErrorCode setPtrFE(ForcesAndSourcesCore *ptr);
182 };
183 
184 const EntityHandle *
185 EdgeElementForcesAndSourcesCore::UserDataOperator::getConn() {
186  return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE)->cOnn;
187 }
188 
189 double EdgeElementForcesAndSourcesCore::UserDataOperator::getLength() {
190  return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE)->lEngth;
191 }
192 
193 VectorDouble &
194 EdgeElementForcesAndSourcesCore::UserDataOperator::getDirection() {
195  return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE)->dIrection;
196 }
197 
198 VectorDouble &
199 EdgeElementForcesAndSourcesCore::UserDataOperator::getCoords() {
200  return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE)->cOords;
201 }
202 
203 MatrixDouble &
204 EdgeElementForcesAndSourcesCore::UserDataOperator::getTangentAtGaussPts() {
205  return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE)
206  ->tangentAtGaussPts;
207 }
208 
210 EdgeElementForcesAndSourcesCore::UserDataOperator::getEdgeFE() {
211  return static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE);
212 }
213 
215 EdgeElementForcesAndSourcesCore::UserDataOperator::getFTensor1Direction() {
216  return FTensor::Tensor1<double, 3>(getDirection()[0], getDirection()[1],
217  getDirection()[2]);
218 }
219 
221 EdgeElementForcesAndSourcesCore::UserDataOperator::getFTensor1Coords() {
222  double *ptr = &*getCoords().data().begin();
223  return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3>(ptr, &ptr[1],
224  &ptr[2]);
225 }
226 
227 template <>
229 EdgeElementForcesAndSourcesCore::UserDataOperator::getFTensor1TangentAtGaussPts<
230  3>() {
231  double *ptr = &*getTangentAtGaussPts().data().begin();
232  return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 3>(ptr, &ptr[1],
233  &ptr[2]);
234 }
235 
236 template <>
238 EdgeElementForcesAndSourcesCore::UserDataOperator::getFTensor1TangentAtGaussPts<
239  2>() {
240  double *ptr = &*getTangentAtGaussPts().data().begin();
241  return FTensor::Tensor1<FTensor::PackPtr<double *, 3>, 2>(ptr, &ptr[1]);
242 }
243 
244 auto EdgeElementForcesAndSourcesCore::UserDataOperator::getFTensor1Normal(
245  const FTensor::Tensor1<double, 3> &vec) {
250  auto t_dir = getFTensor1Direction();
251  t_normal(i) = FTensor::levi_civita(i, j, k) * t_dir(j) * vec(k);
252  return t_normal;
253 }
254 
255 auto EdgeElementForcesAndSourcesCore::UserDataOperator::
256  getFTensor1Normal() {
257  return getFTensor1Normal(tFaceOrientation);
258 }
259 
260 auto EdgeElementForcesAndSourcesCore::UserDataOperator::
261  getFTensor1NormalsAtGaussPts(const FTensor::Tensor1<double, 3> &vec) {
265 
266  auto &normals_at_gauss_pts =
267  static_cast<EdgeElementForcesAndSourcesCore *>(ptrFE)->normalsAtGaussPts;
268  normals_at_gauss_pts.resize(3, getGaussPts().size2(), false);
269  auto t_normal = getFTensor1FromMat<3>(normals_at_gauss_pts);
270  auto t_dir = getFTensor1TangentAtGaussPts<3>();
271  for(auto gg = 0; gg!=getGaussPts().size2(); ++gg) {
272  t_normal(i) = FTensor::levi_civita(i, j, k) * t_dir(j) * vec(k);
273  ++t_normal;
274  ++t_dir;
275  }
276 
277  return getFTensor1FromMat<3>(normals_at_gauss_pts);
278 }
279 
280 auto EdgeElementForcesAndSourcesCore::UserDataOperator::
281  getFTensor1NormalsAtGaussPts() {
282  return getFTensor1NormalsAtGaussPts(tFaceOrientation);
283 }
284 
285 /**
286  * @deprecated use EdgeElementForcesAndSourcesCore
287  */
290 
291 } // namespace MoFEM
292 
293 #endif //__EDGEELEMENTFORCESANDSURCESCORE_HPP__
294 
295 /**
296  * \defgroup mofem_forces_and_sources_edge_element Edge Element
297  *
298  * \brief Implementation of edge element.
299  *
300  * \ingroup mofem_forces_and_sources
301  */
UBlasMatrix< double >
MoFEM::EdgeElementForcesAndSourcesCore
Edge finite element.
Definition: EdgeElementForcesAndSourcesCore.hpp:30
MoFEM::EdgeElementForcesAndSourcesCore::UserDataOperator::getTensor1Coords
DEPRECATED FTensor::Tensor1< FTensor::PackPtr< double *, 3 >, 3 > getTensor1Coords()
Definition: EdgeElementForcesAndSourcesCore.hpp:163
MoFEM::EdgeElementForcesAndSourcesCoreBase
DEPRECATED typedef EdgeElementForcesAndSourcesCore EdgeElementForcesAndSourcesCoreBase
Definition: EdgeElementForcesAndSourcesCore.hpp:289
DEPRECATED
#define DEPRECATED
Definition: definitions.h:17
FTensor::Tensor1< double, 3 >
EntityHandle
MoFEM::Exceptions::MoFEMErrorCode
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:56
MoFEM::EdgeElementForcesAndSourcesCore::tangentAtGaussPts
MatrixDouble tangentAtGaussPts
Definition: EdgeElementForcesAndSourcesCore.hpp:48
MoFEM::EdgeElementForcesAndSourcesCore::cOnn
const EntityHandle * cOnn
Definition: EdgeElementForcesAndSourcesCore.hpp:54
FTensor::levi_civita
constexpr std::enable_if<(Dim0<=2 &&Dim1<=2), Tensor2_Expr< Levi_Civita< T >, T, Dim0, Dim1, i, j > >::type levi_civita(const Index< i, Dim0 > &, const Index< j, Dim1 > &)
levi_civita functions to make for easy adhoc use
Definition: Levi_Civita.hpp:617
MoFEM::EdgeElementForcesAndSourcesCore::lEngth
double & lEngth
Definition: EdgeElementForcesAndSourcesCore.hpp:51
MoFEM::DeprecatedCoreInterface
Deprecated interface functions.
Definition: DeprecatedCoreInterface.hpp:16
MoFEM::EdgeElementForcesAndSourcesCore::numNodes
int numNodes
Definition: EdgeElementForcesAndSourcesCore.hpp:53
FTensor::Tensor1::data
T data[Tensor_Dim]
Definition: Tensor1_value.hpp:10
MoFEM
implementation of Data Operators for Forces and Sources
Definition: Common.hpp:10
MoFEM::ForcesAndSourcesCore::UserDataOperator
Definition: ForcesAndSourcesCore.hpp:549
MoFEM::EdgeElementForcesAndSourcesCore::Switches
Switches
Definition: EdgeElementForcesAndSourcesCore.hpp:40
MoFEM::EdgeElementForcesAndSourcesCore::cOords
VectorDouble cOords
Definition: EdgeElementForcesAndSourcesCore.hpp:56
i
FTensor::Index< 'i', SPACE_DIM > i
Definition: hcurl_divergence_operator_2d.cpp:27
MoFEM::EdgeElementForcesAndSourcesCore::normalsAtGaussPts
MatrixDouble normalsAtGaussPts
Definition: EdgeElementForcesAndSourcesCore.hpp:49
FTensor::Index< 'i', 3 >
MoFEM::EdgeElementForcesAndSourcesCore::UserDataOperator
default operator for EDGE element
Definition: EdgeElementForcesAndSourcesCore.hpp:68
MoFEM::EdgeElementForcesAndSourcesCore::tFaceOrientation
static FTensor::Tensor1< double, 3 > tFaceOrientation
Definition: EdgeElementForcesAndSourcesCore.hpp:44
MoFEM::ForcesAndSourcesCore
structure to get information form mofem into EntitiesFieldData
Definition: ForcesAndSourcesCore.hpp:22
MoFEM::EdgeElementForcesAndSourcesCore::UserDataOperator::getTangetAtGaussPts
DEPRECATED MatrixDouble & getTangetAtGaussPts()
Definition: EdgeElementForcesAndSourcesCore.hpp:130
UserDataOperator
ForcesAndSourcesCore::UserDataOperator UserDataOperator
Definition: HookeElement.hpp:75
j
FTensor::Index< 'j', 3 > j
Definition: matrix_function.cpp:19
UBlasVector< double >
MoFEM::EdgeElementForcesAndSourcesCore::meshPositionsFieldName
std::string meshPositionsFieldName
Definition: EdgeElementForcesAndSourcesCore.hpp:33
k
FTensor::Index< 'k', 3 > k
Definition: matrix_function.cpp:20
MoFEM::EdgeElementForcesAndSourcesCore::dIrection
VectorDouble dIrection
Definition: EdgeElementForcesAndSourcesCore.hpp:55
MoFEM::FaceElementForcesAndSourcesCoreOnSide
Base face element used to integrate on skeleton.
Definition: FaceElementForcesAndSourcesCoreOnSide.hpp:19